非线性微分方程求解

2024-07-22

非线性微分方程求解(精选十篇)

非线性微分方程求解 篇1

(y') n+a1 (x) (y') n-1y+……an-1 (x) y'yn-1+an (x) yn=0 (1)

的方程称为一阶n次非线性变系数微分方程, 其中a1 (x) , a2 (x) ……, an (x) 为连续函数。

定义2:若 (1) 式中a1 (x) , a2 (x) ……, an (x) 都为常数, 即

(y') n+a1 (y') n-1y+……+an-1y'yn-1+anyn (2)

为一阶n次非线性常系数微分方程。

下面我们就来讨论这两类方程的解。

定理1:若y=u (x) 是微分方程 (1) (或 (2) ) 的一个非零特解, 则y=Cu (x) 为方程 (1) (或 (2) ) 的通解。

二、一阶n次非线性常系数微分方程的特征方程法

我们将方程 (2) 的两边同除以yn得:undefined因为a1a2……an都是常数, 所以我们估计方程有形如:y=eλx (λ为常数) 的解, 将y=eλx带入到方程中得特征方程:

f (λ) =λn+a1λn-1+……+an=0 (3)

由上式得到特征根λ1, 则方程有非零特解y=eλ1x

将 (3) 的解代入到y=eλx中得到特解, 因此方程 (2) 的通解为y=Ceλx。

说明:若特征方程有多个解λ1λ2……λn, 那么方程的通解为y=ceλ1x, y=ceλ2x……y=ceλnx, (C为常数) , 且这些解是相互独立的。

三、一阶n次非线性变系数微分方程与常数变异法

分析 (1) 式和 (2) 式的关系, 如果 (1) 式中系数函数ai (x) (i=1…n) 都是常数即为 (2) 式, 也就是说 (2) 式是 (1) 式的一种特殊情况, 那么 (1) 式是否也同样含有形如 (2) 式的解呢?我们考虑将 (2) 的通解中常数C变成函数C (x) , 即y=C (x) eλx, 代入到 (1) 中得:

( (C (x) eλx) ') n+a1 (x) ( (C (x) eλx) ') n-1 (C (x) eλx) +……+an (x) (C (x) eλx) n=0展开后并整理得:

undefined

如果上式中括号里面的为常数, 方程 (1) 就为非线性常系数微分方程, 即undefined解得undefined。

故作变量代换undefined, 把方程 (4) 中undefined变成常数k1, 同样我们可以把 (4) 式中后面的中括号变成相应的常数系数。方程 (1) 就转化为方程 (2) 。

定理2: 一阶n次变系数非线性微分方程

(y') n+a1 (x) (y') n-1y+……an-1 (x) y'yn-1+an (x) yn=0

经变量y=C (x) eλx代换, 可转化为一阶n次常系数非线性微分方程

( (eλx) ') n+k1 ( (eλx) ') n-1eλx+…+kn (eλx) n=0

的充要条件是存在常数k1, 使得

undefined

其中ki (i=1…n) 为常数。

四、应用举例

