0 引言

在对注聚参数进行优化设计时，传统的做法通常是按照各注聚参数的不同取值设计一系列组合方案，然后利用油藏数值模拟方法对所有组合方案进行筛选^[1-3]，但这种做法工作量巨大，需要耗费很长时间^[4-6]。油藏生产优化经历了开发方式由简单到复杂、控制变量由单一到多元、优化解法由梯度算法到无梯度算法的发展过程，其本身也经历了由单次优化到闭合优化的转变^[7-9]。无梯度优化算法^[10-12]对最优化算法理论研究很有意义，被应用于许多实际问题的解决，在很多行业都得到了广泛应用。赵辉等^[13-14]、万琦^[15]、张凯等^[16]将无梯度算法应用于油藏开发生产优化。油藏开发生产优化是通过求解描述油藏生产的最优控制模型，优化油水井的产出和注入控制参数^[17-18]，使油藏开发处于最优状态，该问题属于复杂的大系统动态最优控制问题，可以用无梯度优化算法解决。

针对有限差分随机逼近（FDSA）算法和同时扰动随机逼近（SPSA）算法的原理和优缺点，提出一种有限差分梯度引导的同时扰动随机逼近（SPSA-FDG）算法，并将其应用到聚合物驱注采参数优化中，以期为化学驱注采参数的最优控制提供新的思路。

1 聚合物驱生产优化数学模型

在实际生产中，聚合物是分段塞注入的，因此，聚合物驱最优控制即通过最优化方法和油藏数值模拟技术确定每个段塞的尺寸和最优注入质量浓度，以实现油田开发效果的最佳化。由于聚合物驱投资大、风险高，因此，通常采用增量累计净现值作为性能指标对该类项目的盈利水平进行科学评价，其表达式为

$ J = \sum\limits_{i = 1}^t {\left( {C_{{\rm{in}}}^i - C_{{\rm{out}}}^i} \right){{\left( {1 + {i_{\rm{c}}}} \right)}^{ - i}}} $

(1)

式中：J为增量累计净现值，元；C_inⁱ为第i年增量现金流入量，元；C_outⁱ为第i年增量现金流出量，元；C_inⁱ-C_outⁱ为第i年的增量净现金流量，元；t为评价时间，a；i_c为基准收益率，%。

增量年现金流入量为聚合物驱相对于水驱年增原油产量的销售收入，增量年现金流出量为聚合物驱生产过程中的投资、税费以及生产成本。选取聚合物驱项目的增量累计净现值作为最优控制的性能指标函数，同时选取聚合物注入质量浓度和段塞尺寸作为聚合物驱最优控制研究的优化参数。在实际生产中，参数的选取会受到注入能力和聚合物价格的影响，故对聚合物注入质量浓度和段塞尺寸要加以限制，即控制变量要满足一定的约束条件，而边界约束是最常见的约束形式^[19-21]。结合最优化问题的一般数学表达形式，建立聚合物驱最优控制数学模型，即

$ \begin{array}{l} \max J = \sum\limits_{i = 1}^t {\left\{ {{Q_{{\rm{o}}i}}{P_{\rm{o}}}\alpha - \left[ {{Q_{{\rm{o}}i}}{C_{\rm{m}}} + {Q_{{\rm{p}}i}}{P_{\rm{p}}} + } \right.} \right.} \\ \;\;\;\;\;\;\;\;\;\;\;\;\left. {\left. {\left( {{R_{\rm{s}}} + R{P_{\rm{o}}}} \right){Q_{{\rm{o}}i}}\alpha } \right]} \right\} \times {\left( {1 + {i_{\rm{c}}}} \right)^{ - i}} - n{I_{\rm{s}}}\\ \left( {{\rm{s}}.{\rm{t}}.\;\;\;u_i^{{\rm{low}}} \le {u_i} \le u_i^{{\rm{up}}},i = 1,2, \cdots ,{N_u}} \right) \end{array} $

