《最新在科学计算中的应用9非传统解法选讲PPT课件.ppt》由会员分享,可在线阅读,更多相关《最新在科学计算中的应用9非传统解法选讲PPT课件.ppt(54页珍藏版)》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、在科学计算中的应用在科学计算中的应用9非传非传统解法选讲统解法选讲9.1遗传算法 9.1.1遗传算法及其在最优化问题中的应用 遗传算法是基于进化论,在计算机上模拟生命进化机制而发展起来的一门新学科,它根据适者生存、优胜劣汰等自然进化规则搜索和计算问题的解。 美国 Michigen 大学的 John Holland 于 1975 年提出的。 遗传算法最优化工具箱 MATLAB 7.0的遗传算法与直接搜索工具箱编写函数:function sol,y=c10mga1(sol,options) x=sol(1); y=x.*sin(10*pi*x)+2;调用gaopt( )函数 a,b,c,d=gao
2、pt(-1,2,c10mga1); a,ca = 1.85054746606888 3.85027376676810c = 1.0e+002 * 0.01000000000000 0.01644961385548 0.03624395818177 0.02000000000000 0.01652497353988 0.03647414028140 0.16000000000000 0.01850468596975 0.03850268083951 0.23000000000000 0.01850553961009 0.03850273728228 1.00000000000000 0.0185
3、0547466069 0.03850273766768比较: ff=optimset; ff.Display=iter; x0=1.8; x1=fmincon(f,x0,-1,2,ff); f(x1)ans = -3.85027376676808 f(a(1) % 遗传算法结果ans = -3.85027376676810 ezplot(x*sin(10*pi*x)+2,-1,20) 改变求解区间 a,b,c,d=gaopt(-1,20,c10mga1); a,ca = 19.45005206632863 21.45002604650601c = 1.0e+002 * 0.0100000000
4、0000 0.17243264358456 0.18858649532480 0.02000000000000 0.19253552639304 0.21133759487918 0.25000000000000 0.19450021530572 0.21450017081177 0.27000000000000 0.19450024961756 0.21450018981219 0.29000000000000 0.19450055493368 0.21450025935531 1.00000000000000 0.19450052066329 0.21450026046506 ezplot
5、(x*sin(10*pi*x)+2,12,20) 放大区间 a,b,c,d=gaopt(12,20,c10mga1); a,ca = 19.85005104334383 21.85002552164857c = 1.0e+002 * 0.01000000000000 0.17647930304626 0.19610637643594 0.03000000000000 0.17648091337382 0.19616374074697 0.05000000000000 0.18841858256128 0.202288599115410.21000000000000 0.198500642509
6、44 0.21850023812862 0.23000000000000 0.19850055906254 0.21850025289993 1.00000000000000 0.19850051043344 0.21850025521649 例:求最小值编写函数:function sol,f=c10mga3(sol,options) x=sol(1:4); f=-(x(1)+x(2)2-5*(x(3)-x(4)2-(x(2)-2*x(3)4-10*(x(1)-x(4)4; a,b,c,d=gaopt(-1,1; -1 1; -1 1; -1 1,c10mga3); a,ca = -0.066
7、6 0.0681 -0.0148 -0.0154 -0.0002c =1.0000 -0.3061 0.2075 -0.2235 -0.1206 -0.25805.0000 -0.2294 0.2076 0.0352 -0.1217 -0.125393.0000 -0.0666 0.0682 -0.0148 -0.0154 -0.0002100.0000 -0.0666 0.0681 -0.0148 -0.0154 -0.0002 求解区域太小,有误差GAOT 的最优化函数 其中:p可给目标函数增加附加参数, v为精度及显示控制向量, P0为初始种群, fun1为终止函数的名称,默认值maxG
8、enTerm, n为最大的允许代数。 例:求最小值 tic, xmM=-ones(4,1),ones(4,1)*1000; a,b,c,d=gaopt(xmM,c10mga3,maxGenTerm,2000); a(1:4), dd=c(1:100:end,:); c(end,:), tocans = -0.0049 0.0049 -0.0081 -0.0081dd = 1.0e+009 * 0.0000 0.0000 -0.0000 -0.0000 0.0000 -5.9663 0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000elapsed_time
9、 = 76.5200描述函数:matlab7.0function f=c10mga3a(x) f=(x(1)+x(2)2+5*(x(3)-x(4)2+(x(2)-2*x(3)4+10*(x(1)-x(4)4; x,f=ga(c10mga3a,4) %四个自变量Optimization terminated: maximum number of generations exceeded.x = 0.06976151754582 -0.05491931584170 0.04952579333589 0.06130810339402f = 0.00147647985822 ff=gaoptimset
10、; ff.Generations=2000; ff.PopulationSize=80;ff.CrossoverFcn=crossoverheuristic; x=ga(c10mga3a,4,ff)Optimization terminated: maximum number of generations exceeded.x = -0.00216363106525 0.00216366042770 -0.00039985387788 -0.00039996677375f = 1.739330597649231e-010 f=inline. % 目标函数描述(x(1)+x(2)2+5*(x(3
11、)-x(4)2+(x(2)-2*x(3)4+10*(x(1)-x(4)4,x); 时间少,精度高 ff=optimset; ff.MaxIter=10000; ff.TolX=1e-7; tic, x,f1=fminsearch(f,10*ones(4,1),ff); toc; x,f1Elapsed time is 0.595406 seconds.ans = 1.0e-006 * 0.03039572499758 -0.03039585246164 -0.75343487601326 -0.75343518285272f1 = 9.014052814563438e-024例:求下面的最优化
12、问题 x,y=meshgrid(-1:0.1:3,-3:0.1:3); z=sin(3*x.*y+2)+x.*y+x+y; surf(x,y,z); shading interp % 用光滑曲面表示目标函数函数描述:传统方法 function y=c10mga5(x) y=sin(3*x(1)*x(2)+2)+x(1)*x(2)+x(1)+x(2); x0=1,3; x=fmincon(c10mga5,x0,-1;-3,3;3)x = -1.00000000000000 1.19031291227215函数描述:function sol,y=c10mga6(sol,options) x=sol
13、(1:2); y=-sin(3*x(1)*x(2)+2)-x(1)*x(2)-x(1)-x(2); xmM=-1 3; -3 3; a,b,c,d=gaopt(xmM,c10mga6,maxGenTerm,500); aa = 2.51604948433614 -3.00000000000000 9.00709500762913遗传算法优化中间结果(40代即可,无需500代,可用默认100)9.1.3 遗传算法在有约束最优化问题中的应用 不能直接用于有约束最优化问题求解 需通过变换处理划为无约束最优化问题 对等式约束可通过等式求解将若干个自变量用其它自变量表示。 不等式约束可用惩罚函数方法转移
14、到目标函数中。 仍采用 gaopt() 或 ga() 函数求解 例:描述函数:function sol,y=c10mga4(sol,options) x=sol(1:2); x=x(:); x(3)=(6+4*x(1)-2*x(2)/3; y1=-2 1 1*x; y2=-1 1 0*x; if (y19 | y2-4 | x(3) a,b,c=gaopt(-1000 0; -1000 0,c10mga4,maxGenTerm,1000); c=c(1:15:end,:); c(end,:); a,ca = -6.99981015633155 -10.99962347934527 28.999
15、05078165773c = 1.0e+003 * 0.00100000000000 -0.32769544124065 -0.20423049398177 -0.10000000000000 0.05900000000000 -0.00146223175991 0 0.00131115879955 0.10200000000000 -0.00416116639726 -0.00666729713459 0.01480583198631 0.84900000000000 -0.00689401645967 -0.01080365682806 0.02847008229837 0.8920000
16、0000000 -0.00694511749224 -0.01089232545085 0.02872558746118 0.93200000000000 -0.00698531391213 -0.01097813084259 0.02892656956064 0.96800000000000 -0.00699692906988 -0.01099399300138 0.02898464534940 1.00000000000000 -0.00699981015633 -0.01099962347935 0.02899905078166可用线性规划得出更精确的结果 f=1 2 3; A=-2 1
17、 1; 1 -1 0; B=9; 4; Aeq=4 -2 -3; Beq=-6; x=linprog(f,A,B,Aeq,Beq,-inf;-inf;0,0;0;inf); xOptimization terminated successfully.ans = -6.99999999999967 -10.99999999999935 0.00000000000000 f*xans = -28.99999999999836 建议求解方法:用 GA 找出全局最优解的大致位置,以其为初值调用最优化函数求精确解。9.2神经网络及其在数据拟合中的应用9.2.1神经网络基础知识单个人工神经元的数学表示形式
18、例:常用传输函数曲线 x=-2:0.01:2; y=tansig(x); plot(x,y) x=-2:0.01:2; y=logsig(x); plot(x,y)BP 神经网络结构示意图 其中:xm,xM分别为列向量,为各样本数据的最大最小值。 hi为一行向量,各隐层节点数。fi每层传输函数,同一层应使用相同的传输函数。例: 考虑一个前馈网络,2个隐层,第一个有8个节点,采用Sigmoid传输函数,第二层节点个数应该等于输出信号的路数,故节点数为1,传输函数为对数Sigmoid函数。 net=newff(0,1; -1,5,8,1,tansig,logsig);3个隐层,1层4个点,线性函数
19、;2层6个点, Sigmoid函数;3层1个点, logsig函数。 net=newff(0,1; -1,5,4 6 1,purelin,tansig,logsig);可用下面的语句格式设定其它参数。 net.trainParam.epochs=300; net.trainFcn=trainlm;9.2.2 神经网络的训练与泛化 神经网络训练函数 X为n*M,n为输入变量的路数,M为样本的组数,Y为m*M,m为输出变量的路数。tr为结构体数据,返回训练的相关跟踪信息。Y1和E为输出和误差矩阵。 可多次训练,原加权矩阵为初值。 目标值曲线函数 神经网络泛化例:由前面最小拟合的例子中的数据进行曲线
20、拟合,2个隐层,隐层节点选择为5。 x=0:.5:10; y=0.12*exp(-0.213*x)+0.54*exp(-0.17*x).*sin(1.23*x); x0=0:0.1:10; y0=0.12*exp(-0.213*x0)+0.54*exp(-0.17*x0).*sin(1.23*x0); net=newff(0,10,5,1,tansig,tansig); net.trainParam.epochs=1000; % 设置最大步数 net=train(net,x,y); % 训练神经网络 net.IW1 net.LW2,1 % 隐层权值和输出层权值ans = 0.4765 -1.9
21、076 0.5784 0.9450 -0.2888 -2.7916 0.3052 -2.9388 0.9780 1.1814可改变求解算法 net=newff(0,10,5,1,tansig,tansig); net.trainParam.epochs=100; net.trainFcn=trainlm; net,b1=train(net,x,y); net=newff(0,10,5,1,tansig,tansig); net.trainParam.epochs=100; net.trainFcn=traincgf; net,b2=train(net,x,y); net=newff(0,10,
22、5,1,tansig,tansig); net.trainParam.epochs=100; net.trainFcn=traingdx; net,b3=train(net,x,y);可改变各层传输函数 net=newff(0,10,5,1,tansig,logsig); net.trainParam.epochs=100; net.trainFcn=trainlm; net,b2=train(net,x,y); net=newff(0,10,5,1,logsig,tansig); net,b3=train(net,x,y); net=newff(0,10,5,1,logsig,logsig)
23、; net,b4=train(net,x,y);可改变结构,选择隐层15个节点 net=newff(0,10,15,1,tansig,tansig); net.trainParam.epochs=100; net.trainFcn=trainlm; net,b2=train(net,x,y); figure; y1=sim(net,x0); plot(x0,y0,x0,y1,x,y,o)例:二元函数的拟合 x,y=meshgrid(-3:.6:3, -2:.4:2); x=x(:); y=y(:); z=(x.2-2*x).*exp(-x.2-y.2-x.*y); % 这三个变量均应为行向量
24、net=newff(-3 3; -2 2,10,10,1,tansig,tansig,tansig); net.trainParam.epochs=1000; net.trainFcn=trainlm; net,b=train(net,x; y,z); % 训练神经网络 x2,y2=meshgrid(-3:.1:3, -2:.1:2); x1=x2(:); y1=y2(:); figure; z1=sim(net,x1; y1); z2=reshape(z1,size(x2); surf(x2,y2,z2)改变第二层节点数 net=newff(-3 3; -2 2,10,20,1,tansig
25、,tansig,tansig); net,b=train(net,x; y,z); % 训练神经网络 z1=sim(net,x1; y1); z2=reshape(z1,size(x2); surf(x2,y2,z2) 效果不好给出密集一点的的样本点 x,y=meshgrid(-3:.2:3, -2:.2:2); x=x(:); y=y(:); z=(x.2-2*x).*exp(-x.2-y.2-x.*y); net=newff(-3 3; -2 2,10,10,1,tansig,tansig,tansig); net.trainParam.epochs=100; net.trainFcn=t
26、rainlm; net=train(net,x; y,z); x1,y1=meshgrid(-3:.1:3, -2:.1:2); a=x1; x1=x2(:); y1=y2(:); z1=sim(net,x1; y1); z2=reshape(z1,size(a); surf(x2,y2,z2) net=newff(-3 3; -2 2,10,20,1, tansig,tansig,tansig); net=train(net,x; y,z); % 修改节点个数后的泛化效果 figure; z1=sim(net,x1; y1); z2=reshape(z1,size(a); surf(x2,y
27、2,z2)9.2.3 神经网络界面 启动神经网络界面 nntool 遗传算法的MATLAB编程包括如下几个程序文件: genetic.m (主程序文件) gen_encode.m (二进制数组编码P) gen_decode.m ( 将二进制数组P解码为状态矩阵) crossover.m (两个染色体间的交叉) mutation.m (变异) shuffle.m (打乱染色体次序)function xo,fo = genetic(f,x0,l,u,Np,Nb,Pc,Pm,eta,kmax)% 基因算法求f(x)最小值 s.t. l = x = u%f为待求函数,x0初值,l,u上下限,Np群体大
28、小,Nb每一个变量的基因值(二进制数)%Pc交叉概率,Pm变异概率,eta学习率,kmax最大迭代次数N = length(x0);%确定各变量缺省值if nargin 10 kmax = 100; %最大迭代次数缺省为100endif nargin 1|eta = 0 eta = 1; %学习率eta,(0 eta 1)endif nargin 8 Pm = 0.01; %变异概率缺省0.01endif nargin 7 Pc = 0.5; %交叉概率缺省0.5endif nargin 6 Nb = 8*ones(1,N);%每一变量的基因值(二进制数)endif nargin 5 Np =
29、 10; %群体大小(染色体数)end%生成初始群体NNb = sum(Nb);xo = x0(:); l = l(:); u = u(:);fo = feval(f,xo);X(1,:) = xo;for n = 2:Np X(n,:) = l + rand(size(x0).*(u - l); %初始群体随机数组endP = gen_encode(X,Nb,l,u); %编码为2进制字串for k = 1:kmax X = gen_decode(P,Nb,l,u); %解码为10进制数 for n = 1:Np fX(n) = feval(f,X(n,:); end fxb,nb = mi
30、n(fX); %选择最适合的,函数值最小的 if fxb fo fo = fxb; xo = X(nb,:); end fX1 = max(fxb) - fX; %将函数值转化为非负的适合度值 fXm = fX1(nb); if fXm eps %如果所有的染色体值相同,终止程序 return; end %复制下一代 for n = 1:Np X(n,:) = X(n,:) + eta*(fXm - fX1(n)/fXm*(X(nb,:) - X(n,:); %复制准则 end P = gen_encode(X,Nb,l,u); %对下一代染色体编码 %随机配对/交叉得新的染色体数组 is =
31、 shuffle(1:Np); for n = 1:2:Np - 1 if rand Pc P(is(n:n + 1),:) = crossover(P(is(n:n + 1),:),Nb); end end %变异 P = mutation(P,Nb,Pm);endfunction P = gen_encode(X,Nb,l,u)%将群体X的状态编码为二进制数组PNp=size(X,1); %群体大小N = length(Nb); %变量(状态)维数for n = 1:Np b2 = 0; for m = 1:N b1 = b2+1; b2 = b2 + Nb(m); Xnm =(2Nb(m
32、)- 1)*(X(n,m) - l(m)/(u(m) - l(m); %编码方程 P(n,b1:b2) = dec2bin(Xnm,Nb(m); %10进制转换为2进制 endendfunction X = gen_decode(P,Nb,l,u)% 将二进制数组P解码为群体X的状态矩阵Np = size(P,1); %群体大小N = length(Nb); %变量维数for n = 1:Np b2 = 0; for m = 1:N b1 = b2 + 1; b2 = b1 + Nb(m) - 1; X(n,m) = bin2dec(P(n,b1:b2)*(u(m) - l(m)/(2Nb(m
33、) - 1) + l(m); %解码方程 endendfunction chrms2 = crossover(chrms2,Nb)%两个染色体间的交叉Nbb = length(Nb);b2 = 0;for m = 1:Nbb b1 = b2 + 1; bi = b1 + mod(floor(rand*Nb(m),Nb(m); b2 = b2 + Nb(m); tmp = chrms2(1,bi:b2); chrms2(1,bi:b2) = chrms2(2,bi:b2); chrms2(2,bi:b2) = tmp;endfunction P = mutation(P,Nb,Pm)%变异Nbb
34、 = length(Nb);for n = 1:size(P,1) b2 = 0; for m = 1:Nbb if rand Pm b1 = b2 + 1; bi = b1 + mod(floor(rand*Nb(m),Nb(m); b2 = b2 + Nb(m); P(n,bi) = P(n,bi); end endendfunction is = shuffle(is)%打乱染色体次序N = length(is);for n = N:-1:2 in = ceil(rand*(n - 1); tmp = is(in);is(in) = is(n); is(n) = tmp; %将第n个元素
35、与第in个元素交换end%例:ex.m 遗传算法计算最优解f = inline(x(1)4-16*x(1)2-5*x(1)*x(2)+x(2)4-16*x(2)2-5*x(2),x);l = -5 -5; %下限u = 5 5; %上限x0 = 0 0;Np = 30; %群体大小Nb = 12 12; %代表每个变量的二进制位数Pc = 0.5; %交叉概率Pm = 0.01; %变异概率eta = 0.8; %学习率kmax = 200; %最大迭代次数xos,fos=fminsearch(f,x0)xo_gen,fo_gen = genetic(f,x0,l,u,Np,Nb,Pc,Pm,eta,kmax)