解Duffing方程周期解的Levenberg-Marquardt方法
2013-12-03姜志侠李延忠孟品超尹伟石
姜志侠, 李延忠, 李 军, 孟品超, 尹伟石
(长春理工大学 理学院应用数学系， 长春 130022)
Duffing方程是非线性理论中常用的微分方程, 其模型具有代表性. 关于Duffing方程解的存在性及数值计算方法的研究目前已有许多结果[1-5].
对于二阶微分方程
x″=f(t,x),
(1)
王怀忠等[1]证明了在满足限制共振条件
(2)
下, 方程(1)存在唯一的T周期解. 其中:f(t,x)=f(t+T,x)在 R×R → R上连续, 且关于x二次连续可微;N为一个非负整数.
对无强迫振动项的Duffing方程u″+g(u)=0, 刘淑媛等[3]利用山路引理证明了超线性Duffing方程周期解的存在性, 并给出一种求Duffing方程周期解的Mountain Pass算法； 李勇等[5]应用同伦算法, 构造性地证明了Poincare-Birkhoff定理, 使得用同伦算法求出周期解的数值实现成为可能.
1 问题的转化
设x(t,P)是方程(1)带有初始值
x(0)=P0,x′(0)=P1
(3)
x(T,P*)=x(0,P*),x′(T,P*)=x′(0,P*)
(4)
时, 微分方程(1)有一个T周期解. 因此, 对于最优化问题
(5)
如果存在满足J(P*)=0的一个全局最优解P*, 则方程(1)以P*为始点的解即为一个周期解. 其中X(T,P)=(x(T,P),x′(T,P))T, 且X(0,P)=(x(0,P),x′(0,P))T=P.
无约束最优化问题(5)的Karush-Kuhn-Tucker(KKT)一阶必要性条件[6]为
‖X(T,P*)-P*‖2=0,
(6)
由文献[7]可得:
定理1对满足条件(2)的微分方程(1), 如果存在Pk∈R(k=1,2,…), 使得
(7)
则当k→ +∞时,Pk→P*, 且J(P*)=minJ(Pk)=0.
因此, 求解微分方程(1)的周期解问题可以转化为求解相关的无约束最优化问题(5), 且以问题(5)最优解为始点的方程(1)的解是周期解. 冯子璇等[7]使用拟牛顿法求解转化后的最优化问题, 用得到的全局最优解作为微分方程(1)周期解的初始值得到了周期解.
在实际问题中, 常会碰到目标函数为非线性函数的平方和形式, 即最小二乘问题, Levenberg-Marquardt型方法是求解这类问题最有效的方法之一, 特别是对于残差为零或接近于零的问题. 针对奇异非线性方程组, 文献[8]给出了Levenberg-Marquardt方法的一种新的参数迭代方法, 并证明了在弱于非奇异性条件的局部误差有界条件下, Levenberg-Marquardt方法仍具有局部二次收敛速度. 本文提出一种基于Levenberg-Marquardt型算法求解微分方程的周期解算法.
2 最优化问题的求解方法
对问题(5), 记f(P)=X(T,P)-X(0,P),M(P)为f(P)的Jacobian矩阵, 令F(P)=‖X(T,P)-X(0,P)‖2/2,H(x)为F(P)的Hessian矩阵, 则可得:G(P)=2M(P)Tf(P),H(P)=2M(P)TM(P)+2Q(P), 其中:
随着P逐渐接近最优解, ‖X(T,P)-X(0,P)‖的值逐渐趋于零,Q(P)也趋于零. 因此, 当‖X(T,P)-X(0,P)‖在最优解处很小时, 使用Levenberg-Marquardt方法所用的搜索方向作为优化过程的搜索方向更有效.
Levenberg-Marquardt方法所用的搜索方向dk是一组线性等式的解： (M(Pk)TM(Pk)+λkI)dk=-M(Pk)f(Pk), 其中λk决定搜索方向和幅值大小, 只要λk足够大, 即可保证f(Pk+dk)
算法步骤如下:
1) 计算f(Pk),M(Pk),F(Pk)；
2) 计算搜索方向dk： (M(Pk)TM(Pk)+λkI)dk=-M(Pk)f(Pk);
3) 迭代Pk+1=Pk+λkdk: 如果‖M(Pk+1)‖<ε, 停止； 否则, 转4)；
4) 如果k=n, 则令P1=Pk+1, 返回1)； 否则, 转5)；
5) 按下式计算线性预测平方和Fl(Pk+1):fl(Pk+1)=M(Pk)Tdk+F(Pk),Fl(Pk+1)=fl(Pk+1)T×fl(Pk+1); 通过对F(Pk+1)和F(Pk)做三次内插值, 计算Fl(Pk+1)的极小值Fk(P*), 将最小化的估计步长记为α*;
6) 计算步长： 若Fl(Pk+1)>Fk(P*), 则令λk=λk/(1+α*)； 否则, 令λk=λk+(Fk(P*)-Fl(Pk))/α*, 再令k=k+1, 返回1).
3 算 例
考虑如下周期边值问题:
(8)
方程(8)满足限制共振条件(2), 因此有唯一周期解. 使用Levenberg-Marquardt方法, 在Matlab2011版软件中编程实现, 其中M(P)为f(P)的Jacobian矩阵, 由于f(P)=X(T,P)-X(0,P)没有连续的表达式, 故用差商代替M(P). 图1～图3为得到的周期解相关曲线. 表1列出了以(5.237 935 1,6.650 527 3)为初始点的数值计算结果.
图1 x与t的关系曲线 图2 x′与t的关系曲线 图3 x′与x的关系曲线 Fig.1 Related curve of x vs t Fig.2 Related curve of x′ vs t Fig.3 Related curve of x vs x′
表1 数值结果
续表1
由图1～图3及表1可见, 使用Levenberg-Marquardt方法求解Duffing方程的周期解有效, 本文算法利用最优化方法逐次求解了微分方程数值周期解的初始点, 避免了数值算法对初始点选择依赖性而产生的计算困难.
感谢吉林大学数学学院徐旭和韩月才教授的悉心指导.
[1] WANG Huai-zhong, LI Yong. Existence and Uniquesness of Periodic Solutions for Duffing Equations Across Many Point of Resonance [J]. Journal of Differential Equations, 1994, 108(1): 152-169.
[2] ZHU Wen-zhuang, WAN Jun, ZHANG Jin. Periodic Solution of a Class of Two Order Differential Equations [J]. Journal of Jilin University: Science Edition, 2004, 42(1)： 61-63. (祝文壮, 万军, 张锦. 一类二阶微分方程周期解的存在性 [J]. 吉林大学学报: 理学版, 2004, 42(1)： 61-63.)
[3] LIU Shu-yuan, LÜ Xian-rui, QI Yi. Numerical Mountain Pass Periodic Solutions of a Duffing Equation [J]. Journal of Jilin University: Science Edition, 2007, 45(4)： 519-523. (刘淑媛, 吕显瑞, 齐毅. 求Duffing方程周期解的Mountain Pass方法 [J]. 吉林大学学报: 理学版, 2007, 45(4)： 519-523.)
[4] LIU Wen-bin, LI Yong. Existence of Periodic Solutions of Higher Order Duffing Equation [J]. Acta Mathematica Sinica, 2003, 46(1): 49-56. (刘文斌, 李勇. 高阶Duffing方程周期解的存在性 [J]. 数学学报, 2003, 46(1): 49-56.)
[5] LI Yong, LIN Zheng-hua. A Constructive Proof of the Poincarae-Birkhoff Theorem [J]. Tran of the Amer Math Soc, 1995, 347(6): 2111-2126.
[6] 薛毅. 最优化原理与方法 [M]. 北京： 北京工业大学出版社, 2008.
[7] FENG Zi-xuan, XU Xu, JI Shu-guan. Finding the Periodic Solution of Differential Equation via Solving Optimization Problem [J]. Journal of Optimization Theory and Applications, 2009, 143(1): 75-86.
[8] YANG Liu, CHEN Yan-ping. On the Convergence of a New Levenberg-Marquardt Method [J]. Mathematica Numerica Sinica, 2005, 27(1)： 55-62. (杨柳, 陈艳萍. 一种新的Levenberg-Marquardt算法的收敛性 [J]. 计算数学, 2005, 27(1)： 55-62.)