(2)

式中：Q_oi为年增油量，t；P_o为原油价格，元/t；α为原油商品率，%；Q_pi为年注聚量，t；C_m为吨油增量操作费用，元；P_p为聚合物价格，元/t；R_s为资源税，元/t；R为综合税率，%；n为区块等效井数，口；I_s为单井增量投资费用，元/口；u_i为第i个控制变量；u_i^up为控制变量u_i的上约束；u_i^low为控制变量u_i的下约束；N_u为控制变量数，个。

聚合物驱最优控制问题就是在控制变量满足约束条件的情况下，计算性能指标函数J的最大值及相应的最优控制变量u。

2 聚合物驱无梯度优化求解算法

聚合物驱注采参数最优控制是一类复杂的大规模最优化问题，目前尚无有效方法计算聚合物驱的真实梯度，因此，采用无梯度优化算法求解聚合物驱最优控制模型。有限差分随机逼近算法梯度近似程度高、收敛速度快，但每个迭代步的性能指标函数计算量正比于控制变量维数，当控制变量维数较高时，计算量非常大。SPSA-FDG算法对控制变量进行同时扰动，每个迭代步最少仅须计算一次性能指标函数即可获得迭代搜索方向，减少了每个迭代步的计算量，但同时该算法计算的梯度近似精度低，降低了收敛速度^{[16, 22]}。

针对上述问题，提出了一种由有限差分梯度引导的同时扰动随机逼近算法，该方法在迭代过程中根据性能指标函数对各控制变量的敏感性变化，实时调整各控制变量扰动步长的大小比例，以便使各控制变量所引起的性能指标函数改变量近似相等，保证计算得到的各控制变量SPSA近似梯度大小比例与FDSA算法类似，从而提高收敛速度。

2.1 有限差分随机逼近算法原理

基于有限差分随机逼近算法，利用中心差分格式可以计算得到性能指标函数对控制变量的近似梯度为

$ \begin{array}{*{20}{c}} {\hat d_i^l\left( {u_{{\rm{opt}}}^l} \right) = \frac{{J\left( {u_{{\rm{opt}}}^l + {\varepsilon _{l,i}}{\Delta _{l,i}}} \right) - J\left( {u_{{\rm{opt}}}^l - {\varepsilon _{l,i}}{\Delta _{l,i}}} \right)}}{{2{\varepsilon _{l,i}}{\Delta _{l,i}}}}}\\ {\left( {i = 1,2, \cdots ,{N_u}} \right)} \end{array} $

(3)

式中：l为第l个迭代计算步；u_opt^l为第l个迭代计算步所获得的N_u维最优控制变量列向量，u_opt^l=[u_{opt, 1}^l u_{opt, 2}^l, …, u_{opt, Nu}^l]；ε_{l, i}为第l个迭代计算步第i个控制变量的扰动步长；Δ_{l, i}为第l个迭代计算步第i个控制变量的随机扰动值，取值为1或-1。

2.2 同时扰动随机逼近算法原理

同时扰动随机逼近算法独立于具体的油藏数值模拟器，而且能够对所有的优化控制变量进行同时扰动优化，计算目标函数对所有控制变量的近似梯度信息时仅须进行2次油藏数值模拟计算，而且计算次数不随控制变量的增加而增多，因此，该方法尤其适合求解控制变量参数体系复杂的海上油田化学驱生产优化问题。

同时扰动随机逼近算法在第l个迭代计算步，目标函数J在u_opt^l处的SPSA梯度表达式为

$ \begin{array}{*{20}{c}} {\hat g_i^l\left( {u_{{\rm{opt}}}^l} \right) = \frac{{J\left( {u_{{\rm{opt}}}^l + {\varepsilon _l}{\Delta _l}} \right) - J\left( {u_{{\rm{opt}}}^l} \right)}}{{{\varepsilon _{l,i}}{\Delta _{l,i}}}}}\\ {\left( {i = 1,2, \cdots ,{N_u}} \right)} \end{array} $