例 求解方程 (y') 2+2tanx·y'y+sec2x·y2=0

解:通过观察确定方程为一阶2次变系数微分方程, 其中a1=2tanx, a2=sec2x, 由常数变异法可设y=C (x) eλx, 则y'= (C (x) eλx) 代入到方程中:

( (C (x) eλx) ') 2+2tanx· (C (x) eλx) ' (C (x) eλx) +sec2x· (C (x) eλx) 2=0

由定理2对上式展开后整理得到

undefined为常数

即undefined

取k1=0 k2=1

λ=±i所以y=ccosx (cosx+isinx) 或y=ccosx (cosxisinx)

五、结语

(一) 对于方程 (1) 的解我们也可以直接利用变量代换的方式令y=u (x) v (x) 带入到方程中化简, 得到相应的k1……kn。也可以求出方程的解。这里不作求解。

(二) 非线性方程在研究线性系统及其扰动系统的稳定性、生态学中等都有广泛的应用, 而非线性方程又没有统一的解法, 对于某些特殊的非线性方程研究他们的解法尤为重要。

摘要:用特征方程法求一阶次非线性常系数微分方程的解, 由一阶次非线性常系数微分方程的解去探求一阶次非线性变系数微分方程的解。

关键词:一阶非线性,常系数,变系数,微分方程

参考文献

[1].丁同仁, 李承治.常微分方程教程[M].北京:高等教育出版社, 1991

[2].戴中林.一类一阶高次微分方程的解法[J].大学数学, 2006

[3].龚东山, 牛富俊, 刘岳巍.一类一阶微分方程独立通解的研究[J].长沙大学学报, 2008

非线性微分方程求解 篇2

专业: 信息与计算科学

班级: 13405011 学号: 1340501123 姓名: 实验日期:报告日期:实验地点:邢耀光 数理学院五楼机房

2016.05.09

2015.05.13

超定方程组的求解

邢耀光

(班级:13405011 学号1340501123)

摘要:在实验数据处理和曲线拟合问题中,求解超定方程组非常普遍。比较常用的方法是最小二乘法。形象的说,就是在无法完全满足给定条件的情况下,求一个最接近的解。最小二乘法(又称最小平方法)是一种数学优化技术。它通过最小化误差的平方和寻找数据的最佳函数匹配。利用最小二乘法可以简便地求得未知的数据,并使得这些求得的数据与实际数据之间误差的平方和为最小。

关键字:最小二乘问题,残量,超定方程组,正则化方程组,Cholesky分解定理。

正文:

最小二乘法的背景:

最小二乘法(又称最小平方法)是一种数学优化技术。它通过最小化误差的平方和寻找数据的最佳函数匹配。利用最小二乘法可以简便地求得未知的数据,并使得这些求得的数据与实际数据之间误差的平方和为最小。最小二乘法还可用于曲线拟合。其他一些优化问题也可通过最小化能量或最大化熵用最小二乘法来表达。最小二乘法经常运用在交通运输学中。

交通发生预测的目的是建立分区产生的交通量与分区土地利用、社会经济特征等变量之间的定量关系,推算规划年各分区所产生的交通量。因为一次出行有两个端点,所以我们要分别分析一个区生成的交通和吸引的交通。

最小二乘问题:

最小二乘问题多产生于数据拟合问题。例如,假定给出m个点t1,...,tm和这m个点上的实验或观测数并假定给出在ti上取值的n个已知函数1(t),...,n(t)。考虑i 的线性组合,f(x;t)x11(t)x22(t)...xnn(t),(1)

我们希望在t1,...,tm点上f(x;t)能最佳的逼近y1,...,ym这些数据。为此,若定义残量 据y1,...,ymj1

则问题成为:估计参数x1,...,xn,使残量r1,...,rm尽可能地小。(2)式可用矩阵-向量形式表示为

ri(x)yixjj(ti),i1,...,m,(2)

n r(x)bAx,(3)其中

1(t1)n(t1)y1A, b,(t)(t)y1mnmm

TT)r(x(x,...x,x)(r(x),...,r(x)).1nm1

当mn时,我们可以要求r(x)0,则估计x的问题就可以用第一章中讨论的方法解决。当mn时,一般不可能使所有残量为零,但我们可要求残向量r(x)在某种范数意义下最小。最小二乘问题就是求x使残向量r(x)在2范数意义下最小。

定义1:给定矩阵ARmn及向量bRm,确定xRn,使得

bAx2r(x)2minr(y)2minAyb2.(4)

yRnyRn这就是所谓的最小二乘问题,简称为LS问题,其中的r(x)常常被称为残向量。

在所讨论的最小二乘问题中,若r线性依赖于x,则称其为线性最小二乘问题:若r非线性依赖于x,则称其为非线性最小二乘问题。

最小二乘问题的解x又可称做线性方程组

Axb,ARmn

(5)的最小二乘解,即x在残向量r(x)bAx的2范数最小的意义下满足方程组(5)。当mn时称(5)式为超定方程组。

定理1:(Cholesky分解定理)若ARnn对称正定,则存在一个对角元均为正数的下三角阵LRnn,使得

ALL.(6)(6)式称为Cholesky分解,其中的L称作A的Cholesky因子。

因此,若线性方程组Axb的系数矩阵是对称正定的,则我们自然可按如下的步骤求其解:

(1)计算A的Cholesky分解:ALL ;

(2)求解Lyb得y ;

(3)求解Lxy得x; 简单而实用的方法是直接比较ALL两边的对应元素来计算L。设

TTTTl11ll2122.Llllnnn1n2T比较ALL两边对应的元素,得关系式

aij 首先,由a11l11,得

l11再由ai1l11li1,得

li1ai1l11,i1,...,n.这样便得到了矩阵L的第一列元素。假定已经算出L的前k1列元素,由

akk得 2lp1jipjpl,1jin(7)

a11.lp1k2kp,122 lkkakklkp.(8)

p1k1再由

aikliplkpliklkk,ik1,...,n,p1k1k1 likaikliplkplkk,ik1,...,n.(9)

p1这样便求出了L的第k列元素。这种方法称为平方根法。

记最小二乘解的解集为LS,即

定理 ATAxATb.(10)

方程组(10)常常被称为最小二乘问题的正则化方程组或法方程组,它是一个含有n个变量和n个方程的线性方程组。在A的列向量线性无关的条件下,AA对称正定,故可用平方根法求解方程组(6),这样,我们就得到了求解最小二乘问题最古老的算法———正则化方法,其基本步骤如下:

(1)计算CAA, dAb;

(2)用平方根法计算C的Cholesky分解:CLL;(3)求解三角方程组Lyd和Lxy.TTTTT

2:xLS 当且仅当

LSxRn:x是LS问题(3)的解,实验 :

一:超定方程组的求解

原理:设A是mn阶矩阵mn,则线性方程组Axb为超定方程组,这里xR,bR。如

mm果A的秩为n,则称A为列满秩矩阵。超定方程组的解满足法方程AAxAb,该解使得

TTbAx 22min,称之为最小二乘解。

11 题目: 111TT1.11.12121.21.2221.31.3x3

21.41.441.51.525

用正则化方法求解,要求:(1)BLL 不得使用MathCAD指令Cholesky;(2)BLL使用MathCAD指令Cholesky。

11解:(1)A111 1.11.126.58.5551.21.22T8.5511.375 1.31.32 , 则 BAA6.51.41.428.5511.37515.2981.51.5215B21,20.5L2.907 , gATbLB2.236, 211111L1128.25

L31 B3123.824 , L22B22L210.316 , L11 4 LB32L31L2132L0.822 , LL2233B3331L320.037 , 22

2.236002.2362.9073.824 即 L2.9070.31600.037 , LT00.3160.822, 3.8240.822000.037

6.70810 则yL1g3.162xLT1y10 ,  , 9.27310132.4781011

x即为所求的最小二乘解。

11.11.1211.21.22(2)A2.2360011.31.32cholesky(B)2.9070.316011.41.423.8240.8220.03711.51.52,2.236002.2362.9073.824 则 L2.9070.3160,LT00.3160.8223.8240.8220.037000.037,6.70810 则yL1g3.162xLT1y10 , 

9.27310132.4781011,x即为所求的最小二乘解。

二:已知如下数据: xi0.00.20.40.60.81.01.2yi0.91.92.83.34.05.76.5 利用最小二乘法拟合曲线 ya1xa2.0.00.90.21.9解:令B0.00.20.40.60.81.01.20.42.80.91.92.83.34.05.76.5 ,x0.6,y3.3 0.81.04.05.71.26.510.010.210.41 则A10.6XATAATy0.8434.571,即p(x)0.8434.571x, , 10.811.011.2故最小二乘法拟合曲线为y4.571x0.843.程序附录: 一;

11.11.12111.21.22256.58.55A11.31.32b3AxbBATA gATb, B6.58.5511.37511.41.428.5511.37515.298,45, ,11.51.52 f(B)nrows(B)

Lidentityn()fork1nLkkBkkifk1k1LkkBkkLkp2otherwisep1forik1nBLkikiLifk1kk(break)ifknk1BikLipLkpLikp1LotherwisekkL

15g20.5, , 28.25,002.23602.23602.2362.9073.824Tf(B)2.9070.3160L2.9070.3160L00.3160.82200.0373.8240.8220.037,Lf(B),3.8240.8220.037,0,6.708103.16210yx119.27310132.4781011TyLg , xLy ,, ,二;

120.00.20.40.60.81.01.2TT B , xB , yB, 0.91.92.83.34.05.76.5

00.90.21.90.42.8x0.6 ,y3.3 , x0, x0.2,nrows(x),n120.8415.71.26.507,i1n , Ai11,.111Ax,AXy,A1i2i111p(s)0.8434.571x.0.20.410.843,p(s)XT1

TTAy ,X0.6,XAA4.571s0.811.2

心得体会:

通过本次的课程设计,让我学会了很多,学会了简单的MathCAD 软件的用法。让我更加深刻了解最小二乘问题,和以往对知识的疏忽得以补充。不仅掌握了学习的知识,而其还可以培养和熟练使用资料,运用工具书的能力,把我们所学的课本知识与实践结合起来,起到温故而知新的作用。

参考文献:

1,《数值线性代数》(第二版)北京大学出版社 徐树方,高立,张平方,编著。2013.01

email:974671870qq.com

1340501123:邢耀光

非线性微分方程求解 篇3

[关键词]:二阶变系数齐次线性微分方程 二阶变系数非齐次线性微分方程 一阶线性微分方程 通解 特解

一、引言

在微分方程的理论中,二阶线性微分方程占有十分重要的位置.国内外现行 《高等数学》中的方程[1],只是对常系数微分方程的情况做了详细的讨论,《常微分方程》也未对二阶变系数微分方程的解作进一步的阐述.一般的变系数微分方程的通解没有普遍的求法.在高等数学中,只有对欧拉方程这种特殊的变系数微分方程研究了它的解法[2].本文通过研究几种特殊的二阶变系数线性微分方程通解的求法,推导二阶变系数线性微分方程的通解的一般求法,包括二阶变系数齐次线性微分方程和二阶变系数非齐次线性微分方程.诣在解决二阶变系数线性微分方程的求解问题,并得出成规的求解的方法与结论,以便适应在工程技术的实际领域或学生在学习相关专业中的需要.

二、一些特殊型二阶变系数线性微分方程的求解方法

考虑如下的二阶变系数微分方程

(1)

其中 都是连续函数,当 满足一定条件时,通过适当的变量代换,方程(1)可化为欧拉方程,进而求出其通解.

引理1[3] 假如方程

中 .

只需引进变量 ,则方程可化为

.

定理1当 时,方程(1)可通过变量代换 化成欧拉方程且通解为

(Ⅰ) 时, .

(Ⅱ) 时, .

(Ⅲ) 时, .

证明:设 ,这里 为待定的连续可微函数,此时有

将 代入方程(1),得

(2)因为

所以(2)式化为

. (3)

若設 ,则(3)为欧拉方程,将其代入(3)得

. (4)

. (5)

令 ,则 ,(5)式化为

(6)

(6)式变成为常系数线性微分方程未知函数为 ,自变量换成 .求解(6),再由 代回 ,得到 ,从而得到方程(1)的解.具体解法如下:

方程(6)对应的特征方程为

判别式 .分以下三种情况讨论:

(Ⅰ)当 即 时,特征方程有两个相异实根

, .

方程(6)的通解为

.

将 代入上式,得

.

所以方程(1)的通解为

.

若 ,即 时,方程(1)的通解为

.

(Ⅱ)当 即 时,特征方程有两个二重根

.

方程(6)的通解为

.

将 代入上式有

.

所以方程(1)的通解为

.

(Ⅲ)当 即 时,特征方程有一对共轭复根

, .

方程(6)的通解为

.

将 代入上式,有

.

所以方程(1)得通解

.

例1 求方程 的通解.

解:这里 , ,且 .此处 .符合情况(Ⅰ).将 , 代入通解公式并化简得方程的通解为

.

例2 求方程 的通解.

解:这里 , ,此处 .符合情况(Ⅲ),将 , 代入通解公式并化简得方程的通解为

.

可见,在方程(1)的系数满足条件 时,便可由定理的公式直接求出这类方程的通解,避免了繁琐的变量代换求通解的过程.

设 阶变系数线性非齐次微分方程

(1)

对应的齐次方程为

(2)

定理2 若方程(2)的特解为 则方程(1)的通解为

引理2[4] 对于 阶变系数线性微分方程

有以下结论:

(1)若 ,则特解为 ;

(2)若 ,则特解为 ;

(3)若 ,则特解为 ;

(4)若 ,则特解为 .

上述寻找特解的方法要求系数要满足一定的条件,有时并不好实现.但对于一些二阶变系数线性微分方程可通过变量代换化为常系数方程,从而很容易求解.

考虑二阶变系数齐次线性微分方程

(3)

其中 都是连续函数,当 , 满足一定条件时,通过适当的变量代换,方程(3)可化为常系数微分方程,进而求出其通解.

定理4 若存在常数 ,使得方程(3)的系数 满足:

其中 为常数,则方程(3)可以化为二阶常系数线性齐次方程

, (4)

进而可求得原方程的通解.

证明:令 ,则 , ,代入方程(3),得

. (5)

又因为 ,所以 且 ,代入

(5),得 ,再将 代入此二阶常系数线性齐次方程的通解,便得原方程的通解.

nlc202309011050

注:一般情况下,可令 .

例3 求方程 的通解.

解:因为 ,则方程有特解 .于是,方程有形如

的通解.将 代入方程得 .

. (4)

解方程(4)得

.

(5)

解方程(5)得

.

于是,原方程的通解为

.

例4 求方程 的通解.

解:因为 , .显然,存在常 数,使得 .令 ,则 , ,代入方程得

(6)

方程(6)为二阶常系数线性齐次微分方程,其通解为

.

故原方程的通解为

.

例5 求方程 的通解.

解: , .显然,存在常数 ,使得 .令 ,则

, ,

代入方程得

. (7)

方程(7)为二阶常系数线性非齐次微分方程,其通解为

.

故原方程的通解为 .

定理5 對于二阶变系数线性非齐次微分方程

, (1)

若 ,和实数 满足

, (2)

则(1)的通解为

(3)

其中积分 , 和 都表示一个原函数, 和 为任意常数.

证明:设方程(1)的解为 ,求导得

, ,

将 代入(1)化简得

(4)

在(4)中,不妨令

(5)

方程(5)为二阶变系数线性齐次微分方程,文献[4]给出一类若 满足(2),则方程(5)的一个特解必为

(6)

将(5)代入(4)整理得

(7)

显然(7)为可降阶的微分方程.利用可降阶的微分方程的求解方法可求得(7)的通解(即求得 )为

.

其中积分 , 和 都表示一个原函数, 和 为任意常数.由此得(1)的通解为

例6 求方程 的通解.

解:文献[4]中给出 =-3,所以对应齐次方程 的一个特解为 .又 , ,代入(7)得方程的通解为

= =

=

例7 求方程 的通解.

解:文献[4]中给出 =-1,所以对应齐次方程 的一个特解为 .又 , 代入(7)得方程的通解为

= =

注意:定理中的条件 非常苛刻,只有在特殊条件下才能满足.但通过分析可以发现,定理中的这个条件并不一定必须满足,只要能知道方程(1)对应齐次方程的一个非零特解 ,则由公式(7)同样可以求出方程(1)的通解.

例8 求方程 的通解.

解:已知对应齐次方程 的一个特解为

.

又 , ,代入(7)得方程的通解为

=

=

= 1

三、一般型二阶变系数线性微分方程的求解方法

定义[5] 若 、 为连续非常数的函数,则称方程

(1.1)

为二阶变系数线性微分方程. 如果 恒等于零,那么该方程称为二阶变系数齐次线性微分方程;如果 不恒等于零,那么该方程称为二阶变系数非齐次线性微分方程.

假设二阶变系数非齐次线性微分方程中 具有一阶连续的导数、 连

续.令

, (1.2)

, (1.3)

则方程(1.1)变形为

(1.4)

令 (1.5)

那么原方程(1.1)就化简为

.

解之,得 ,将之代入(1.5)式,则方程(1.5)通过上述变换可降阶为

此一阶线性非齐次微分方程的解就是我们所要求的二阶变系数非齐次线性微分方程的解[6],而方程

的解为

(1.6)

故式(1.6)为二阶变系数非齐次线性微分方程

的通解公式.

另外由(1.2)式又得

将其代入式(1.6)可得二阶变系数非齐次线性微分方程(1.1)通解的另两种形式为

(1.7)

或 (1.8)

特别地1 当 时,方程(1.1)就转化为二阶变系数齐次线性微分方程,而式(1.6)、(1.7)、(1.8)分别为

, (1.9)

, (1.10)

, (1.11)

它们是对应的二阶变系数齐次线性微分方程 的通解公式. 以上的求解过程或方式就是二阶变系数线性微分方程的求解方法,(1.6)、

(1.7)、(1.8)均为二阶变系数非齐次线性微分方程

的通解公式. 公式(1.9)、(1.10)、(1.11)均为二阶变系数齐次线性微分方程

的通解公式,.在具体应用时,应依据问题灵活使用.

特别地2 形如

型的方程可化为伯努利方程.

原方程变形为

则原方程就化为伯努利方程

即可求得其解[7] .

四、举例

nlc202309011050

运用二阶变系数线性微分方程的一般求解法求二阶变系数线性微分方程的解时,其重点是构造(1.2)和(1.3)式,难点或关键点是从(1.2)式和(1.3)式,求出 和 . 或由(1.2)式和(1.3)式变形得

(1.12)

和 , (1.13)

再从中求得 和 ,然后用上述方法或上述公式,可求得二阶变系数线性微分方程的解.

注:方程(1.12)或(1.13)是黎卡提(Riccati)方程,见《常微分方程》教[8].求二阶变系数线性微分方程解时,必须观察二阶变系数线性微分方程的特

征.如果是上述特殊类型的二阶变系数线性微分方程,就用特殊类型的二阶變系数线性微分方程的求解方法求之;如果不是上述特殊类型的二阶变系数线性微分方程,就用二阶变系数线性微分方程的一般求解方法求之.

二阶变系数线性微分方程的一般求解步骤:

第一步:构造(1.2)式和(1.3)式;

第二步:计算出 , ;

第三步:将第二步的结果代入上述公式求出通解来.

例9 求 的通解

解:由方程特征可知,

, ,

则 的通解为

.

例10 求 的通解.

解:令 ,解之得

由以上公式,所求方程的通解为

.

五、总结

对一般的二阶变系数线性微分方程而言,由《常微分方程》 教材[5] 知,只要能求出二阶变系数齐次线性微分方程的一个特解,则二阶变系数线性齐次或非齐次微分方程的解即可求得.尽管专家学者目前的研究给出了一些特殊类型的二阶变系数线性微分方程的求解方法,然而,如何求出其中的某一特解是无法可循的.通过研究特殊类型的二阶变系数线性微分方程的求解方法,深入研究了一般二阶变系数线性微分方程的求解方法,包括二阶变系数齐次线性微分方程和二阶变系数非齐次线性微分方程.

参考文献:

[1]同济大学应用数学系.高等数学[M] .北京:高等教育出版社,2003:259-317.

[2]丁同仁,李承治.常微分方程教程[M].北京:高等教育出版社,1991:183-185.

[3]陈惠汝,刘红超.二阶变系数线性微分方程的变量代换解法[J].高等函授学报自然科学版,2008,21(3):21-22.

[4]王高雄,周之铭.常微分方程[M].第2版.北京:高等教育出版社,1983.

[5]王高雄.常微分方程[M].北京:高等教育出版社,1985:101-163.

[6]龚东山,刘岳巍,贾筱景.计算一类常微分方程特解的新方法[J].河北北方学院学报:自然科学版,2008,24(06);1-3.

[7]常秀芳,李高.伯努利方程的几种新解法[J].雁北师范学院学报,2007,23:(02).89-91.

[8]王高雄.常微分方程[M].北京:高等教育出版社,1985:62.

求解微分方程的线性多步法 篇4

其中f为t和u的已知函数,u0为给定的初值。我们假设函数f (t, u)在区域:t0≤t≤T,|u|<∞内连续,并且u满足Lipschitz条件,即存在常数L,对所有t∈[t, T]和u1, u2,有|f (t, u1)-f (t, u2)|≤L u1-u2|。在上述基本假设下,初值问题(1)在区间[t, T]上有唯一解,并且u (t)为连续可微的。

将区间[0, T]作N等分,小区间的长度h=T/N称为步长,点列tn=nh (n=0, 1,…,N)称为节点,t0=0。由已知初值u (t0)=u0,可算出u (t)在t=t0的导数值u′(t0)=f (t0, u (t0))=f (t0, u0)。利用Taylor展式得u1=u0+hf (t0, u0), u1就是u (t1)的近似值。利用u1又可算出u2,如此下去可算出u在所有结点上的值,一般递推公式为un+1=un+hf (tn, un) , n=0, 1, …, N-1。由文献知求解该类方程

由文献[1]知,求解该类方程分为单步法和多步法。只要利用h, tm和um即可算出um+1的算法称为单步方法,单步法包括Euler法、梯形法、Runge-kutta法,其一般形式为um+1=um+hφ(tm, um h),函数φ称为增量函数,它依赖于所给的微分方程。用单步法求出{φm},m=1, 2,…,N只需要一个初值u0。

需要用到h, tm, tm+1,…,tm+k-1和h, um, um+1,…,um+k-1才能求出um+k (k>1)的算法称为多步方法,线性多步法包括数值积分法和待定系数法,其一般形式为,这里αj,βj (j=0, 1, K, k)是常数,它需要有k个初值u0, u1,…,uk-1才能得到整个序列{um},m=1, 2,…,N。

2. 线性多步法的改进

常微分方程初值问题只给我们提供了一个初值,但线性k步方法需要k个初始值u0, u1,…,uk-1,所以这其余的k-1个初始值就要通过其他方法来先期得到。最容易想到的是采用简单的单步方法,如Euler方法,或改进的Euler法,但它们的收敛阶较低。为保证整个计算的精度,我们应使初始值的计算精度至少与所用的线性多步方法同阶。

因此,若使用的线性多步方法为q阶[2],[3]的,可以用q阶Taylor展开方法或Runge-kutta方法来求u0, u1,…,uk-1,然后进入正式的求解过程。用Taylor展开确定初始值的同时还能得到关于选择h的信息。如要求计算误差不超过ε,则应使

同时成立。当微商次数增大时,计算量也将急速增大,因此利用低阶微商构造高精度的方法对实际计算是很有意义的,对u (tn+1)=u (tn)+蘩ntn+1u′(t) dt利用u′(t)在点tn+1, tn处的函数值及其各阶导数值构造Hernite型插值多项式,再代积分,舍去误差项即可得到带导数值的计算公式。

对隐式的计算格式一般采用迭代解法,将隐式线性多步方法改写成

由于Fm是已知的,于是(2)式可用迭代法求解。

无疑是十分重要的,很自然的想法是用显格式算出un+k (0)———这称为预测,再用隐格式进行迭代求un+k———这称为校正,因此整个过程称为预测—校正算法[4],简称PC算法。显然,预测的P式与校正的C式的阶一般应取成相同的。

预测—校正算法求un+k的过程一般如下:

由于使用了同阶的预测公式求得初始近似un+k (0),因此一般说来校正次数不会太多,大约两三次就够了。若校正步数太多则往往提示步长过大,应适当减小步长后再进行计算。实际计算中经常利用算式与算式同阶的特点,对(3)式再作一些修正。用线性组合的方法消去每个格式中的截断误差主项。仅花费很小的代价就使方法的精度再提高一阶,从而不需要再通过迭代来改善计算的近似值。

设预测算式和校正算式的局部截断误差分别为

整理后分别得到:

这两个式子告诉我们,在求(un+k (c)和un+k (p))后,再利用它们的局部截断误差[5],[6]的主项系数的某种组合对它们进行修正,可使得修正后的局部截断误差比原来高一阶,即达到o (hq+2)。这里被称为修正系数。

直接按照(5)式对un+k (p)进行修正是不行的,因为此时的un+k (c)还没有算出。但是,在计算un+k的近似值时,un+k-1已经算出来了,因此u (c) n+k-1和u (p) n+k-1也都已经有了,这时可以用u (c) n+k-1-u (p) n+k-1去代替(5)中的un+k (c)-un+k (p)。

综上所述,我们可以将带修正的预测—校正算法的一般形式表示如下:

在多步方法中,改善步长将带来一个新问题,就是结点变了。在预期结点处已经算出的近似值,对于计算当前结点说,很多都用不上了,而在以新步长为单位的结点处,往往没有计算过相应的近似函数值。在使用的方法中,一般采用值的方法去补上所需的这些点上的近似值再进行计算。

下面给出几种常用的带修正的预测—校正算法。

(1) Milne预校算法。Hamming针对Milne算法绝对不稳提出了如下修正:

P和C均为四阶方法,M1和M2中的修正系数分别为

(2) Hamming预校算法。该预校算法的P式和C式分别为:P∶与Milne预校算法相同。

P和C也均为四阶方法,M1和M2中的修正系数分别为

(3) Adams四阶预校算法。此算法中的P算式和C算式分别为:P∶Adams四步四阶外插公式。P∶Adams三步四阶内插公式。

M1和M2中的修正系数分别为270251和-19270。理论研究表明,这个算法的稳定性比前两个算法好,是常用的预校算法。

3. 结语

针对微分方程初值问题的线性多步法,我们利用局部截断误差的主项系数的某种组合对其进行修正,使修正后的局部截断误差比原来高一阶,即达到o (hq+2),并给出了带修正的预测—校正算法的一般形式,见文中(7)式。

摘要:微分方程初值问题的求解一直是人们所关注的热点问题, 本文针对微分方程初值问题的线性多步法, 利用它们的局部截断误差的主项系数的某种组合对它们进行修正, 提出了预测—校正算法, 并给出了几种常用的预测—校正算法。

关键词:微分方程,初值问题,线性多步法

参考文献

[1]李立康, 於崇华, 朱政华编著.微分方程数值解法[M].复旦大学出版社, 2005.

[2]胡健伟.微分方程数值方法[M].科学出版社, 1999.

[3]李荣华.微分方程数值解法 (第三版) [M].高等教育出版社, 1996.

[4]Dennis G.Zill编著.微分方程与边界值问题[M].机械工业出版社, 2003.

[5]黄振侃.编著.数值计算——微分方程数值解[M].北京工业大学出版社, 2006.

非线性微分方程求解 篇5

中立型延迟微分方程一般线性方法的非线性稳定性

本文研究一类非线性中立型延迟微分方程一般线性方法的数值稳定性.证明了一般线性方法为(k,p,O)-代数稳定时,在一定的`约束条件下,其数值解保持微分方程理论解的稳定性质,特别是证明了在约束网格情形代数靛的-般线性方法能无条件保持解析解的稳定性.

作 者:李超群 LI CHAOQUN 作者单位:中国地质大学数学系,武汉,430074刊 名:应用数学学报 ISTIC PKU英文刊名:ACTA MATHEMATICAE APPLICATAE SINICA年,卷(期):31(2)分类号:O241.8关键词:中立型延迟微分方程 一般线性方法 (k,P,O)-代数稳定 代数稳定

洋物流求解中国方程 篇6

独资?合资?

2004年,外国物流企业所面临的一个核心问题是“应该采用什么样的运营模式”:是建立全资子公司,还是合资公司?按照中国加入世贸组织的协定:到2005年年底,中国应修改相关法规,取消对外商投资物流企业所有权的限制。因此,届时是否建立合资公司将成为一个商业决定,而不必考虑法规限制。

今天,我想和大家分享一下选择物流合资公司中国合作伙伴的思考。了解外国物流企业的思考方式并确定评估外国合作伙伴的流程,或许将帮助您们在与外国企业建立合作关系时提供一些指导。

需要强调的是:一些外国物流企业已经决定采取全资子公司的形式。瑞士的K&N(Kuehne & Nagel)公司通过两种途径确立了市场地位:首先是在中国18个沿海城市建立起了代表处网络;其次是利用与胜科物流(SembCorp Logistics)的交叉持股关系,通过新科安达(ST Anda)加强在中国市场的渗透。通过内地与香港的CEPA合作,K&N已经获得了在2005年年底全面开放之前建立全资子公司的许可,并在2004年4月获得了交通部的甲级许可。以在上海宝山和外高桥的两个物流中心为基础,K&N还计划在内地如武汉、成都以及长春开设新的办事处和配送中心。这表明K&N将大量投资,从而实现迅速扩大市场覆盖并获得竞争优势的目标。

然而,并不是所有的外国物流企业都有足够的经验或承担风险的意愿来设立全资子公司。一些企业更倾向于寻找中国合作伙伴建立合资公司的方式,更快速、更简洁地进入市场。但是寻找、评估合作伙伴以及与之进行谈判,在世界上任何地方都是很复杂的过程。在中国,由于市场瞬息万变以及中外企业在寻求合作伙伴问题上的文化差异,这个过程尤其具有挑战性。

合资图什么?

外国物流企业在中国建立合资公司,通常会希望中国合作伙伴能够带来以下方面的优势:

整体上更为快速地进入市场

与地方政府和中央政府的关系

在中国的网络能力

在中国的客户关系

人才

同样,外国合作伙伴通常会希望在以下几方面做出贡献, 以期和中国企业相辅相成:

整合的物流解决方案

对某一具体行业领域的深入了解

国际网络能力

跨国公司客户关系

管理技巧

谁将携手洋物流?

美世常常应一些外国企业的要求,帮助发现中方合作伙伴并对之进行评估。在对大量中国企业进行筛选后,我们会对少数入围企业进行详细评估。我们在评估入围企业是否适合与外方建立合资公司的过程中,主要关注以下五个方面。

一、合作方资质

我们通常采用四个标准:

组织效能: 不同的中国物流企业,其组织结构差别迥异。一些企业采用分散式结构,各地区子公司占主导地位;而其它一些企业的总部则对大客户关系、解决方案设计及服务定价等实施更为集中化的控制,这些企业能够对网络有明确的控制,并更便于提供整合解决方案。

管理团队素质:在一些中国企业,管理权力高度集中,控制在一个领导人手中;而在其它一些企业中则存在着范围更广、更深的管理模式。我们对中国企业内部所有重要经理人的教育背景也是非常关注的。

市场声誉:我们通过货主和行业专家了解情况,将候选合作企业与其它中国物流企业从服务质量以及行业道德标准的角度进行比较。

财务优势: 这一点不但从对合资企业初期投资角度来看非常重要,而且从长期合作过程中对合资公司持续投资以保证业务增长来看也非常重要。如果合作企业有母公司,那么我们还需要密切关注其母公司的优势及其对子公司的承诺。

二、网络能力

我们将该条标准分为以下几个方面:

产品和服务的广度和深度/网络覆盖: 对于外国物流企业来说,能够获得合作方在国内市场的网络至关重要。进入中国市场的外国物流企业很难自行建立起国内网络,尤其是在中国的内地。因此,我们会花费大量时间,参观合作方在国内网络的设备,确认其网络覆盖程度。

设备/管理质量: 通过参观国内网络中具有代表性的设备,如仓储设备和地区配送中心(RDC),我们将对服务质量以及设备的管理和运营状况进行评估。

供应商管理质量:我们将了解中国第三方物流行业在管理运输和仓储服务提供商方面的经验,这对于合资公司的未来发展是至关重要的。

三、关系网

外国企业希望中方合作者能够带来政府关系和客户关系。

政府关系: 外国企业非常重视中方合作者与中央政府及地方政府的关系,不仅仅因为这能够促进政府对合资项目的初期审批,而且因为政府关系能够确保合资公司在未来能够继续获得必要的运营许可,尤其是在运输方面。

客户网络的广度和深度:尽管外国企业可能与跨国企业货主之间已经有良好的关系,它仍然也很注重中方合作者所能带来的本地客户关系。这也能够帮助合资项目在初期获得更多收入。然而,由于定价和关系等考虑,把中国企业的客户关系转移到合资公司具有一定的挑战性。

在某一行业领域内的客户群覆盖:中国面临着向行业专业化发展的趋势。中方合作者在高科技、汽车或消费产品等具体行业领域内的优势具有重大价值。

四、战略吻合度

如果双方之间缺乏战略一致性,那么将很难成立合资公司。这条标准可以总结为以下三方面:

愿景和承诺:中方在成立合资公司事宜上与外方是否愿景相同,双方是否真正地承诺全力合作。

第三方物流战略一致性:中方合作伙伴在客户选择、价值定位以及增长途径等方面,是否与外国企业有共同的战略目标?

对于合作关系的观点: 中方是否对合作关系有长期展望?合作关系是否平等?这是一个双赢的合作关系,还是中方企业想利用这个机会向外方企业学习经验然后独立运营?

五、达成一致的可能性

如果中方企业不能高度符合这条标准,那么其它方面的标准都无异于空谈。这个问题就是:“交易是否能够达成?”

股份结构的公开:这是一个较难的问题,需要在早期明确。强大的中国企业通常很难放弃多数股权。

运营控制的公开:在运营层面上,双方应该在合作早期明确某些重要职位如总经理、财务总监等应由哪方人员担任?这些职位的确定可以体现出持股情况,但也可以有例外。

业务范围的灵活性以及对独家性的承诺: 有时,一些实力强大的中方企业可能会被要求将某一领域的业务放到合资公司中并承诺不会在合资公司外继续发展该领域业务。这意味着,为了获得长期利益,不得不以一定的短期利益作为牺牲。

以夷人之道度夷!

综上所述,从2005年开始,一些实力强大的外国物流企业将会以全资子公司的形式进入中国市场;尽管政策法规会有所放松,但一些企业还会寻求中方合作者,来帮助它们迅速、高效地进入市场。

我们看到,这些希望成立合资公司的外国企业已经意识到了做好充分准备工作的重要性,以及对中方候选合作伙伴进行全面评估的重要意义。当中国企业需要对外国合作伙伴进行评估并建立合资企业的时候,可以使用相同的分析技巧,建立更为成功的合资企业,使双方获得收益。管理

遗传算法求解非线性方程组研究综述 篇7

为解决以上问题,众多研究者将目光转向了遗传算法。遗传算法研究兴起于20 世纪80 年代末,是模仿自然选择和遗传机制的一种智能优化算法,隐含并行性和群体全局搜索是其两个显著特征,具有较强的鲁棒性,对于一些大型、复杂非线性系统求解具有独特的优越性能[3 -4]。但遗传算法也存在着收敛速度慢, 精度差等缺点。遗传算法一般用于求解优化问题,但近年来人们也将其用于求解非线性方程组。

1 遗传算法求解非线性方程组基本思想

非线性方程组是指有n个变量m个方程所组成的方程组

式中至少含有一个非线性方程。其求解是指在其定义域D内找出一组数x*= ( x1*,x2*,…,xn*) 能满足方程组中的每一个方程。遗传算法的思路是将求解方程组( 1) 转化为如下优化问题

则求解方程组就转化为求一组值x*= ( x1*,x2*,…, xn*) 使得F( x*) = 0 成立,即求使函数F( x1,x2,…, xn) 取得最小值0 的一组数。其中F( x1,x2,…,xn) 也被称为适应度函数。

遗传算法将定义域D分为一个或多个独立的搜索空间。每个搜索空间构成了各个种群的范围。每个种群由经过编码的一定数目的个体组成,这里的编码常用二进制编码,也可采用格雷码或实数编码。编码的长度,通常由算法所要求的精度所决定。若采用二进制编码,并设目标解的第i个分量xi∈[ai,bi],其精度要求为 εi,则此分量编码的长度Li必须满足

当算法开始后,根据适应度函数计算种群中每一个体的适应度,根据一定的优化准则,开始产生新一代的种群。为产生下一代,按照适应度选择个体,父代要求基因重组而产生子代。这个过程通常由选择、交叉和变异3 步构成。

选择算子: 通常采用赌轮法或锦标赛法。在赌轮法中,各个体的选择概率与其适应度成正比。

交叉算子: 通常在编码中随机产生一个交叉点,对于以交叉概率Pc选中的配对个体,将交叉点右侧的编码进行简单交叉。例如: 若采用二进制,编码长度为8,则交叉点有1 ~ 7 共7 个位置,若设交叉点为5,选中的配对个体A为11001000,个体B为00110111,则交叉后的新个体A'为11001111,B'为00110000。在某些改进的算法中,也有人采用了多点交叉。

变异算子: 通常采用基因位等概率变异,对于种群中的每一个体,以变异概率Pm改变某一个或某一些基因座上的基因值为其他的等位基因。

经过以上3 步所产生子代的适应度又被重新计算,子代被插入到种群中将父代取而代之,构成新的一代。这一过程周而复始。在遗传算法求解非线性方程组的过程中,通常会设置种群的进化代数,当达到预定的代数或预设的精度要求时,一般认为已求得非线性方程组的解。

Goldberg[5]将采用二进制编码,选择算子使用赌轮法,交叉算子使用单点交叉,变异算子使用基因位等概率变异的遗传算法称为基本遗传算法。

在上述求解非线性方程组的过程中,遗传算法的并行性有效提高了求解效率,择优机制使得求解过程具有良好的全局优化性和稳健性,目标函数解释为编码化个体的适应度使求解过程具有良好的可操作性和简单性。但在求解过程中,遗传算法收敛速度慢、局部搜索能力较差、容易出现早熟等问题均影响了求解效率。

2 遗传算法求解非线性方程组的改进

2. 1 搜索空间方面的改进

搜索空间是遗传算法求解非线性方程组中由定义域D决定的一个或多个范围。搜索空间的大小,在一定程度上决定了求解的效率。

曾毅[6]提出了一种将种群中的个体按照适应值从大到小排列,根据前N位个体变量的对应二进制编码最高位与最优解对应的二进制编码的最高位是否相同的情况,对搜索空间进行缩小和移动的方法。该方法使得每个变量对应的二进制编码长度无需过长,便可求出高精度的解。

成媛媛等[8]提出了根据历代最优解将搜索区间不断划分成较小的区间,不断在各个小的搜索区间递归地调用动态自适应遗传算法,直到达到精度要求或小区间宽度为零停止的方法。这一方法较好地克服了传统方法的局限性,能够在较大的初始区间上以概率为1 搜索到区间内的全部解,并能较好地实现并行计算, 大幅提高了一般遗传算法求解非线性方程组的收敛速度和稳定性。

排新颖等[9]提出了一种梯度的变异步长,即在迭代开始时使用大的步长来进行全局寻优,在后期基本确定真正解的区间时,采用较小的步长进行局部细化搜索。

从上述文献来看,在搜索空间方面的改进,主要以动态缩小搜索空间等为主。缩小了搜索空间,意味这提高了求解的速度。因此,这些改进将是对基本遗传算法求解非线性方程组的一种有效补充。

2. 2 编码方面的改进

综上所述,在基本遗传算法求解非线性方程组的过程中,采用了二进制编码。其在程序上易于实现,但实践发现,该编码也带来了编码和解码误差等问题。

因此,曾毅[7],彭灵翔[10]均提出了在求解过程的编码环节使用实数编码来替代二进制编码,从而省略了复杂的编码和解码过程。实数编码虽使得设计和分析难度增加,但通用性较好,可方便地表示范围较大的数,便于在较大的空间进行搜索,同时也改善了遗传算法的计算复杂性,提高了运算效率。

传统二进制编码存在Hamming悬崖问题。对于二进制编码方法表示的个体,变异操作有时虽只是一个基因座的差异,而对应的参数值却相差较大。但若使用格雷码来对个体进行编码,则编码串之间的一位差异,对应的参数值也只有微小的差别。这就相当于增强了遗传算法的局部搜索能力,便于对连续函数进行局部空间搜索。因此Angel Fernando Kuri - Morales[11]采用格雷码,取得了较好的使用效果。

2. 3 传算子方面的改进

在使用遗传算法求解非线性方程组时,一般会使用到选择、交叉和变异3 个算子。选择算子一般通过赌轮法,交叉和变异操作一般采用根据经验所得的固定概率。

罗亚中等[12]在算子设计时以联赛竞争算子为选择算子,该算子可直接通过比较个体的目标函数值完成操作,因此该文的遗传算法适应度函数直接选择为优化目标函数。文中采用了最优保存策略,即当前群体中最优个体不参与交叉运算和变异运算,而是用其来替代本代群体中经过交叉、变异操作后所产生的最差个体。

胡小兵等[13]通过设计特定的交叉概率和变异概率函数,实现了交叉概率Pc和变异概率Pm的动态改变,从而增加求解过程中算法的进化能力和收敛速度。

田巧玉等[14]采用了启发式交叉的方式,使用参数 “Ratio”指定子辈离较好适应度的父辈有多远。如: 父辈是parent1 和parent2,而parent1 有较好的适应度,则启发式交叉生成的子辈如下

该文中还采用了Gaussian变异算子。此变异算子在进行变异操作时,将一正态分布、具有均值的随机数加入到父向量的每一项,替换原有的基因值。该Gaussian变异算子由两个参数“Scale ”和“Shrink ” 决定。

针对采用非线性排序,排序的结果通常易陷入局部解。徐红[15]提出了设置参数r,将适应值差别因素以r的比例加入到选择中,能改善搜索性能。参数r的设计及由此得到的交叉概率p'如下

交叉前排序的作用使相邻个体的适应度差别最小,特性相对集中。即采用先排序后交叉可能带来两种效果: ( 1) 加快收敛速度。( 2) 形成未成熟收敛。希望出现的是前一种情况,因而在搜索初期采用交叉前排序的操作以期提高搜索效率。

彭灵翔等[10]在变异算子上使用了最佳个体保留算子和灭绝与移民算子。采用最佳个体保留算子是一种常见的做法,可将有最好适应度的染色体群作为种子传到下一代。灭绝与移民算子在当交配池中的染色体几乎相同时,或最佳个体的适应度在一定代数内的增幅小于门槛值时起作用。灭绝与移民过程一般分为两部分,开始将最佳染色体之外的染色群全部灭绝。 然后移民产生Np- 1 个新的染色体。移民过程是: 从第2 个到Np/2 个染色体用产生初始种群的方法产生父辈染色体。其余染色体用最佳适应度的染色体通过变异产生。

综上所述,对于算子的改进,主要采用动态变化来替代固定概率的方式。无论哪种改进,均以加快收敛速度和防止早熟为主要目的。从各文献提供的实验数据来看,均取得了一定的效果。

3 混合遗传算法求解非线性方程组

混合遗传算法求解非线性方程组主要是通过遗传算法与其他算法相结合等手段来求解的。

3. 1 与经典迭代方法的混合

遗传算法中的杂交算子和变异算子能在全变量空间内搜索方程组的解,经典迭代方法能在收敛点附近局部细致地搜索解。若能将两者混合,有可能发挥好两类算法各自的优势而取得更好的求解效果。目前文献中主要提出了两类混合策略。

( 1) 将经典迭代方法作为算法的一个遗传算子。 赵明旺[16]采用牛顿迭代算法与遗传算法结合的混合算法求解非线性方程组,其混合策略是: 对子代群体中每个个体以概率Pn进行牛顿算子迭代搜索,产生的新迭代值取代父代加入子代群体。该文认为确定牛顿算子概率Pn的大小需考虑的因素为: 1) 所求解的非线性方程组牛顿法迭代过程的收敛性。若其迭代过程收敛较快,则相应的Pn可取较小值。否则,则取较大值。 2) 预定的繁殖代数。若预定的繁殖代数较大,则相应的Pn可取小一些。否则,取大一些。

进一步地,赵明旺在文献[17]中将上述思想推广运用至求解相容非线性方程组问题,即求解方程数与变量数可能不一致的情况。

罗亚中等[12 -18],屈颖爽[19]均提出了将拟牛顿法与遗传算法相结合的混合算法求解非线性方程组。文献[12]依据自适应混合算子概率Pn对群体利用拟牛顿法进行局部搜索,此文设计了两类混合算子: 拟牛顿迭代混合算子和Powell混合算子,并进行了分析比较。 拟牛顿混合算子的构造方法是: 以被选中个体的染色体为初始点进行拟牛顿迭代操作,若算法收敛则以收敛点作为个体新的染色体,若不收敛则不改变个体的染色体。Powell混合算子的构造方法是: 以个体的染色体为初始点进行Powell寻优计算,以优化结果作为个体新的染色体。而文献[18]的中心思想是只对每一代的最优个体利用拟牛顿法进行精细局部搜索。其文中认为自适应混合算子概率Pn应该随着进化的增加而变大,最后趋近于一个常数P0,因此提出了如下公式

其中,T为遗传算法中设置的最大代数; t为当前进化代数; 常数P0∈( 0,1 ]反映了局部强搜索算子对每个个体的最大可能作用程度。P0大则混合算子对遗传算法搜索空间的局部开采愈充分,但相应的计算成本愈大,此文选择P0= 0. 10。

屈颖爽[19]将知名的拟牛顿法—BFGS方法与遗传算法进行混合,其混合策略与赵明旺[16]的一致,即对子代群体中每个个体以概率Pn进行拟牛顿算子迭代搜索产生的新迭代值取代父代加入子代群体。值得注意的是此文利用有限Markov理论建立了基于基本遗传算法与BFGS算法的混合算法的全局收敛性定理, 即证明了保留当前最好解的混合算法收敛到全局最优解。

王鹏[20]进一步改进了屈颖爽[19]的工作,主要体现在两个方面: ( 1) 拟牛顿迭代步骤中,首先计算每个体的适应值,从中选择性能好的个体按一定的概率进行拟牛顿迭代。( 2) 对个体的拟牛顿迭代仅进行一次,在取得良好效果的同时节省了计算量。

排新颖等[9]考虑将梯度法与遗传算法混合,其主要思想是用梯度对最优个体进行局部寻优。选择的下降方向是待优化函数负梯度方向,设置一个步长L,在当前点进行局部寻优: 若寻优后的函数值小于原来函数值的 λ( λ <1) 倍,则寻优成功,更新当前个体,否则L = Lμ,直到L小于一个设定的值。

曹春红等[21]将几何约束问题转化为非线性方程组的形式,并用牛顿迭代算法与遗传算法结合的混合算法求解,取得了较好的效果,其混合策略与赵明旺[16]一致。

周丽等[22]提出了一种混合小生境遗传算法与拟牛顿算法求解非线性方程组的新方法,其混合策略与罗亚中等[12]一致。对于由小生境遗传算法产生的一代种群,除最优个体外以一定的概率Pn选择其中的个体采用拟牛顿法进行局部寻优,如果寻优后的个体适应度值高于寻优前的适应度值,以优化的结果替换寻优前的个体,否则说明此个体不收敛,不进行替换操作。

4 采用经典迭代方法进行求解

C. L. Karr等[23]提出了一种求解非线性方程组的混合算法,此法中遗传算法主要为牛顿迭代算法提供好的初值。此文以求取Gauss - Legendre求积公式节点与系数形成的非线性方程组为典型例子,说明了牛顿迭代算法对求解此类方程组时对迭代初始值高度敏感,而混合算法的求解效果相当不错。

田巧玉等[14]在实现遗传算法与拟牛顿法相结合的思路是: 在应用遗传算法进行优化设计的基础上,采用拟牛顿法进行二次优化,将获得的结果作为最优的设计结果,即首先运行遗传算法,找到最优点附近的一个点,将它作为拟牛顿法的初始点,以找到目标函数的最小值。此文给出了遗传算法迭代终止条件有两个: 一是进化到指定的最大代数; 另一个是适应度限,即当前代的最佳适应度值小于或等于规定的值时就停止。

吴国辉等[24]考虑到遗传算法全局群体搜索能力及该算法在起始搜索阶段速度非常快的特点,提出先用遗传算法进行n步迭代,直至其收敛结果满足预期要求,进入拟牛顿法的收敛域之内,然后将遗传算法迭代结果作为拟牛顿法迭代初始值继续迭代至精确解。 此文对遗传算法迭代的控制办法是当前群体平均适应度小于等于预先设定的精度值。

4. 1 与其他智能算法的混合

遗传算法是模仿自然选择和遗传机制的一种智能优化算法,研究者已将其与其他智能算法如模拟退火算法等进行混合,或者融入如量子计算原理、生物免疫思想等以构造出更具优势的新的智能算法。

4. 2 遗传退火算法

模拟退火算法是由Metropolis等基于热力学的退火机制提出来的一种对退火过程进行模拟的算法[25]。 遗传算法具有较强的全局搜索能力但容易陷入早熟, 而模拟退火算法具较强的局部把握能力,但全局搜索能力不足。因此对两者取长补短的遗传退火算法应运而生,并成功用于全局优化问题的求解[25 -26]。付振岳等[27]将遗传退火算法进一步加以改进并应用到含有超越函数的非线性方程组的求解中,胡斐等[28]将自适应遗传算法与模拟退火算法结合提出了自适应的遗传退火算法,均取得了较好的求解效果。

4. 3 量子遗传算法

量子遗传算法由Narayanan等人[29]最早提出的, 并经Han等人[30]发展的一种基于量子计算原理的概率优化方法,其用量子位编码来表示染色体,用量子门作用和量子门更新来完成进化搜索,具有种群规模小、 全局寻优能力强的特点,已成功用于求解TSP与0 -1 背包等问题,获得了比传统遗传算法更好的效果。杜娟等[31]提出了一种基于混合量子遗传算法的非线性方程组求解方法。为克服量子遗传算法的求解精度低和加强局部搜索能力,采用拟牛顿法作为一个强局部搜索算子,把量子遗传算法的寻优结果作为拟牛顿法的初值,取得了较好的效果。徐红[15]提出了一种改进量子遗传算法求解非线性方程组的新方法,结果表明此法具有较高的收敛可靠性及较高的精度,对于非线性方程组求解问题具有良好的适应性,特别是一些很难求解的方程组,利用此方法求解更有意义。

4. 3 免疫遗传算法

免疫遗传算法是近年来基于生物免疫机制提出的一种新的智能计算算法[32],其将生物免疫系统的学习、记忆和多样性的特点引入遗传算法。由于兼顾了搜索速度、全局搜索能力和局部搜索能力,免疫遗传算法正成为优化设计领域的研究热点之一。王景芳等[33]提出了一种烃类蒸汽转化反应器物料衡算的非线性方程组的免疫遗传算法求解方法,取得了良好的效果。该算法由于在遗传算法的基础上引入了高斯变异和基于抗体浓度的更新策略调节机制,能有效地保持抗体的多样性,从而避免了遗传算法中存在的早熟收敛问题。

5 并行遗传算法求解非线性方程组

遗传算法具有并行性,其按并行方式搜索一个种群数目的点,而不是单点。其的并行性表现在以下两个方面: 一是内在并行性,本身适应大规模并行。最简单的并行方式是让若干台甚至是数千台计算机各自进行独立种群的演化计算,运行过程中甚至不进行任何通信,等到运算结束时才通信比较,选取最佳个体。这种并行处理方式对并行系统结构并无任何限制和要求,即遗传算法适合在目前所有的并行机或分布式系统上进行并行处理,尤其适用现在PC机多线程多进程的并行计算。二是遗传算法的内含并行性。由于遗传算法采用种群的方式组织搜索,因而可同时搜索解空间的多个区域,并相互交流信息。使用这种搜索方式,虽然每次只执行与种群规模n成比例的计算,但实质上已进行了大约O( n3) 次有效搜索,这就使遗传算法能以较少的计算获得较大的收益。

目前并行遗传算法的实现方案分3 类: 主从式模式( Master - Slave Model) 、粗粒度模型( Coarse - Grained Model) 和细粒度模型( Fine - Grained Model) 。

并行遗传算法正广泛应用于各个领域,已有人将并行遗传算法用在非线性方程组求解上。刘灿文等[34]提出一种自适应并行遗传算法: 采用粗粒度模型,将初始种群分为3 个子种群,并让这3 个子群体按照特点互补的演化规律在并行机的3 个进程上分别进化,定期移植或交换部分优秀个体。其中进程0 为主进程,主要任务之一是保存自身进化而得的优秀个体并吸收另外两个进程进化而得的优秀个体,使不遭受破坏,目的在于保持优秀个体的多样性。对交叉和变异概率针对不同进程、不同阶段设置不同的值,来提高进化效率,如进程0 中的子群体在进化中具有较低的交叉率和变异率,使得优秀个体不易被破坏,而对进程1 和进程2 用较高的交叉率和变异率,尤其在进化后期, 用以探测新的超平面上的优秀个体,不断为进程0 提供新的超平面上的优秀个体,防止进程0“早熟”并加快其收敛到全局最优的速度。

付振岳等[35]针对一些复杂的含有三角函数或是对数的非线性方程组用遗传算法进行求解,由于该类问题涉及的求解空间大,遗传算法需要更大的初始化种群,为提高算法的性能,引入了并发机制和最大堆来提高硬件利用率和降低某些关键步骤的时间复杂度。 根据CPU的利用率和计算等待时间的比率以及CPU的个数,设置合适的线程数,采用细粒度模型,将种群根据线程数分割成相应的组群,每个线程对组群先后进行堆初始化→交叉→堆初始化→变异→堆初始化→ 淘汰等操作,设置一个全局变量记录全体线程在执行遗传退火算法中产生的最优个体,每个线程在完成一轮遗传退火行为后,将本组群的最优值与全局最优值进行比较,若本组群最优值优于全局最优值,则用本组群最优值替换全局最优值,反之,则用全局最优值替换本组群最优值。该算法有效地提高了遗传退火算法的性能,加快了求解的速度。

6 结束语

文中分析了遗传算法求解非线性方程组的研究现状。从所列的文献看,大部分的研究工作主要集中在对遗传算法求解线性方程组各个环节的改进中,也有不少学者提出了遗传算法与其他智能算法混合使用来求解的方法,比如模拟退火算法、量子计算原理等。还有众多研究者结合了计算机多处理器技术和并行运算技术,提出了并行遗传算法这一新概念。从现有遇到的问题看,对该问题的研究尚有不足,需对其进行深一步的研究。

未来的工作,可包括两个方面: 一在现有模型基础上,研究如何将遗传算法与局部优化方法结合,以便有效提高解的质量。另一方面朝并行模型混合化方向研究。可以预期,随着计算机技术的进步和生物学研究的深入,遗传算法求解非线性方程组将更具通用性和有效性。

摘要:针对近年来国内外基于遗传算法求解非线性方程组的研究情况进行了系统总结,揭示了各种遗传策略的综合运用,遗传算法与传统迭代方法或其他智能优化方法的有效结合,并行遗传算法是基于遗传算法方法求解非线性方程组成功的关键,同时对未来的研究方向进行了展望。

非线性微分方程求解 篇8

常用的积分变换有两类:Fourier变换和Laplace变换。当自变量x的取值区间是无穷区间 (-∞, +∞) 时, 可考虑对方程中的x采用Fourier变换, 这时被变换的函数还应满足当x→±∞时趋于0, 而且对二阶的方程还要求被变换的函数的一阶导数也趋于0。当自变量x的取值区间是半无穷区间[0, +∞) , 而且在x=0时给定了函数值以及直到低于方程阶数的各阶导数值时, 可考虑对方程中的x采用Laplace变换, 这时被变换函数还应是指数级的。对一个偏微分方程相继实施这两种变换, 可将偏微分方程化为一个代数方程, 然后再相继实施逆变换, 就可以得到偏微分方程的解的表达式。

1.1 像原函数的微分性质:

1.2 卷积定理

设F1 (s) =L[f1 (t) ], F2 (s) =L[f2 (t) ], 则

Fourier变换的性质

(1) 像原函数的微分性质:设函数f (t) 在 (-∞, +∞) 上连续或只有有限个可去间断点, 当|t|→+∞时, f (n) (t) →0, 则F[f (n) (t) ]= (iω) nF[f (t) ]

(2) 筛选性质

对任意的连续函数f (t) , 有

2 一维热传导方程的求解

undefined

其中M为充分大的正数, 且假定undefined

对t采用Laplace变换, 由 (1) 及 (2) 式结合拉氏变换的微分性质, 得

undefined

其中undefined是U (x, t) 的Laplace变换,

再对 (4) 式中的x采用Fourier变换, 用undefined表示undefined (x, s) 的Fourier变换, 结合付氏变换的微分性质, 得

undefined

(5) 式是一个代数方程, 解得

undefined

对 (6) 式进行拉氏逆变换L-1, 得

undefined

对上式再进行逆变换F-1, 结合付氏变换的筛选性质, 得

undefined

3 一维波动方程的求解

undefined

假定undefined

对t采用Laplace变换, 由 (7) 、 (8) 及 (9) 式结合拉氏变换的微分性质, 得

undefined

其中undefined,

再对 (10) 式中的x采用Fourier变换, 令undefined, 结合付氏变换的微分性质, 得undefined

(11) 式是一个代数方程, 解得

undefined

对 (12) 式进行拉氏逆变换L-1, 结合拉氏变换的卷积性质, 得

undefined

对上式再进行逆变换F-1, 得

undefined

由逆变换的定义, 结合付氏变换的筛选性质, 得

undefined

undefined

从而undefined

4 结语

用积分变换法求解线性偏微分方程, 可将复杂繁冗的计算简化为代数运算, 而且解的表达式比较紧凑, 便于分析各元素之间的相互依赖关系, 而且方法易掌握.

参考文献

[1]张元林.积分变换 (第四版) [M].北京:高等教育出版社, 2003, (5) .

非线性微分方程求解 篇9

关键词:Jordon标准型,特征多项式,广义特征向量,基解矩阵

一、问题的提出

一阶常系数线性齐次微分方程组的一般形式是BX′=DX, 其中B= (bij) n×n, X′= (x1′, x2′, …, xn′) , D= (dij) n×n, X= (x1, x2, …, xn) T.我们仅把X′当作未知向量, 通过求解代数方程组, 可将微分方程组BX′=DX转化为X′=AX的形式, 其中A= (aij) n×n, 当矩阵B是非奇异矩阵时, A=B-1D.本文将对微分X′=AX组的解法进行探讨.

在求解线性微分方程组时, 至关重要的问题是:求出微分方程组的基解矩阵.许多《常微分方程》的教材对用约当 (Jordon) 标准型来求基解矩阵几乎都不作介绍, 更没有从理论上加以讨论.很多代数学方面的文献也只是在如何写出一般方阵的约当 (Jordon) 标准型做了介绍, 而对确定非奇异矩阵T, 使得T-1AT=E的过程却轻描淡写或一带而过.为了简便、通俗地解决这个问题, 本文给出了一种广义特征向量的矩阵变换迭代法.

对任意给定的微分方程组:

X′=AX (1)

其中, X′= (x1′, x2′, …, xn′) T, A= (aij) n×n, X= (x1, x2, …, xn) T.

设其系数矩阵A的特征根为:λ1 (τ1重) , λ2 (τ2重) , …, λr (τr重) , 即矩阵A的特征多项式为:

二、几个概念及相关性质

定义1形如的方阵称为Si阶约当 (Jordon) 块.

定义2由若干个Jordon块组成的分块对角阵, 其中Ji (i=1, 2, …, m) 为Si阶Jordon块, 且当时称为n阶Jordon标准型.

定义3[1]对任意的n阶常数矩阵A, 定义其指数

性质如果T是非奇异矩阵, 则

定义4[2]设λi为矩阵A的一个特征值, 且有

则称τi为λi的代数重数, 一般简称重数;同时, 若kank (λiE-A) =n-ti, 则称ti为λi的几何重数.

定义5[2]若对非零向量αi有:则称αi是矩阵A的属于特征根λi的k阶广义特征向量.

定义6λ是矩阵A的特征值, α是一个非零向量, 则向量组

称为向量α的长度k的广义特征向量链.

三、引理及问题的分析

引理1每一个n阶矩阵A都与一个Jordon型矩阵相似, 这个Jordon型矩阵在不考虑其中的Jordon块的排列次序外是被矩阵A唯一确定的.

与矩阵A相似的Jordon型矩阵称为矩阵A的Jordon标准型.

对于方程组X′=AX中的矩阵A, 由引理1可知必存在非奇异的矩阵T, 使得T-1AT=J (8)

其中J为定义2的Jordon标准型, 即

由于矩阵J与Ji (i=1, 2, …, m) 的特殊形式, 利用定义3很容易得到

这就是微分方程组X′=AX的基解矩阵.

现在关键就是寻求上述的非奇异变换矩阵T, 由 (9) 式可以联想, 关键是寻找Jordan块.

将 (10) 式中的m个Jordan块作适当的分组, 组成r个小Jordan标准型J (τi) 即

不失一般性, 现只对矩阵A的特征值λi进行讨论, 设λi的重数为τi, 几何重数为ti.由于λi的几何重数是ti, 则其线性无关特征向量必有ti个, 设为

引理2设α是矩阵A的属于特征值的λ的k阶广义特征向量, 则α的长度为k的广义特征向量链:

必定是线性无关.

四、定理及其证明

定理1设βi1, βi2, …, βiti为矩阵A属于特征值λi线性无关特征向量, 则以下ti条广义特征向量链中的所有向量都是线性无关的.

为了表达简明, 可对问题分情形讨论.

情形一若li1=li2=…=liti.

现将 (18) 式代入 (16) 式, 然后对 (16) 式两边同时施行矩阵变换 (λiE-A) li1-2, 同理可得k (2) i1=k (2) i2=…=k (2) iti=0.

如此依次操作, 即可得到:ki1 (j) =ki2 (j) =…=kiti (j) =0 (j=3, 4, …, li1) .

这即是说 (16) 式左边各项的系数全为零, 所以, 当li1=li2=…=liti时, 命题结论正确.

情形二若li1, li2, …, liti不全相等.

设在li1, li2, …, liti中从右边到左边的第一个不连等的系数为lih-1, 即lih-1

则对 (16) 式两边同时施行矩阵变换 (λiE-A) lih-1, 得

现设从 (15) 式的右边到左边第二个不连等系数为liq-1, 即liq-1

同样可设liq=liq-1+p, 则可类似前面证明可得:

如此下去, 对从 (16) 左边到右边不连等的序号进行讨论直至结束, 容易得到 (18) 左边的组合系数全为零.

从而, 即 (14) 式定义的向量组线性无关.

综合情形一和情形二, 定理1结论正确.证毕.

要使T-1AT=J成立, 即AT=TJ, 亦即

现对 (20) 第i个块进行分析, 有ATi=TiJτi. (21)

则 (21) 可以分解成ti个等式:

所以是满足 (23) 的, 其中αij (s) = (-1) s-1βij (lij-s+1) (s=1, 2, …, lij) .

别向量改变正负号并不能影响向量组的线性关系性, 则可由定理1以及属于不同特征值的广义特征向量必线性无关知:当时, 得到的矩阵T= (T1, T2, …, Tr) 就是变换T-1AT=J式中的非奇异变换矩阵T.证毕.

五、微分方程组X′=AX的求解算法

综上讨论可以得到用化Jordon标准型法求解微分方程组X′=AX的一般算法步骤:

1. 求出微分方程组X′=AX中的矩阵A的特征多项式并分解成一次因式方幂的乘积 (形如 (2) ) , 而对于一般的多项式想对其进行因式分解是比较困难的, 不妨借助数学软件Mathematica, 使用命令:Roots[g (λ) =0, λ], 求出g (λ) 所有的根, 然后将其分解成:g (λ) = (λ-λ1) τ1 (λ-λ2) τ2… (λ-λr) τr.

2. 对特征值λi的重数τi分成几何重数ti份, 分别为li1, li2, …, liti, 即li1+li2+…+liti=τi, 并且解出ti个属于特征值λi的无关特征向量如 (12) .

3. 分别解ti个线性方程组αij (lij) = (λiE-A) lij-1γij (j=1, 2, …, ti) 的一个非零解.

4. 根据 (14) 式可求出Ti, 进而求出T= (T1, T2, …, Tr) .

5. 根据 (10) 式可以求出微分方程组X′=AX的基解矩阵exp At.

6. 微分方程组X′=AX的解为:X=C exp At, C为任意常列向量.

上面给出的微分方程组的求解算法的步骤明确, 程序性强, 容易通过编程在计算机上实现.

参考文献

[1]王高雄等.常微分方程[M].广州:高等教育出版社, 1983.

[2]方保, 矩阵论[M].北京:清华大学出版社, 2004.

[3]王萼芳等.高等代数[M].北京:高等教育出版社, 2003.

非线性微分方程求解 篇10

人们把碳氢化合物的集合称之为烃, 按相态又将其区分为气态烃、液态烃和固态烃。烃类在一定条件下能与水蒸气进行化学反应, 生成含CO、CO2、CH4、H2等组分的转化气, 人们习惯称该反应为烃类蒸气转化反应。转化气是基本原料, 在合成氨、羟基合成、制氢等方面获得广泛应用, 因而烃类蒸气转化在化肥工业、石油化工等领域占有特定地位。

用于蒸气转化的气态烃有天然气、油田气、炼厂气、焦炉气等;液态烃有石脑油、抽余油等轻质馏分, 为扩大原料来源, 近来已着手研究重油蒸气转化。

烃类蒸气转化热力学及物料、热量计算有两个目的:一是考察烃类蒸气转化反应的限度及其随工艺条件的变化, 为工艺条件选择提供热力学依据;二是做该系统的物料、热量计算, 求出单耗指标, 为方案评比时提供依据。

随着对烃类蒸气转化反应特征及其传递过程特征认识的深化, 借鉴相关学科, 开发转化过程数学模拟系统, 并在实践中加以检验, 得出结构参数及操作参数, 诸如催化剂活性、水碳比、操作负荷、空气过剩系数对反应器状态、操作弹性、能耗、温度分布、浓度分布的影响态势。

采用烃类蒸气转化工艺, 具有能耗低、运转周期长、经济效益高的显著优点, 成为化工、化肥工业的骨干力量。

2 物料衡算的非线性二元方程组

对烃类蒸气转化反应系统, 人们感兴趣的是三个反应, 即一氧化碳变换反应、甲烷转化反应、甲烷裂解化反应。

在进行物料衡算时, 给定两个平衡温距, 增加了两个限制方程, 相当于给定了两个关键组分, 即:

ΚΡms=0.974Ρ2 (3a+4b-4-m2) 3×10-10nΤ2 (1-a-b) (R-a-2b) (1)

ΚΡwgs=b (3a+4b-4-m2) 1a (R-a-2b) (2)

nΤ=1+R+2 (a+b) -4-m2+D+E

式中:nT——达到平衡时系统的总摩尔数, mol;P——总压 (平衡时系统常数) , Pa;M——原料烃中的氢碳比, 无因次;R——原料烃中的水碳摩尔比, 无因次;D——原料烃中的氮碳比;E——原料烃中的氩碳比;KPms——一氧化碳变换反应平衡常数;KPwgs——甲烷转化反应平衡常数。

平衡常数的表达式通式为:

lnΚΡ=Δa0RΤ+Δa1lnΤR+Δa2ΤR+Δa3Τ22R+Δa4Τ33R

+J

式中:T——温度, K;J——反应分商。

在焓多项式系数及标准生成自由焓已知的条件下, 求得:

式 (1) 、式 (2) 中, 仅ab两个未知数, 过去多用试差法、代入法、迭代法求解, 很繁琐。随着计算机的普及, 有的用拟牛顿法求解, 但也费事;在这里我们采用免疫遗传算法, 求解快速而准确且回避陷入局部极值点。

我们在制氢物料衡算、流程优化、装置节能方面, 遇到上述非线性二元方程组 (1) 、 (2) 的求解。在求得ab值后, 可获得:

3 免疫遗传算法

免疫遗传算法 (IGA) 是近年来基于生物免疫机制提出的一种改进遗传算法, 是一种新型的计算智能方法, 它是生命科学中免疫原理与传统遗传算法的结合。生物免疫系统具有抗体多样性、自我调节、免疫记忆功能等特征[1], 免疫遗传算法就是在遗传算法的基础上引入生物免疫系统的基本特性。现在研究和应用都表明, 免疫遗传算法兼顾了搜索速度、全局搜索能力和局部搜索能力, 正成为优化设计领域的研究热点之一[2]。

免疫遗传算法把待求解的问题对应为抗原, 问题的解对应为抗体, 用抗原和抗体的亲和度描述可行解与最优解的逼近程度[3,4]。算法首先接收一个抗原 (对应特定问题) , 然后随机产生一组初始抗体 (对应初始候选解) ;接着计算每一个抗体的适应度 (亲和度) , 对抗体进行交叉和变异;再通过基于浓度的群体更新策略生成下一代抗体群, 直至满足终止条件, 算法结束。免疫遗传算法的流程图见图1。

其基本步骤如下:

(1) 算法初始化。抗原输入及参数设定:输入目标函数及约束条件, 作为抗原的输入;设定种群规模Popsize、选择概率Ps、交叉概率Pc、变异方式等参数;

(2) 初始抗体产生。在第一次迭代时, 抗体通常在解空间中用随机的方法产生;

(3) 亲和度及浓度的计算。计算各抗体和抗原的适应度并计算各抗体的浓度;

(4) 终止条件判断。判断是否满足终止条件, 如果是, 则将与抗原适应度最高的抗体加入免疫记忆数据库中, 然后终止;否则继续;

(5) 选择、交叉、变异操作。根据设置的选择概率Ps、交叉概率Pc和变异方式选择抗体进行选择、交叉和变异操作;

(6) 根据以上的操作更新群体后转到步骤 (3) 。

4 问题求解设计

4.1 亲和度函数设计

在温度T=705 ℃时, 算得:一氧化碳变换反应平衡常数KPms=14.503 9。

在温度T=740 ℃时, 算得:甲烷转化反应平衡常数KPwgs=1.321 1。

在实例中输入:P=9.84e5, m=3.668 3, D=0.011 33, R=3.861 8;代入式 (1) 、式 (2) 得a, b的二元非线性方程组:

g1 (ab) =94.3081 (3a+4b-0.1659) 3 (2a+2b+4.7073) 2 (1-a-b) (3.8618-a-2b) -14.5039=0 (4)

g2 (ab) =b (3a+4b-0.1659) a (3.8618-a-2b) -1.3211=0 (5)

设亲和度函数:

f (ab) =11+ (g1 (ab) ) 2+ (g2 (ab) ) 2 (0f (ab) 1)

亲和度函数图象见图2, 原方程组的求解转变为求亲和度函数最大值f (a, b) =1的解。

4.2 基于浓度的群体更新

为了保证抗体的多样性, 提高算法的全局搜索能力, 采用了一种基于抗体间欧氏距离和适应度来计算抗体相似度和浓度的方法。记抗体xixj的欧氏距离为D (xi, xj) , 适应度分别为f (xi) 和f (xj) , 给定适当常数δ>0, ε>0, 如满足下式:

D (xixj) δ|f (xi) -f (xj) |ε

则称抗体xi与抗体xj相似, 与抗体xi相似的抗体的个数称为抗体xi的浓度, 记为Ci;抗体xi被选择的几率为p (xi) , 即:

p (xi) =αCi[1-f (xi) Μ (x) ]+βf (xi) Μ (x)

式中:α, β—— (0, 1) 之间的可调参数;M (x) ——所有抗体的最大适应度值;Ci——抗体xi的浓度。

可见, 当抗体浓度高时, 适应度高的抗体被选中的几率就小;当抗体浓度不高时, 适应度高的抗体被选中的几率就大。这样既保留了优秀个体, 又可减少相似抗体的选择, 确保了个体的多样性。

4.3 遗传操作

免疫遗传算法能使抗体保持多样性并且最终能够收敛到最优解的主要操作, 就是在算法中有选择、交叉和变异算子的存在, 从而使整个抗体群沿着适应度较好的方向搜索。

(1) 选择算子。采用如下的选择算子:

Ρs (xi) =αρ (xi) Σi=1nρ (xi) + (1-α) 1Νe-Ciβ

式中:ρ (xi) ——适应度函数的矢量距;Ci——抗体xi的浓度;αβ——常数调节因子;N——该种群内抗体的总数。

(2) 两点交叉法。

X1l=[x11, x21, …, xn1], X2l=[x12, x22, …, xn2]是l代的两个抗体, 在第i个点和第d个点实施两点算术交叉, 产生的下一代抗体是:

X1l+1=[x11, …, xi, …, xj, xj+11, …, xn1]

X2l+1=[x12, …, xi, …, xj, xj+12, …, xn2]

其中:xkxk (ikj) 由如下线性组合产生:

xk=ζxk1+ (1-ζ) xk2

xk=ζxk2+ (1-ζ) xk1

式中:ζ∈[0, 1]——比例系数。

(3) 高斯变异法。

采用该法时先将各抗体解码为相应的网络结构, 按下式改变网络的所有权值:

xim=xi+γe-f (xi) μ (0, 1)

式中:xim——变异后的抗体;xi——变异前的抗体;μ (0, 1) ——均值为0, 方差为1的正态分布随机变量;γ∈ (-1, 1) ——个体的变异率;f (xi) ——抗体xi的适应度, 即目标函数的适应值。

上式说明抗体的变异程度和适应度成反比, 即适应度越低 (目标函数的适应值越小) , 个体的变异率越高, 反之亦然。变异后, 重新将组成一个新抗体。

取初始种群规模N=20等参数情况下, 算得:a=0.205 2, b=0.359 1, 适应度f (a, b) =0.999 8, 进而得出系数总摩尔数nT=5.835 9, 并由此算得各组分含量如表1所示。

5 结 论

本文研究基于免疫遗传算法求解烃类转化非线性二元方程组, 取得了良好的效果。免疫遗传算法是借鉴生物免疫系统的自适应识别和排除侵入机体的抗原性异物的功能, 将生物免疫系统的学习、记忆和多样性的特点引入遗传算法。该算法由于在遗传算法的基础上引入了高斯变异和基于抗体浓度的更新策略调节机制, 能有效地保持抗体的多样性, 从而避免了遗传算法中存在的早熟收敛问题[1]。在解决实际问题时, 目标函数和约束条件作为抗原输入, 随后产生初始抗体群, 并通过一系列遗传操作及抗体亲和度的计算, 在保持抗体多样性的情况下, 找出针对该抗原的抗体, 即问题的解。免疫遗传算法的基本特征包括:提高了算法全局搜索能力, 避免陷入到局部最优解;具有最优个体记忆功能;具有快速的全局收敛性能。

摘要:利用免疫遗传算法将烃类蒸气转化反应器物料衡算的非线性二元方程组的求解问题转化为非线性最优化问题, 从而快速地获得生产技术上的有效解。在解决实际问题时, 目标函数和约束条件作为抗原输入, 随后产生初始抗体群, 并通过一系列遗传操作及抗体亲和度的计算, 在保持抗体多样性的情况下, 找出针对该抗原的抗体, 即问题的解。用此方法来解决生产过程中的实际问题, 能优化生产, 提高经济效益。

关键词:免疫遗传算法,烃类蒸气转化,物料衡算,非线性二元方程组,求解

参考文献

[1]谢克明, 郭红波, 谢刚.人工免疫算法及其应用[J].计算机工程与应用, 2005, 29 (18) :77-80.

[2]TIMMIS J, NEALT, HUNTJ.Artificial Immune Systemfor DataAnalysis[J].Biosystems, 2000, 55 (1-3) :143.

[3]DASGUPTA D, FORREST S.Artificial Immune Systems in In-dustrial Applications[C]//Proc 2nd International Conferenceon Intelligent Processing and Manufacturing of Materials.Hono-lulu, 1999:257.

上一篇:排球技术训练下一篇:数据资源安全