张霖森, 程 兰, 张守贵
(重庆师范大学 数学科学学院, 重庆 401331)
迄今为止,已有很多关于各种变分不等式问题的研究[1-6].但相对而言,关于曲率障碍的四阶椭圆变分不等式问题的研究还不多[7-9].在这类障碍问题中,弹性膜在边界固定和垂直方向受到外力作用与曲率障碍条件下,需要确定弹性膜的平衡位置.对于这类问题的数值解法,有间断Galerkin方法[1]、有限元方法[3,7-8]、有限差分法[6].交替方向乘子法(ADMM)在结构优化问题中有着广泛的应用,例如二维的变分不等式[10]、接触问题[11-12]和Stokes问题[13].ADMM的每一次迭代,只需要求解一个线性问题,而且辅助未知量和Lagrange乘子是显式计算的.对于任意的正参数,ADMM都是全局收敛的.但是该方法对罚参数非常敏感,很难根据具体问题选择合适的罚参数.
本文重点分析了ADMM和罚参数的自适应法则求解关于单侧曲率障碍四阶变分不等式的组合算法[14-15].首先将ADMM应用于求解曲率障碍问题,由于ADMM的收敛速度严重依赖于罚参数,为了改进算法性能,我们提出了一种自适应法则,通过平衡原理和迭代函数来选择合适的罚参数[16-17],然后得到了一种自适应交替方向乘子算法(SADMM),并分析了该方法的收敛性.最后给出了一些数值结果来验证该方法的可行性和有效性.
本文考虑以下四阶椭圆变分不等式问题:给定f∈L2(Ω),求解u∈K,使得
(1)
其中Ω是平面区域R2中的有界闭区域,其边界为Г=∂Ω,闭凸集
(2)
定义以下符号:
则问题(1)可以写成如下变分不等式:求解u∈K,使得
a(u,v-u)≥l(v-u), ∀v∈K.
(3)
由式(3)可知式(1)是一个变分不等式问题.则式(3)有如下性质:
其中C>0;
由于双线性形式a(·,·)是对称的,故变分不等式问题(1)等价于如下变分形式的极小值问题:求解u∈K,使得
J(u)≤J(v), ∀v∈K,
(4)
其中
(5)
由式(1)、(4)和文献[8]定义如下Lagrange函数L:H2(Ω)×(L2(Ω))2→R,
(6)
和函数空间
(7)
(8)
(9)
L(u,μ)≤L(u,λ)≤L(v,λ).
(10)
由式(10)中L(u,μ)≤L(u,λ)直接得到的结论(8),可以推出
Δu≤0, a.e. inΩ,
(11)
(12)
则
L(u,λ)=J(u).
(13)
在此基础上,考虑增广Lagrange乘子法.我们观察到,问题(1)和(4)都等价于:求解{u,p}∈H,使得
j(u,p)≤j(v,q), ∀{v,q}∈H,
(14)
其中
(15)
(16)
对于ρ>0,由增广Lagrange乘子法
定义
(17)
接着定义
(18)
求解pn+1∈L2(Ω),使得
Lρ(un,pn+1,λn)≤Lρ(un,q,λn), ∀q∈L2(Ω).
(19)
(20)
更新Lagrange乘子,得到
λn+1=λn+ρ(Δun+1-pn+1).
(21)
在这种方法中,最小化问题(19)中的pn+1可以由un和λn显式求解.最小化问题(20)形成一个变分问题,对于给定的λn,pn+1和ρ>0,该变分问题具有唯一的解.因此,我们得到下面的ADMM[10,16].
算法1(ADMM)
第二步 计算辅助变量pn+1∈L2(Ω),
(22)
(23)
第四步 更新Lagrange乘子
λn+1=λn+ρ(Δun+1-pn+1).
(24)
对于算法1,我们注意到它对于任何固定参数ρ>0都是无条件收敛的.然而,如果参数太小或太大,该方法的效率将大大降低.为了改进ADMM的性能,我们提出了一个可变参数ρn的自适应规则.下面我们假设一个非负序列{sn}满足
(25)
在此基础上,我们提出SADMM,具体算法过程如下[11].
算法2(SADMM)
第二步 计算辅助变量pn+1∈L2(Ω),
(26)
(27)
第四步 更新Lagrange乘子
λn+1=λn+ρn(Δun+1-pn+1).
(28)
第五步 选取罚参数,使得
(29)
关于算法2中第五步罚参数的选取将在算例分析部分详细说明.为了证明算法2的收敛性,需要如下的基本结论.
定义以下符号:
根据问题的性质和算法的原理,可得算法2的收敛性结果.
定理1 设{un,pn,λn}是算法2产生的序列,则
(30)
其中p=Δu.
证明令δun=un-u,δΔun=Δun-Δu,δρn=ρn-ρ,δλn=λn-λ.由式(9)和(23),有
在以上两式中取v=δΔun+1,将它们相减得到
a(δun+1,δun+1)=(δλn+ρn(Δun+1-pn+1),-δΔun+1),
(31)
在式(31)中,用p=Δu,我们得到
a(δun+1,δun+1)=(δλn+ρn(δΔun+1-δpn+1),-δΔun+1).
(32)
根据式(28),有
δλn+1=δλn+ρn(δΔun+1-δpn+1).
(33)
从引理1的式(11)和(12),有λ≥0,-Δu≥0,且在Ω上(-Δu,λ)=0,因此可以得到
ρn(-Δu,(λn+ρnΔun)+-λ)≥0,
(34)
又有
(λn+ρnΔun-(λn+ρnΔun)+,(λn+ρnΔun)+-λ)≥0.
(35)
结合式(34)和(35),有
(λn+ρnδΔun-(λn+ρnΔun)+,(λn+ρnΔun)+-λ)≥0.
再根据式(26),可得到
(-ρn(Δun-pn+1)+ρnδΔun,δλn+ρn(Δun-pn+1))≥0,
因此
(36)
根据式(28),可得到
δλn+1-ρnδΔun+1=δλn-ρnδΔun+ρn(Δun-pn+1).
最后,由式(36)得到
因此,该定理得证.
证明根据sn≥0,0≤ρn+1≤(1+sn)ρn和式(33),有
(37)
将式(30)代入式(37),有
(38)
由式(33)和(38),可得
(39)
即存在一个常数C>0,使得对于任意的n≥0,有
(40)
根据式(38),我们也可以得出
(41)
由式(33)、(40)和(41),有
因此
(42)
(43)
因为p=Δu,再结合式(43),可以得到
(44)
因此,在L2(Ω)上,pn→p.
对于∀n≥1,由式(39),有
(45)
根据式(33)和(45)可知{λn}是有界的,可得λn→λ.
注1 在算法2的第五步中令sn=0,有ρn=ρ0,则算法2中罚参数为固定参数,即为算法1,从而可以得算法1也是收敛的.
本文利用自适应交替方向乘子算法求解四阶单侧曲率障碍的变分不等式问题,迭代过程中通过迭代函数自动调整罚参数,用变参数ρn代替ρ,从而达到提高算法效率的目的[2,10,13].下面具体考虑算法2(SADMM)中选择罚参数ρn的自适应法则.根据算法的收敛性证明式(30)知下式成立:
‖λn+1-λ-ρn(Δun+1-Δu)‖≤‖λn-λ-ρn(Δun-Δu)‖.
利用平衡加快收敛,我们希望
‖λn-λ‖≈‖ρn(Δun-Δu)‖.
使用以下方法来确定合适的参数ρn,给定常数μ>0,如果
‖λn-λ‖>(1+μ)‖ρn(Δun-Δu)‖,
那么在下次迭代的时候要增加ρn;如果
用cn计算ρn改变的次数,令c0=0,则
采用如下方法得到非负序列{sn}:
为了检验本文算法的可行性,这里根据具体的算例进行数值分析.考虑一个具体的障碍问题:求解u∈K,使得
采用本文的ADMM 和SADMM两种方法进行求解,取步长h=1/40和参数ρ=1,u和-Δu的数值解结果如图1和图2所示,解析解如图3所示,图4为逐点误差,可见数值解和精确解是吻合的.
图1 u的数值解 图2 -Δu的数值解
图3 精确解u 图4 逐点误差
对于不同的初始罚参数取值ρ,表1分别对步长h=1/10,1/20,1/40,1/80时ADMM和SADMM所需的迭代次数进行了比较,表中*表示迭代达到了最大迭代次数200后即停止.表2分别给出了迭代所需的CPU运行时间,单位为s.表1和表2的结果表明,SADMM不仅有效减少了迭代次数,收敛速度更快,而且很稳定.SADMM的收敛速度和CPU运行时间受初始参数ρ的影响很小.
表1 算法随步长变化所需迭代次数的情况
表2 算法随步长变化所需CPU时间情况
本文提出了求解四阶单侧曲率障碍变分不等式问题的SADMM.先引入辅助函数和Lagrange函数,将问题转化为鞍点问题,使用ADMM求解.为了提高计算效率,通过平衡原理和迭代函数,迭代过程使用自适应法则选取合适的罚参数.数值算例的结果表明,参数对SADMM影响较小,该方法加快了收敛速度,且较为稳定.该方法还可以推广到四阶双侧曲率障碍问题的数值求解中.
致谢本文作者衷心感谢重庆师范大学教学改革研究项目(xyjg007;sz202114)对本文的资助.