(4)

式中：ε_l为第l个迭代计算步N_u维控制变量对应的扰动步长对角矩阵，其中ε_{l, i}为第l个迭代计算步第i个控制变量的扰动步长，ε_l=$ \left[{\begin{array}{*{20}{c}} {{\varepsilon _{l, 1}}}&{}&{}&{}&{}\\ {}& \ddots &{}&{}&{}\\ {}&{}&{{\varepsilon _{l, \mathit{i}}}}&{}&{}\\ {}&{}&{}& \ddots &{}\\ {}&{}&{}&{}&{{\varepsilon _{l, {\mathit{N}_\mathit{u}}}}} \end{array}} \right]$；u_opt^l为在第l个迭代计算步所获得的最优控制变量；Δ_l为N_u维随机扰动向量，其中所包含元素Δ_{l, i}(i=1, 2, …, N_u)为服从参数为±1的对称Bernoulli分布。

在获得随机扰动梯度后，便可采用迭代法进行优化求解，在第l+1迭代计算步所获得的控制变量为

$ u_{{\rm{opt}}}^{l + 1} = u_{{\rm{opt}}}^l + {\alpha _l}{{\hat g}^l}\left( {u_{{\rm{opt}}}^l} \right) $

(5)

式中：α_l为迭代步长，为了保证当l→∞时，u_opt^l能够收敛到局部最优解。

2.3 有限差分梯度引导的同时扰动随机逼近算法

由式（2）可知，性能指标函数对控制变量u_i和u_j的FDSA近似梯度比值为

$ \begin{array}{*{20}{c}} {\frac{{\hat d_i^l\left( {u_{{\rm{opt}}}^l} \right)}}{{\hat d_j^l\left( {u_{{\rm{opt}}}^l} \right)}} = \frac{{J\left( {u_{{\rm{opt}}}^l + {\varepsilon _{l,i}}{\Delta _{l,i}}} \right) - J\left( {u_{{\rm{opt}}}^l} \right)}}{{J\left( {u_{{\rm{opt}}}^l + {\varepsilon _{l,j}}{\Delta _{l,j}}} \right) - J\left( {u_{{\rm{opt}}}^l} \right)}} \cdot \frac{{{\varepsilon _{l,j}}{\Delta _{l,j}}}}{{{\varepsilon _{l,i}}{\Delta _{l,i}}}}}\\ {\left( {i,j = 1,2, \cdots ,{N_u}} \right)} \end{array} $

(6)

假设SPSA算法中各控制变量采用的扰动步长取值与FDSA算法相同，则性能指标函数对控制变量u_i的SPSA近似梯度与FDSA近似梯度之间有如下关系

$ \begin{array}{*{20}{c}} {\hat g_i^l\left( {u_{{\rm{opt}}}^l} \right) = \frac{{J\left( {u_{{\rm{opt}}}^l + {\varepsilon _l}{\Delta _l}} \right) - J\left( {u_{{\rm{opt}}}^l} \right)}}{{J\left( {u_{{\rm{opt}}}^l + {\varepsilon _{l,i}}{\Delta _{l,i}}} \right) - J\left( {u_{{\rm{opt}}}^l} \right)}} \cdot \hat d_i^l\left( {u_{{\rm{opt}}}^l} \right)}\\ {\left( {i = 1,2, \cdots ,{N_u}} \right)} \end{array} $

(7)

由式（7）可知，性能指标函数对控制变量u_i和u_j的SPSA近似梯度比值为

$ \begin{array}{*{20}{c}} {\frac{{\hat g_i^l\left( {u_{{\rm{opt}}}^l} \right)}}{{\hat g_j^l\left( {u_{{\rm{opt}}}^l} \right)}} = \frac{{J\left( {u_{{\rm{opt}}}^l + {\varepsilon _{l,j}}{\Delta _{l,j}}} \right) - J\left( {u_{{\rm{opt}}}^l} \right)}}{{J\left( {u_{{\rm{opt}}}^l + {\varepsilon _{l,i}}{\Delta _{l,i}}} \right) - J\left( {u_{{\rm{opt}}}^l} \right)}} \cdot \frac{{\hat d_i^l\left( {u_{{\rm{opt}}}^l} \right)}}{{\hat d_j^l\left( {u_{{\rm{opt}}}^l} \right)}}}\\ {\left( {i,j = 1,2, \cdots ,{N_u}} \right)} \end{array} $

(8)

若要使性能指标函数对控制变量u_i和u_j的SPSA近似梯度大小比例与FDSA近似梯度大小比例相同，则要满足

$ \begin{array}{*{20}{c}} {\frac{{J\left( {u_{{\rm{opt}}}^l + {\varepsilon _{l,j}}{\Delta _{l,j}}} \right) - J\left( {u_{{\rm{opt}}}^l} \right)}}{{J\left( {u_{{\rm{opt}}}^l + {\varepsilon _{l,i}}{\Delta _{l,i}}} \right) - J\left( {u_{{\rm{opt}}}^l} \right)}} = 1}\\ {\left( {i,j = 1,2, \cdots ,{N_u}} \right)} \end{array} $

(9)

由于控制变量扰动步长取值为非负数，将式（9）带入式（6）可得，控制变量u_i和u_j的扰动步长应满足

$ \begin{array}{*{20}{c}} {\frac{{{\varepsilon _{l,i}}}}{{{\varepsilon _{l,j}}}} = \frac{{\hat d_j^l\left( {u_{{\rm{opt}}}^l} \right){\Delta _{l,j}}}}{{\hat d_i^l\left( {u_{{\rm{opt}}}^l} \right){\Delta _{l,i}}}} = \left| {\frac{{\hat d_j^l\left( {u_{{\rm{opt}}}^l} \right)}}{{\hat d_i^l\left( {u_{{\rm{opt}}}^l} \right)}}} \right|}\\ {\left( {i,j = 1,2, \cdots ,{N_u}} \right)} \end{array} $

(10)

在第l个迭代计算步建立SPSA扰动步长修正矩阵为

$ {M_l} = \hat d_{\max }^l\left( {u_{{\rm{opt}}}^l} \right)\left[ {\begin{array}{*{20}{c}} {\frac{1}{{\left| {\hat d_1^l\left( {u_{{\rm{opt}}}^l} \right)} \right|}}}&{}&{}&{}&{}\\ {}& \ddots &{}&{}&{}\\ {}&{}&{\frac{1}{{\left| {\hat d_i^l\left( {u_{{\rm{opt}}}^l} \right)} \right|}}}&{}&{}\\ {}&{}&{}& \ddots &{}\\ {}&{}&{}&{}&{\frac{1}{{\left| {\hat d_{{N_u}}^l\left( {u_{{\rm{opt}}}^l} \right)} \right|}}} \end{array}} \right] $

(11)

式中：$ {\hat d^l}\left( {u_{{\rm{opt}}}^l} \right) = {\rm{max}}\left[{\left| {\hat d_i^l\left( {u_{{\rm{opt}}}^l} \right)} \right|} \right]$为第l个迭代计算步各控制变量FDSA近似梯度绝对值的最大值。

则第l个迭代计算步各控制变量的SPSA-FDG近似梯度为

$ \begin{array}{*{20}{c}} {\hat f_i^l\left( {u_{{\rm{opt}}}^l} \right) = \frac{{J\left( {u_{{\rm{opt}}}^l + {\varepsilon _{l,\max }}{M_l}{\Delta _l}} \right) - J\left( {u_{{\rm{opt}}}^l} \right)}}{{{\varepsilon _{l,\max }}{M_{l,i}}{\Delta _{l,i}}}}}\\ {\left( {i = 1,2, \cdots ,{N_u}} \right)} \end{array} $

(12)

式中：ε_{l, max}为根据能够取得$ \hat d_{{\rm{max}}}^l\left( {u_{{\rm{opt}}}^l} \right)$的控制变量所选取的适当扰动步长，常数；M_{l, i}为第l个迭代计算步第i个控制变量的扰动步长修正系数。

由式（12）可知，当采用SPSA-FDG算法时，如果每个迭代计算步都进行扰动步长修正，则类似于FDSA算法；如果从不进行扰动步长修正，则为SPSA算法。由此可见，SPSA-FDG算法是该类算法中相对一般的算法，FDSA算法和SPSA算法则为该类算法的2个特例。

3 聚合物驱生产优化应用实例

基于CMG油藏数值模拟软件建立了油藏数值模拟模型，再利用SPSA-FDG算法开展了聚合物驱最优控制研究，并将其优化结果与FDSA算法和SPSA算法的优化结果进行了对比。

所建油藏数值模拟模型划分为7个模拟层（图 1），采用五点法井网，共包括4口注入井和9口生产井。各注采井组纵向非均质性差异大，其中INJ-1井和INJ-3井所在井组纵向非均质性弱，变异系数为0.18，INJ-2井和INJ-4井所在井组纵向非均质性强，变异系数为0.86。油藏的平均渗透率为1 500 mD，地下原油黏度为60 mPa·s，孔隙体积为116.89万m³，聚合物溶液的残余阻力系数为1.8，可及孔隙体积为0.82，最大吸附体积质量为0.3 kg/m³。

聚合物驱过程采用四段塞的注入方式，包括前置段塞、主段塞、副主段塞和后置段塞，本次聚合物驱最优控制研究即通过最优化方法确定每口井每个段塞的最佳注入质量浓度和段塞尺寸，共包括32个控制变量。在优化过程中，对各控制变量仅考虑边界约束条件，其中注入质量浓度考虑矿场聚合物驱实践经验和聚合物实际注入能力，上、下边界分别设置为0 mg/L和2 300 mg/L，各阶段段塞尺寸上、下边界分别设置为0 PV和0.1 PV。

分别采用FDSA，SPSA和SPSA-FDG等3种算法求解聚合物驱最优控制问题，均基于相同的初始估计。4口注聚井的初始状态估计相同，其中前置段塞质量浓度和尺寸分别为1 700 mg/L和0.012 5 PV；主段塞质量浓度和尺寸分别为1 400 mg/L和0.037 5 PV；副主段塞质量浓度和尺寸分别为1 200 mg/L和0.037 5 PV；后置段塞质量浓度和尺寸分别为1 000 mg/L和0.012 5 PV。对于SPSA-FDG算法，当基本SPSA算法进行到第13，30，50次迭代计算时，分别引入FDG梯度对扰动步长进行修正。

从各优化方法在不同迭代次数或数值模拟次数下的净现值变化（图 2、表 1）可以看出，优化后最终净现值均显著增加。虽然SPSA-FDG算法梯度估计精度没有FDSA算法高，但却高于SPSA算法，同时又保留了SPSA算法每次迭代仅须进行2次油藏数值模拟计算的优点，因此，在达到收敛时仅须进行248次油藏数值模拟计算，是3种方法中总体上计算量最小的，且收敛到的净现值也略高于其他2种方法。综上所述，SPSA-FDG算法收敛性要好于FDSA算法和SPSA算法，而且易于和任何油藏数值模拟器相结合，能够用来求解聚合物驱最优控制问题。

图 3为利用SPSA-FDG算法优化后的聚合物驱调控图。从图 3可以看出，各注入井聚合物驱前置段塞聚合物质量浓度最高，用量中等，主要起到前缘调剖作用；主段塞及副主段塞聚合物质量浓度较高，用量最多，是聚合物驱改善油水流度比、降水增油的关键部分；后置段塞为流度保护段塞，聚合物质量浓度最低，用量最少，主要是为了防止后续水驱指进破坏主段塞而影响开发效果。

对利用各种优化算法计算得到的各注入井聚合物配注量优化结果（表 2）进行分析可知，在聚合物注入总量一定的情况下，注入井所在井组非均质性强，聚合物配注量大；相反，若注入井所在井组非均质性弱，则聚合物配注量小。为了进一步讨论各单井聚合物优化配注对聚合物驱总体开发效果的影响，设计了各单井均匀注聚方案。均匀注聚方案的聚合物注入总量与SPSA-FDG算法优化结果相同，为1 123 PV·mg/L，采用单段塞方式注入，段塞尺寸取SPSA-FDG算法优化得到的各单井段塞尺寸的平均值，即0.652 PV，则聚合物注入质量浓度为1 722.4 mg/L。

图 4为分别采用SPSA-FDG算法优化注聚方案与均匀注聚方案计算得到的油藏含水率及累积采出程度变化曲线。从图 4可看出，采用由SPSA-FDG算法计算得到的优化方案使开发效果得到明显改善。这主要是因为经SPSA-FDG算法优化后，位于强非均质性井组的注入井聚合物配注量增大，且由于聚合物驱提高采收率幅度随着油藏非均质性的增强呈上升趋势，因此，非均质性较强的井组因聚合物配注量增大而增加的产油量大于非均质性较弱的井组因相同聚合物配注量减小而减少的产油量，从而使得油藏总体采出程度得到提高。

图 5为分别利用SPSA-FDG算法优化注聚方案与均匀注聚方案计算得到的第4层油藏模型剩余油饱和度分布。从图 5可以看出，采用由SPSA-FDG算法计算得到的优化方案进行聚合物驱开发时，在聚合物注入总量一定的情况下，油藏波及系数增大，强非均质性井组的剩余油饱和度明显降低，弱非均质性井组的剩余油饱和度变化不大，油藏剩余油分布更加均匀。这是由于重力作用和油层的正韵律性使得注入水向油层下部高渗层段窜流，以致随着纵向变异系数增大而产生暴性水淹，导致水驱采收率急剧下降。

另外，对采用SPSA-FDG算法优化注聚方案与均匀注聚方案计算得到的净现值随开发时间变化的关系曲线进行对比和分析发现，在聚合物注入总量一定的情况下，采用SPSA-FDG算法优化注聚方案，经济净现值比均匀注聚方案提高了11.64%，且投资回收期比均匀注聚方案缩短。

4 结论

（1）本文提出的由有限差分梯度引导的同时扰动随机逼近（SPSA-FDG）算法与有限差分随机逼近（FDSA）算法相比，数值模拟次数可以减少一半以上，与同时扰动随机逼近（SPSA）算法相比，迭代次数可以减少一半以上。

（2）注入井聚合物驱前置段塞聚合物质量浓度最高，用量中等，主要起到前缘调剖作用；主段塞及副主段塞聚合物质量浓度较高，用量最多，是聚合物驱改善油水流度比、降水增油的关键部分；后置段塞为流度保护段塞，聚合物质量浓度最低，用量最少，主要是为了防止后续水驱指进破坏主段塞而影响开发效果。

（3）与均匀注聚方案对比发现，采用SPSA-FDG算法计算得到的优化方案，聚合物驱见效时间提前，含水下降漏斗加深，累积采出程度提高，油藏总体开发效果更好。

（4）SPSA-FDG算法能够解决海上油田聚合物驱开发阶段的注采最优控制问题，可以最大限度地提高聚合物驱开发阶段的经济效益，从而改善海上油田不同开发阶段开发方式组合模式的总体开发效果。

致谢: 在项目完成过程中，中海油研究总院教授级高级工程师康晓东和中国石油大学（华东）侯健教授、杜庆军老师均给予了帮助，在此表示感谢！