数学实验“线性方程组的最速下降法与共轭梯度法解法”实验报告(内含matlab程序代码)

合集下载

最速下降法和共轭梯度法

最速下降法和共轭梯度法

v b Ax v, v t v, Av x x tv
output k+1, x enddo
算法中也可以用残差作为循环控制,我们可以利用以下的事后误差估计式
|| x x* || || r || ( A) || x || || b ||
证明 按照残差的定义,并且由于假设 x 是方程组的准确解,我们有
(0)
for k=0 to M-1 do
v ( k 1) b Ax ( k ) v( k ) , v( k ) t v ( k ) , Av ( k ) x ( k 1) x ( k ) tk v ( k )
output k+1, x enddo
( k 1)
为了节省空间,实际的算法为 input x ,A, b, M output 0, x for k=0 to M-1 do
定理 1 (A-标准正交系定理)设 A 对称正定, u
(1)
,u
(2)
,…, u
(n)
是一组关于内积
x, y A 的标准正交系。定义 x ( i ) x ( i 1) b Ax (i 1) , u ( i ) u ( i )
其中 x
(0)
是 R 中的任意点,则 Ax
n
( k 1) r (k ) , v(k ) (k ) (k ) x x v (k ) (k ) v Av , ( k 1) b Ax ( k 1) r ( k 1) (i ) k A (i ) v ( k 1) r ( k 1) r (i ) , v v (i ) A i 0 v , v
2
v, b Ax q( x)

matlab数学实验

matlab数学实验

《管理数学实验》实验报告班级姓名实验1:MATLAB的数值运算【实验目的】(1)掌握MATLAB变量的使用(2)掌握MATLAB数组的创建,(3)掌握MA TLAB数组和矩阵的运算。

(4)熟悉MATLAB多项式的运用【实验原理】矩阵运算和数组运算在MA TLAB中属于两种不同类型的运算,数组的运算是从数组元素出发,针对每个元素进行运算,矩阵的运算是从矩阵的整体出发,依照线性代数的运算规则进行。

【实验步骤】(1)使用冒号生成法和定数线性采样法生成一维数组。

(2)使用MA TLAB提供的库函数reshape,将一维数组转换为二维和三维数组。

(3)使用逐个元素输入法生成给定变量,并对变量进行指定的算术运算、关系运算、逻辑运算。

(4)使用MA TLAB绘制指定函数的曲线图,将所有输入的指令保存为M文件。

【实验内容】(1)在[0,2*pi]上产生50个等距采样数据的一维数组,用两种不同的指令实现。

0:(2*pi-0)/(50-1):2*pi 或linspace(0,2*pi,50)(2)将一维数组A=1:18,转换为2×9数组和2×3×3数组。

reshape(A,2,9)ans =Columns 1 through 71 3 5 7 9 11 132 4 6 8 10 12 14Columns 8 through 915 1716 18reshape(A,2,3,3)ans(:,:,1) =1 3 52 4 6ans(:,:,2) =7 9 118 10 12 ans(:,:,3) =13 15 17 14 16 18(3)A=[0 2 3 4 ;1 3 5 0],B=[1 0 5 3;1 5 0 5],计算数组A 、B 乘积,计算A&B,A|B,~A,A= =B,A>B 。

A.*Bans=0 0 15 121 15 0 0 A&Bans =0 0 1 11 1 0 0 A|Bans =1 1 1 11 1 1 1~Aans =1 0 0 00 0 0 1A==Bans =0 0 0 01 0 0 0A>=Bans =0 1 0 11 0 1 0(4)绘制y= 0.53t e -t*t*sin(t),t=[0,pi]并标注峰值和峰值时间,添加标题y= 0.53t e -t*t*sint ,将所有输入的指令保存为M 文件。

用MATLAB实现最速下降法_牛顿法和共轭梯度法求解实例——张小强

用MATLAB实现最速下降法_牛顿法和共轭梯度法求解实例——张小强

机电产品优化设计课程设计报告姓名:张小强学号:201222080633学院:机械电子工程学院实验的题目和要求一.课程名称:最优化设计方法二.实验日期:2013年6月27日三.实验目的:掌握最速下降法,牛顿法和共轭梯度法的算法思想,并能上机编程实现相应的算法。

四.实验要求:用MATLAB 实现最速下降法,牛顿法和共轭梯度法求解实例。

五.实验原理:最速下降法是以负梯度方向最为下降方向的极小化算法,相邻两次的搜索方向是互相直交的。

牛顿法是利用目标函数)(x f 在迭代点k x 处的Taylor 展开式作为模型函数,并利用这个二次模型函数的极小点序列去逼近目标函数的极小点。

共轭梯度法它的每一个搜索方向是互相共轭的,而这些搜索方向k d 仅仅是负梯度方向k g -与上一次搜索方向1-k d 的组合。

五.运行结果如下: 题目:f=(x-2)^2+(y-4)^2①.最速下降法:M 文件:function [R,n]=steel(x0,y0,eps)syms x ;syms y ;f=(x-2)^2+(y-4)^2;v=[x,y];j=jacobian(f,v);T=[subs(j(1),x,x0),subs(j(2),y,y0)];temp=sqrt((T(1))^2+(T(2))^2);x1=x0;y1=y0;n=0;syms kk ;while (temp>eps)d=-T;f1=x1+kk*d(1);f2=y1+kk*d(2);fT=[subs(j(1),x,f1),subs(j(2),y,f2)];fun=sqrt((fT(1))^2+(fT(2))^2);Mini=Gold(fun,0,1,0.00001);x0=x1+Mini*d(1);y0=y1+Mini*d(2);T=[subs(j(1),x,x0),subs(j(2),y,y0)];temp=sqrt((T(1))^2+(T(2))^2);x1=x0;y1=y0;n=n+1;endR=[x0,y0];调用黄金分割法:M文件:function Mini=Gold(f,a0,b0,eps)syms x;format long;syms kk;u=a0+0.382*(b0-a0);v=a0+0.618*(b0-a0);k=0;a=a0;b=b0;array(k+1,1)=a;array(k+1,2)=b;while((b-a)/(b0-a0)>=eps)Fu=subs(f,kk,u);Fv=subs(f,kk,v);if(Fu<=Fv)b=v;v=u;u=a+0.382*(b-a);k=k+1;elseif(Fu>Fv)a=u;u=v;v=a+0.618*(b-a);k=k+1;endarray(k+1,1)=a;array(k+1,2)=b;endMini=(a+b)/2;输入:[R,n]=steel(0,1,0.0001)输出:R = 1.99999413667642 3.99999120501463n = 1②.牛顿法:M文件:syms x1x2;f=(x1-2)^2+(x2-4)^2;v=[x1,x2];df=jacobian(f,v);df=df.';G=jacobian(df,v);epson=1e-12;x0=[0,0]';g1=subs(df,{x1,x2},{x0(1,1),x0(2,1)});G1=subs(G,{x1,x2 },{x0(1,1),x0(2,1)});k=0;mul_count=0;sum_count=0;mul_count=mul_count+12;sum_count=sum_count+6; while(norm(g1)>epson)p=-G1\g1;x0=x0+p;g1=subs(df,{x1,x2},{x0(1,1),x0(2,1)});G1=subs(G,{x1,x2},{x0(1,1),x0(2,1)});k=k+1;mul_count=mul_count+16;sum_count=sum_count+11; end;kx0mul_countsum_count结果:k = 1x0 =24mul_count = 28sum_count = 17③.共轭梯度法:M文件:function f=conjugate_grad_2d(x0,t)x=x0;syms xi yi af=(xi-2)^2+(yi-4)^2;fx=diff(f,xi);fy=diff(f,yi);fx=subs(fx,{xi,yi},x0);fy=subs(fy,{xi,yi},x0);fi=[fx,fy];count=0;while double(sqrt(fx^2+fy^2))>ts=-fi;if count<=0s=-fi;elses=s1;endx=x+a*s;f=subs(f,{xi,yi},x);f1=diff(f);f1=solve(f1);if f1~=0ai=double(f1);elsebreakx,f=subs(f,{xi,yi},x),countendx=subs(x,a,ai);f=xi-xi^2+2*xi*yi+yi^2;fxi=diff(f,xi);fyi=diff(f,yi);fxi=subs(fxi,{xi,yi},x);fyi=subs(fyi,{xi,yi},x);fii=[fxi,fyi];d=(fxi^2+fyi^2)/(fx^2+fy^2);s1=-fii+d*s;count=count+1;fx=fxi;fy=fyi;endx,f=subs(f,{xi,yi},x),count输入:conjugate_grad_2d([0,0],0.0001)结果:x = 0.24998825499785 -0.24999998741273f = 0.12499999986176count = 10ans = 0.12499999986176六.结论如下:最速下降法越接近极小值,步长越小,前进越慢。

最优化方法实验报告(2)

最优化方法实验报告(2)

最优化方法实验报告Numerical Linear Algebra And ItsApplications学生所在学院:理学院学生所在班级:计算数学10-1学生姓名:甘纯指导教师:单锐教务处2013年5月实验三实验名称:无约束最优化方法的MATLAB实现实验时间: 2013年05月10日星期三实验成绩:一、实验目的:通过本次实验的学习,进一步熟悉掌握使用MATLAB软件,并能利用该软件进行无约束最优化方法的计算。

二、实验背景:(一)最速下降法1、算法原理最速下降法的搜索方向是目标函数的负梯度方向,最速下降法从目标函数的负梯度方向一直前进,直到到达目标函数的最低点。

2、算法步骤用最速下降法求无约束问题n R()min的算法步骤如下:xxf,a )给定初始点)0(x ,精度0>ε,并令k=0;b )计算搜索方向)()()(k k x f v -∇=,其中)()(k x f ∇表示函数)(x f 在点)(k x 处的梯度;c )若ε≤)(k v ,则停止计算;否则,从)(k x 出发,沿)(k v 进行一维搜索,即求k λ,使得)(min )()()(0)()(k k k k v x f v x f λλλ+=+≥; d )令1,)()()1(+=+=+k k v x x k k k k λ,转b )。

(二)牛顿法1、算法原理牛顿法是基于多元函数的泰勒展开而来的,它将)()]([-)(1)(2k k x f x f ∇∇-作为搜索方向,因此它的迭代公式可直接写出来:)()]([)(1)(2)()(k k k k x f x f x x ∇∇-=-2、算法步骤用牛顿法求无约束问题n R x x f ∈),(min 的算法步骤如下:a )给定初始点)0(x ,精度0>ε,并令k=0;b )若ε≤∇)()(k x f ,停止,极小点为)(k x ,否则转c );c )计算)()]([,)]([)(1)(2)(1)(2k k k k x f x f p x f ∇∇-=∇--令;d )令1,)()()1(+=+=+k k p x x k k k ,转b )。

MATLAB实现最速下降法_和牛顿法和共轭梯度法

MATLAB实现最速下降法_和牛顿法和共轭梯度法

MATLAB实现最速下降法_和牛顿法和共轭梯度法最速下降法:题目:f=(x-2)^2+(y-4)^2M文件:function [R,n]=steel(x0,y0,eps) syms x;syms y;f=(x-2)^2+(y-4)^2;v=[x,y];j=jacobian(f,v);T=[subs(j(1),x,x0),subs(j(2),y,y0)]; temp=sqrt((T(1))^2+(T(2))^2); x1=x0;y1=y0;n=0;syms kk;while (temp>eps)d=-T;f1=x1+kk*d(1);f2=y1+kk*d(2);fT=[subs(j(1),x,f1),subs(j(2),y,f2)];fun=sqrt((fT(1))^2+(fT(2))^2);Mini=Gold(fun,0,1,0.00001);x0=x1+Mini*d(1);y0=y1+Mini*d(2);T=[subs(j(1),x,x0),subs(j(2),y,y0)];temp=sqrt((T(1))^2+(T(2))^2);x1=x0;y1=y0;n=n+1;endR=[x0,y0]调用黄金分割法:M文件:function Mini=Gold(f,a0,b0,eps) syms x;format long; syms kk;u=a0+0.382*(b0-a0);v=a0+0.618*(b0-a0);k=0;a=a0;b=b0;array(k+1,1)=a;array(k+1,2)=b; while((b-a)/(b0-a0)>=eps) Fu=subs(f,kk,u);Fv=subs(f,kk,v);if(Fu<=Fv)b=v;v=u;u=a+0.382*(b-a);k=k+1;elseif(Fu>Fv)a=u;u=v;v=a+0.618*(b-a);k=k+1;endarray(k+1,1)=a;array(k+1,2)=b; endMini=(a+b)/2;输入:[R,n]=steel(0,1,0.0001)R = 1.99999413667642 3.99999120501463 R = 1.999994136676423.99999120501463 n = 1牛顿法:题目:f=(x-2)^2+(y-4)^2M文件:syms x1 x2;f=(x1-2)^2+(x2-4)^2;v=[x1,x2];df=jacobian(f,v);df=df.';G=jacobian(df,v);epson=1e-12;x0=[0,0]';g1=subs(df,{x1,x2},{x0(1,1),x0(2,1)});G1=subs (G,{x1,x2},{x0(1,1),x0(2,1)});k=0;mul_count=0;sum_count=0;mul_count=mul_count+12;sum_count=sum_count+6; while(norm(g1)>epson) p=-G1\g1;x0=x0+p;g1=subs(df,{x1,x2},{x0(1,1),x0(2,1)});G1=subs(G,{x1,x2},{x0(1,1),x0(2,1)});k=k+1;mul_count=mul_count+16;sum_count=sum_count+11;end;kx0mul_countsum_count结果::k = 1x0 =24mul_count = 28sum_count = 17 共轭梯度法:题目:f=(x-2)^2+(y-4)^2M文件:function f=conjugate_grad_2d(x0,t)x=x0;syms xi yi af=(xi-2)^2+(yi-4)^2; fx=diff(f,xi);fy=diff(f,yi);fx=subs(fx,{xi,yi},x0); fy=subs(fy,{xi,yi},x0); fi=[fx,fy]; count=0;while double(sqrt(fx^2+fy^2))>ts=-fi;if count<=0s=-fi;elses=s1;endx=x+a*s;f=subs(f,{xi,yi},x);f1=diff(f);f1=solve(f1);if f1~=0ai=double(f1);elsebreakx,f=subs(f,{xi,yi},x),count endx=subs(x,a,ai);f=xi-xi^2+2*xi*yi+yi^2;fxi=diff(f,xi);fyi=diff(f,yi);fxi=subs(fxi,{xi,yi},x);fyi=subs(fyi,{xi,yi},x);fii=[fxi,fyi];d=(fxi^2+fyi^2)/(fx^2+fy^2); s1=-fii+d*s;count=count+1;fx=fxi;fy=fyi;endx,f=subs(f,{xi,yi},x),count 输入:conjugate_grad_2d([0,0],0.0001) 结果:x = 0.24998825499785 -0.24999998741273f = 0.12499999986176count = 10ans = 0.12499999986176。

数学实验“线性方程组高斯消去法”实验报告(内含matlab程序)

数学实验“线性方程组高斯消去法”实验报告(内含matlab程序)

西京学院数学软件实验任务书实验一实验报告一、实验名称:线性方程组高斯消去法。

二、实验目的:进一步熟悉理解Guass 消元法解法思路,提高matlab 编程能力。

三、实验要求:已知线性方程矩阵,利用软件求解线性方程组的解。

四、实验原理:消元过程:设0)0(11≠a ,令乘数)0(11)0(11/a a m i i -=,做(消去第i 个方程组的i x )操作1i m ×第1个方程+第i 个方程(i=2,3,.....n )则第i 个方程变为1)1(2)1(2...i n in i b x a x a =++ 这样消去第2,3,。

,n 个方程的变元i x 后。

原线性方程组变为:⎪⎪⎪⎩⎪⎪⎪⎨⎧=++=++=++)1()1(2)1(2)1(2)1(22)1(22)0(1)0(11)0(11... . .... ...n n nn n n n n n b x a x a b x a x a b x a x a 这样就完成了第1步消元。

回代过程:在最后的一方程中解出n x ,得:)1()1(/--=n nn n n n a b x再将n x 的值代入倒数第二个方程,解出1-n x ,依次往上反推,即可求出方程组的解:其通项为3,...1-n 2,-n k /)()1(1)1()1(=-=-+=--∑k kk n k j j k kj k k k a x a bx五、实验内容:function maintest2clcclear allA=[1 3 4;2 4 5;1 4 6];%系数矩阵 b=[1 7 6]'%常数项num=length(b)for k=1:num-1for i=k+1:numif A(k,k)~=0l=A(i,k)/A(k,k); A(i,:)=A(i,:)-A(k,:).*l; b(i)=b(i)-b(k)*l; endendendAb%回代求xx(num)=b(num)/A(num,num);for i=num-1:-1:1sum=0;for j=i+1:numsum=sum+A(i,j)*x(j);endx(i)=(b(i)-sum)/A(i,i);endxEnd六、实验结果:A =1.0000 3.0000 4.0000 0 -2.0000 -3.00000 0 0.5000b =1.00005.00007.5000x =16 -25 15。

用MATLAB实现最速下降法-牛顿法和共轭梯度法求解实例

用MATLAB实现最速下降法-牛顿法和共轭梯度法求解实例

题目和要求最速下降法是以负梯度方向最为下降方向的极小化算法,相邻两次的搜索方向是互相直交的。

牛顿法是利用目标函数)(x f 在迭代点k x 处的Taylor 展开式作为模型函数,并利用这个二次模型函数的极小点序列去逼近目标函数的极小点。

共轭梯度法它的每一个搜索方向是互相共轭的,而这些搜索方向k d 仅仅是负梯度方向k g -与上一次接待的搜索方向1-k d 的组合。

运行及结果如下:最速下降法:题目:f=(x-2)^2+(y-4)^2M 文件:function [R,n]=steel(x0,y0,eps)syms x ;syms y ;f=(x-2)^2+(y-4)^2;v=[x,y];j=jacobian(f,v);T=[subs(j(1),x,x0),subs(j(2),y,y0)];temp=sqrt((T(1))^2+(T(2))^2);x1=x0;y1=y0;n=0;syms kk ;while (temp>eps)d=-T;f1=x1+kk*d(1);f2=y1+kk*d(2);fT=[subs(j(1),x,f1),subs(j(2),y,f2)];fun=sqrt((fT(1))^2+(fT(2))^2);Mini=Gold(fun,0,1,0.00001);x0=x1+Mini*d(1);y0=y1+Mini*d(2);T=[subs(j(1),x,x0),subs(j(2),y,y0)];temp=sqrt((T(1))^2+(T(2))^2);x1=x0;y1=y0;n=n+1;endR=[x0,y0]调用黄金分割法:M文件:function Mini=Gold(f,a0,b0,eps)syms x;format long;syms kk;u=a0+0.382*(b0-a0);v=a0+0.618*(b0-a0);k=0;a=a0;b=b0;array(k+1,1)=a;array(k+1,2)=b;while((b-a)/(b0-a0)>=eps)Fu=subs(f,kk,u);Fv=subs(f,kk,v);if(Fu<=Fv)b=v;v=u;u=a+0.382*(b-a);k=k+1;elseif(Fu>Fv)a=u;u=v;v=a+0.618*(b-a);k=k+1;endarray(k+1,1)=a;array(k+1,2)=b;endMini=(a+b)/2;输入:[R,n]=steel(0,1,0.0001)R = 1.99999413667642 3.99999120501463 R = 1.99999413667642 3.99999120501463 n = 1牛顿法:题目:f=(x-2)^2+(y-4)^2M文件:f=(x1-2)^2+(x2-4)^2;v=[x1,x2];df=jacobian(f,v);df=df.';G=jacobian(df,v);epson=1e-12;x0=[0,0]';g1=subs(df,{x1,x2},{x0(1,1),x0(2,1)});G1=subs(G,{x1,x2},{x0(1,1),x0(2,1)});k=0;mul_count= 0;sum_count=0;mul_count=mul_count+12;sum_count=sum_count+6;while(norm(g1)>epson)p=-G1\g1;x0=x0+p;g1=subs(df,{x1,x2},{x0(1,1),x0(2,1)});G1=subs(G,{x1,x2},{x0(1,1),x0(2,1)});k=k+1;mul_count=mul_count+16;sum_count=sum_count+11;end;kx0mul_countsum_count结果::k = 1x0 =24mul_count = 28sum_count = 17共轭梯度法:题目:f=(x-2)^2+(y-4)^2M文件:function f=conjugate_grad_2d(x0,t)x=x0;syms xi yi af=(xi-2)^2+(yi-4)^2;fx=diff(f,xi);fy=diff(f,yi);fx=subs(fx,{xi,yi},x0);fy=subs(fy,{xi,yi},x0);count=0;while double(sqrt(fx^2+fy^2))>ts=-fi;if count<=0s=-fi;elses=s1;endx=x+a*s;f=subs(f,{xi,yi},x);f1=diff(f);f1=solve(f1);if f1~=0ai=double(f1);elsebreakx,f=subs(f,{xi,yi},x),countendx=subs(x,a,ai);f=xi-xi^2+2*xi*yi+yi^2;fxi=diff(f,xi);fyi=diff(f,yi);fxi=subs(fxi,{xi,yi},x);fyi=subs(fyi,{xi,yi},x);fii=[fxi,fyi];d=(fxi^2+fyi^2)/(fx^2+fy^2);s1=-fii+d*s;count=count+1;fx=fxi;fy=fyi;endx,f=subs(f,{xi,yi},x),count输入:conjugate_grad_2d([0,0],0.0001)结果:x = 0.24998825499785 -0.24999998741273 f = 0.12499999986176count = 10ans = 0.12499999986176结论如下:最速下降法越接近极小值,步长越小,前进越慢。

研究报告用matlab解线性方程组

研究报告用matlab解线性方程组

用matlab解线性方程组2008-04-12,17:00一。

高斯消去法1.顺序高斯消去法直接编写命令文件a=[]d=[]'[n,n]=size(a);c=n+1a(:,c)=d;,for,k=1:n-1a(k+1:n,,k:c)=a(k+1:n,,k:c)-(a(k+1:n,k)/,a(k,k))*a(k,,k:c);,,,,,%消去endx=[0,0,0,0]',,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,%回带x(n)=a(n,c)/a(n,n);for,g=n-1:-1:1x(g)=(a(g,c)-a(g,g+1:n)*x(g+1:n))/a(g,g)end2.列主高斯消去法*由于“[r,m]=max(abs(a(k:n,k)))”返回的行是“k:n,k”内的第几行,所以要通过修正来把m,改成真正的行的值。

该程序只是演示程序,真正机器计算不需要算主元素所在列以下各行应为零的值。

直接编写命令文件a=[]d=[],'[n,n]=size(a);c=n+1a(:,c)=d;,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,%(增广)for,k=1:n-1[r,m]=max(abs(a(k:n,k)));,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,%选主m=m+k-1;,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,%(修正操作行的值),,,if(a(m,k)~=0),,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,if(m~=k)a([k,m],:)=a([m,k],:);,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,%换行enda(k+1:n,,k:c)=a(k+1:n,,k:c)-(a(k+1:n,k)/,a(k,k))*a(k,,k:c);,,,,,,%消去endendx=[0,0,0,0]',,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,%回带x(n)=a(n,c)/a(n,n);for,g=n-1:-1:1x(g)=(a(g,c)-a(g,g+1:n)*x(g+1:n))/a(g,g)end3.分别用顺序高斯消去法和列主高斯消去法解方程组a*x=d,并比较结果a=[0,1,2,3;9,11,23,34;62.5,23.4,15.5,17.2;192.01,124,25.1,59.3],d=[1;1;1;1]顺序高斯消去法:提示“Warning:,Divide,by,zero.”,x,=NaN,NaN,NaN,NaN 列主高斯消去法:x,=-1.2460,2.8796,5.5206,-4.3069由此可见列主高斯消去法可以解决顺序高斯消去法所不能解决的问题。

数学实验线性方程组最速下降法与共轭梯度法解法实验报告范文内含matlab程序代码_

数学实验线性方程组最速下降法与共轭梯度法解法实验报告范文内含matlab程序代码_

数学实验线性方程组最速下降法与共轭梯度法解法实验报告范文内含matlab程序代码_西京学院数学软件实验任务书课程名称数学软件实验班级数0901For学号07姓名实验课题线性方程组的最速下降法与共轭梯度法Forperonalueonlyintudyandreearch;notfor实验目的熟悉线性方程组的最速下降法与共轭梯度法实验要求等其中Matlab/C/C++/Java/Maple/Mathematica运用一种语言完成实验内容线性方程组的最速下降法线性方程组的共轭梯度法成绩教师实验五实验报告一、实验名称:最速下降法与共轭梯度法解线性方程组。

二、实验目的:进一步熟悉理解掌握最速下降法与共轭梯度法解法思路,提高matlab编程能力。

三、实验要求:已知线性方程矩阵,应用最速下降与共轭梯度法在相关软件编程求解线性方程组的解。

四、实验原理:1.最速下降法:从某个初始点出发,沿在点处的负梯度方向)(0)(0某某)(某f)0)((0)(0A某)f(某br求得的极小值点,即)(1某)某f()(0(0))rfmin(某0然后从出发,重复上面的过程得到。

如此下去,))2((1某某得到序列{}为A的最小,最大特征值。

最速下降法迭代格式:给定初值,)0(某按如下方法决定:)(k某(k)k(k)A某)f(某rb(k)T(k)r,r(k1)(k)(k)r某某kk(k)T)(k,Arr2.共轭梯度法其基本步骤是在点处选取搜索方向,使其与前)(k)k(d某一次的搜索方向关于共轭,即)1(kdA1)(k1)(k)(k0dd,Ad然后从点出发,沿方向求得的极小值点)(k)k(d某)(某f,即)(k1某)k((k1))f(某(某d)minf0如此下去,得到序列{}。

不难求得的)1k(k)()(k0Ad,d某解为)1k(k)((k)d)df(某1k由共轭的定义可得:)1(k)k(0Ad,d(k)(k1),rAd1k(k1)(k1)d,Ad共轭梯度法的计算过程如下:第一步:取初始向量,计算)(0某(0)(0)(0)(0)某Abf(某dr)(0)(0),Adr0(0)(0)d,Ad(1)(0)(0)d某某0第步:计算1k(k)(k)(k)某A)rbf(某(k)(k1)r,Adk1(k1)(k1),dAd)1(k)(k)k(ddr1k(k)(k)r,Adk(k)(k),Add(k1 )(k)(k)d某某0五、实验内容:%最速下降法function[某,k]=fatet(A,b,ep);某0=zero(ize(b),1);某=某0;k=0;m=1000;tol=1;whiletol>=epr=b-A某某0;q=dot(r,r)/dot(A某r,r);某=某0+q某r;k=k+1;tol=norm(某-某0);某0=某;ifk>=mdip('迭代次数太多,可能不收敛!'); ;returnendend某k%共轭梯度法function[k,某]=gong_e(A,b)ep=input('请输入允许误差ep=');某0=input('请输入初始值某0=');k=0;r0=b-A某某0;%求出dangqian梯度whilenorm(r0)>epr0=b-A某某0;k=k+1;ifk==1p0=r0;elelamda=(r0'某r0)/(p0'某A某p0);r1=r0-lamda某A某p0;p0=r0+(r0'某r0)/(r1'某r1)某p0;某1=某0+lamda某p0;某0=某1;r0=r1;endend某=r0;k;end六、实验结果:A=[520;641;125];b=[1018-14]';ep=1.0e-6;某=-0.87507.1875-5.5000k=60仅供个人用于学习、研究;不得用于商业用途。

MATLAB实现最速下降法_和牛顿法和共轭梯度法

MATLAB实现最速下降法_和牛顿法和共轭梯度法

MATLAB实现最速下降法_和牛顿法和共轭梯度法最速下降法:题目:f=(x-2)^2+(y-4)^2M文件:function [R,n]=steel(x0,y0,eps) syms x;syms y;f=(x-2)^2+(y-4)^2;v=[x,y];j=jacobian(f,v);T=[subs(j(1),x,x0),subs(j(2),y,y0)]; temp=sqrt((T(1))^2+(T(2))^2); x1=x0;y1=y0;n=0;syms kk;while (temp>eps)d=-T;f1=x1+kk*d(1);f2=y1+kk*d(2);fT=[subs(j(1),x,f1),subs(j(2),y,f2)];fun=sqrt((fT(1))^2+(fT(2))^2);Mini=Gold(fun,0,1,0.00001);x0=x1+Mini*d(1);y0=y1+Mini*d(2);T=[subs(j(1),x,x0),subs(j(2),y,y0)];temp=sqrt((T(1))^2+(T(2))^2);x1=x0;y1=y0;n=n+1;endR=[x0,y0]调用黄金分割法:M文件:function Mini=Gold(f,a0,b0,eps) syms x;format long; syms kk;u=a0+0.382*(b0-a0);v=a0+0.618*(b0-a0);k=0;a=a0;b=b0;array(k+1,1)=a;array(k+1,2)=b; while((b-a)/(b0-a0)>=eps) Fu=subs(f,kk,u);Fv=subs(f,kk,v);if(Fu<=Fv)b=v;v=u;u=a+0.382*(b-a);k=k+1;elseif(Fu>Fv)a=u;u=v;v=a+0.618*(b-a);k=k+1;endarray(k+1,1)=a;array(k+1,2)=b; endMini=(a+b)/2;输入:[R,n]=steel(0,1,0.0001)R = 1.99999413667642 3.99999120501463 R = 1.999994136676423.99999120501463 n = 1牛顿法:题目:f=(x-2)^2+(y-4)^2M文件:syms x1 x2;f=(x1-2)^2+(x2-4)^2;v=[x1,x2];df=jacobian(f,v);df=df.';G=jacobian(df,v);epson=1e-12;x0=[0,0]';g1=subs(df,{x1,x2},{x0(1,1),x0(2,1)});G1=subs (G,{x1,x2},{x0(1,1),x0(2,1)});k=0;mul_count=0;sum_count=0;mul_count=mul_count+12;sum_count=sum_count+6; while(norm(g1)>epson) p=-G1\g1;x0=x0+p;g1=subs(df,{x1,x2},{x0(1,1),x0(2,1)});G1=subs(G,{x1,x2},{x0(1,1),x0(2,1)});k=k+1;mul_count=mul_count+16;sum_count=sum_count+11;end;kx0mul_countsum_count结果::k = 1x0 =24mul_count = 28sum_count = 17 共轭梯度法:题目:f=(x-2)^2+(y-4)^2M文件:function f=conjugate_grad_2d(x0,t)x=x0;syms xi yi af=(xi-2)^2+(yi-4)^2; fx=diff(f,xi);fy=diff(f,yi);fx=subs(fx,{xi,yi},x0); fy=subs(fy,{xi,yi},x0); fi=[fx,fy]; count=0;while double(sqrt(fx^2+fy^2))>ts=-fi;if count<=0s=-fi;elses=s1;endx=x+a*s;f=subs(f,{xi,yi},x);f1=diff(f);f1=solve(f1);if f1~=0ai=double(f1);elsebreakx,f=subs(f,{xi,yi},x),count endx=subs(x,a,ai);f=xi-xi^2+2*xi*yi+yi^2;fxi=diff(f,xi);fyi=diff(f,yi);fxi=subs(fxi,{xi,yi},x);fyi=subs(fyi,{xi,yi},x);fii=[fxi,fyi];d=(fxi^2+fyi^2)/(fx^2+fy^2); s1=-fii+d*s;count=count+1;fx=fxi;fy=fyi;endx,f=subs(f,{xi,yi},x),count 输入:conjugate_grad_2d([0,0],0.0001) 结果:x = 0.24998825499785 -0.24999998741273f = 0.12499999986176count = 10ans = 0.12499999986176。

实验一用matlab求解线性方程组

实验一用matlab求解线性方程组

实验1.1 用matlab 求解线性方程组第一节 线性方程组的求解 一、齐次方程组的求解rref (A ) %将矩阵A 化为阶梯形的最简式null (A ) %求满足AX =0的解空间的一组基,即齐次线性方程组的基础解系【例1】 求下列齐次线性方程组的一个基础解系,并写出通解:我们可以通过两种方法来解: 解法1:>> A=[1 -1 1 -1;1 -1 -1 1;1 -1 -2 2]; >> rref(A) 执行后可得结果: ans=1 -1 0 0 0 0 -1 1 0 0 0 0 由最简行阶梯型矩阵,得化简后的方程⎪⎩⎪⎨⎧=+--=+--=-+-02200432143214321x x x x x x x x x x x x取x2,x4为自由未知量,扩充方程组为即提取自由未知量系数形成的列向量为基础解系,记所以齐次方程组的通解为解法2: clearA=[1 -1 1 -1;1 -1 -1 1;1 -1 -2 2];B=null(A, 'r') % help null 看看加个‘r ’是什么作用,若去掉r ,是什么结果?执行后可得结果: B=1 0 1 0 0 1 0 1⎩⎨⎧=-=-004321x x x x ⎪⎪⎩⎪⎪⎨⎧====44432221x x x x x x x x ⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡+⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡11000011424321x x x x x x ,00111⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=ε,11002⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=ε2211εεk k x +=易见,可直接得基础解系所以齐次方程组的通解为二、非齐次线性方程组的求解 Matlab 命令的基本格式:X =A\b %系数阵A 满秩时,用左除法求线性方程组AX =b 的解注意:A/B 即为AB -1, 而A\B 即为A -1B.C =[A,b];D =rref(C) % 求线性方程组AX =b 的特解,即D 的最后一列元素【例2】 求下列非齐次线性方程组的解:,00111⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=ε,11002⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=ε⎪⎪⎪⎩⎪⎪⎪⎨⎧=+=++=++=++=+150650650651655454354332121x x x x x x x x x x x x x 2211εεk k x +=解: clearA=[5 6 0 0 0;1 5 6 0 0;0 1 5 6 0;0 0 1 5 6;0 0 0 1 5]; b=[1;0;0;0;1];format rational %采用有理数近似输出格式,比较format short 看看x=A\b执行后可得所求方程组的解. 作业:【第一题】 求下列非齐次线性方程组的通解.A=[1 2 3 1;1 4 6 2;2 9 8 3;3 7 7 2] B=[3;2;7;12] format rational x=A\B x =⎪⎪⎩⎪⎪⎨⎧=+++=+++=+++=+++1227737389222643324321432143214321x x x x x x x x x x x x x x x x42/31/2684838239393950-7/3【第二题】计算工资问题一个木工,一个电工,一个油漆工,三个人相互同意彼此装修他们自己的房子。

MATLAB实现最速下降法(梯度)程序

MATLAB实现最速下降法(梯度)程序

matlab最速下降法2010-08-18 17:13function x=fsxsteep(f,e,a,b)% fsxsteep函数最速下降法% x=fsxsteep(f,e,a,b)为输入函数 f为函数 e为允许误差 (a,b)为初始点;% fsx TJPU 2008.6.15x1=a;x2=b;Q=fsxhesse(f,x1,x2);x0=[x1 x2]';fx1=diff(f,'x1'); %对x1求偏导数fx2=diff(f,'x2'); %对x2求偏导数g=[fx1 fx2]'; %梯度g1=subs(g); %把符号变量转为数值d=-g1;while (abs(norm(g1))>=e)t=(-d)'*d/((-d)'*Q*d);t=(-d)'*d/((-d)'*Q*d); %求搜索方向x0=x0-t*g1; %搜索到的点v=x0;a=[1 0]*x0;b=[0 1]*x0;x1=a;x2=b;g1=subs(g);d=-g1;end;x=v;function x=fsxhesse(f,a,b)% fsxhesse函数求函数的hesse矩阵;% 本程序仅是简单的求二次函数的hesse矩阵!;% x=fsxhesse(f)为输入函数 f为二次函数 x1,x2为自变量;% fsx TJPU 2008.6.15x1=a;x2=b;fx=diff(f,'x1'); %求f对x1偏导数fy=diff(f,'x2'); %求f对x2偏导数fxx=diff(fx,'x1'); %求二阶偏导数对x1再对x1fxy=diff(fx,'x2'); %求二阶偏导数对x1再对x2fyx=diff(fy,'x1'); %求二阶偏导数对x2再对x1fyy=diff(fy,'x2'); %求二阶偏导数对x2再对x2fxx=subs(fxx); %将符号变量转化为数值fxy=subs(fxy);fyx=subs(fyx);fyy=subs(fyy);x=[fxx,fxy;fyx,fyy]; %求hesse矩阵syms x1 x2;X=[x1,x2];fx=X(1)^2+2*X(2)^2;z=fsxsteep(fx,0.001,1,1)。

计算方法——共轭梯度法求解线性方程组

计算方法——共轭梯度法求解线性方程组
k T
n-1 do
k
r r ; k T k d Ad
x
k 1
x k d ;
k k
1
1 共轭梯度法求解线性方程组

r
k 1
b Ax
k 1
k 1

④ 若 || r k 1 || 或 k 1 n ,则输出近似解 x(k+1),停止;否则,转⑤; ⑤ ⑥
18 9 27 9 1 18 45 0 45 2 ,b A 9 16 0 126 9 27 45 9 135 8
由于该方程组的系数矩阵以及右端项都比较简单,因此采用从 matlab 命令窗口手 动输入的方式来输入数据,取计算精度为 10-6,运行过程及结果如图 2 和图 3(由于迭 代的初始值为随机产生,因此每次得到的残量图会有所不同,但最终都趋于 0) :
|| r ||2 k k 22 ; || r ||2
d
k 1
r
k 1
k d ;
k
end do
图 1 共轭梯度法求解线性方程组程序框图
1.2 程序使用说明 共轭梯度法求解线性方程组的 matlab 程序见附录 1,该程序可以求解系数矩阵为 对称正定矩阵的线性方程组。在使用该程序时,可将程序复制到 matlab 命令窗口中直 接运行或者复制到编辑窗口中保存运行,运行时刻根据提示输入,直至得到结果。 开始运行程序时,会出现提示“请选择系数矩阵、右端项以及系数矩阵阶数的输
图 4 命令窗口显示的运行结果
1.6 1.4 1.2 1
残量
0.8 0.6 0.4 0.2 0 0 20 40 60 80 100

共轭梯度法实验报告

共轭梯度法实验报告

数值代数实验报告一、实验名称:用共轭梯度法解线性方程组。

二、实验目的:进一步熟悉理解掌握共轭梯度法解法思路,提高matlab 编程能力。

三、实验要求:已知线性方程矩阵,应用共轭梯度法在相关软件编程求解线性方程组的解。

四、实验原理:1.共轭梯度法:考虑线性方程组Ax b =的求解问题,其中A 是给定的n 阶对称正定矩阵,b 是给定的n 维向量,x 是待求解的n 维向量.为此,定义二次泛函()2T T x x Ax b x ϕ=-.定理1 设A 对称正定,求方程组Ax b =的解,等价于求二次泛函()x ϕ的极小值点. 定理1表明,求解线性方程组问题就转化为求二次泛函()x ϕ的极小值点问题.求解二次函数极小值问题,通常好像盲人下山那样,先给定一个初始向量0x ,确定一个下山方向0p ,沿着经过点0x 而方向为0p 的直线00x x p α=+找一个点1000x x p α=+,使得对所有实数α有()()00000x p x p ϕαϕα+≤+,即在这条直线上1x 使()x ϕ达到极小.然后从1x 出发,再确定一个下山的方向1p ,沿着直线11x x p α=+再跨出一步,即找到1α使得()x ϕ在2111x x p α=+达到极小:()()11111x p x p ϕαϕα+≤+.重复此步骤,得到一串012,,,ααα 和 012,,,p p p ,称k p 为搜索方向,k α为步长.一般情况下,先在k x 点找下山方向k p ,再在直线k k x x p α=+上确定步长k α使()(),k k k k k x p x p ϕαϕα+≤+最后求出1k k k k x x p α+=+.然而对不同的搜索方向和步长,得到各种不同的算法.由此,先考虑如何确定k α.设从k x 出发,已经选定下山方向k p .令 ()()k k f x p αϕα=+()()()2TT k k k k k k x p A x p b x p ααα=++-+()22TT k k k k k p Ap r p x ααϕ=-+,其中k k r b Ap =-.由一元函数极值存在的必要条件有()220TT k k k k f p Ap r p αα'=-=所确定的α即为所求步长k α,即 T k kk Tk kr p p Ap α=. 步长确定后,即可算出1k k k k x x p α+=+.此时,只要0T k k r p ≠,就有()()()()1k k k k k k x x x p x ϕϕϕαϕ+-=+-()2220T k k TT k kk k k k Tk kr p p Ap r p p Ap αα=-=-<即()()1k k x x ϕϕ+<.再考虑如何确定下山方向k p .易知负梯度方向是()x ϕ减小最快的方向,但简单分析就会发现负梯度方向只是局部最佳的下山方向,而从整体来看并非最佳.故采用新的方法寻求更好的下山方向——共轭梯度法. 下面给出共轭梯度法的具体计算过程:给定初始向量0x ,第一步仍选用负梯度方向为下山方向,即00p r =,于是有00010001000,,T T r r x x p r b Ax p Ap αα==+=-.对以后各步,例如第k+1步(k ≥1),下山方向不再取k r ,而是在过点由向量k r 和1k p -所张成的二维平面21{|,,}k k k x x x r p R πξηξη-==++∈内找出使函数ϕ下降最快的方向作为新的下山方向k p .考虑ϕ在2π上的限制:()1,()k k k x r p ψξηϕξη-=++11()()T k k k k k k x r p A x r p ξηξη--=++++12()T k k k b x r p ξη--++. 计算ψ关于,ξη的偏导得: ()()11112,2,T T T k k k k k kT T k k k k r Ar r Ap r r r Ap p Ap ψξηξψξηη----∂=+-∂∂=+∂其中最后一式用到了10T k k r p -=,这可由k r 的定义直接验证.令 0ψψξη∂∂==∂∂, 即知ϕ在2π内有唯一的极小值点001k k k x x r p ξη-=++,其中0ξ和0η满足 00101011,0.T T T k k k k k k T Tk k k k r Ar r Ap r r r Ap p Ap ξηξη----⎧+=⎨+=⎩由于0k r ≠必有00ξ≠,所以可取()01001k k k k p x x r p ηξξ-=-=+作为新的下山方向.显然,这是在平面2π内可得的最佳下山方向.令010k ηβξ-=,则可得1111.T k k k T k k r Ap p Ap β----=-注:这样确定的k p 满足10Tkk p Ap -=,即k p 与1k p -是相互共轭的. 总结上面的讨论,可得如下的计算公式:T k kk Tk kr p p Ap α= , 1k k k k x x p α+=+, 11k k r b Ax ++=-,1T k kk Tk kr Ap p Ap β+=-, 11k k k k p r p β++=+. 在实际计算中,常将上述公式进一步简化,从而得到一个形式上更为简单而且对称的计算公式.首先来简化1k r +的计算公式:11()k k k k k k k k r b Ax b A x p r Ap αα++=-=-+=-.因为k Ap 在计算k α是已经求出,所以计算1k r +时可以不必将1k x +代入方程计算,而是从递推关系1k k k r b Ap α+=-得到.再来简化k α和k β的计算公式.此处需要用到关系式1110,T T T k k k k k k r r r p r p +-+=== 1,2,k =.从而可导出1111,T T k k k k r r r α+++=-, ()111TT T k k k k k k k k k p Ap p r r p r αα+=-=()1111T Tk k k k k k k kr r p r r βαα--=+=.由此可得,T k k k T k k r r p Ap α=, 11.T k k k T k kr r r r β++=.从而有求解对称正定方程组的共轭梯度法算法如下:0x =初值00r b Ax =-;0k =while 0k r ≠1k k =+ if 1k = 00p r =else21122T T k k k k k r r r r β-----= 1122k k k k p r p β----=+ end11111T Tk k k k k r r p Ap α-----=111k k k k x x p α---=+111k k k k r r Ap α---=-endk x x =注:该算法每迭代一次仅需要使用系数矩阵A 做一次矩阵向量积运算. 定理2 由共轭梯度法得到的向量组{}i r 和{}i p 具有如下基本性质: (1)0T i j p r =, 0;i j k ≤<≤ (2)0T i j r r =, i j ≠,0,;i j k ≤≤ (3)0T i j p Ap =, i j ≠,0,;i j k ≤≤ (4)000{,,}{,,}(,,1)k k span r r span p p A r k κ==+,其中0000(,,1){,,,}k A r k span r Ar A r κ+=,通常称之为Krylov 子空间.下面给出共轭梯度法全局最优性定理:定理3 用共轭梯度法计算得到的近似解k x 满足()(){}00min :(,,)k x x x x A r k ϕϕκ=∈+或{}**00min :(,,)k AAx x x x x x A r k κ-=-∈+,其中Ax=*x 是方程组Ax b =的解,0(,,)A r k κ是由所定义的Krylov 子空间. 定理2表明,向量组0,,k r r 和0,,k p p 分别是Krylov 子空间0(,,1)A r k κ+的正交基和共轭正交基.由此可知,共轭梯度法最多n 步便可得到方程组的解*x .因此,理论上来讲,共轭梯度法是直接法.然而实际使用时,由于误差的出现,使k r 之间的正交性很快损失,以致于其有限步终止性已不再成立.此外,在实际应用共轭梯度法时,由于一般n 很大,以至于迭代()O n 次所耗费的计算时间就已经使用户无法接受了.因此,实际上将共轭梯度法作为一种迭代法使用,而且通常是k r 是否已经很小及迭代次数是否已经达到最大允许的迭代次数max k 来终止迭代.从而得到解对称正定线性方程组的实用共轭梯度法,其算法如下:x =初值0;k =;r b Ax =-T r r ρ=while)()max2band k kε><1k k =+if 1k = p r = else;p r p βρρβ==+ end;;T Ap p x x p ωαρωα===+ ;;T r r r r αωρρρ=-== end算法中,系数矩阵A 的作用仅仅是用来由已知向量p 产生向量Ap ω=,这不仅可以充分利用A 的稀疏性,而且对某些提供矩阵A 较为困难而由已知向量p 产生向量Ap ω=又十分方便的应用问题是十分有益的。

最优化牛顿法最速下降法共轭梯度法matlab代码

最优化牛顿法最速下降法共轭梯度法matlab代码

最优化⽜顿法最速下降法共轭梯度法matlab代码⽜顿法迭代公式:(1)2()1()[()]()k k k k x x f x f x +-=-??Matlab 代码:function [x1,k] =newton(x1,eps)hs=inline('(x-1)^4+y^2'); 写⼊函数ezcontour(hs,[-10 10 -10 10]); 建⽴坐标系hold on; 显⽰图像syms x y 定义变量f=(x-1)^4+y^2; 定义函数grad1=jacobian(f,[x,y]); 求f 的⼀阶梯度grad2=jacobian(grad1,[x,y]); 求f 的⼆阶梯度k=0; 迭代初始值while 1 循环grad1z=subs(subs(grad1,x,x1(1)),y,x1(2)); 给f ⼀阶梯度赋初值 grad2z=subs(subs(grad2,x,x1(1)),y,x1(2)); 给f ⼆阶梯度赋初值x2=x1-inv(grad2z)*(grad1z)'; 核⼼迭代公式if norm(x1-x2)break;elseplot([x1(1),x2(1)],[x1(2),x2(2)],'-r*'); 画图k=k+1; 迭代继续x1=x2; 赋值endendend优点:在极⼩点附近收敛快缺点:但是要计算⽬标函数的hesse 矩阵最速下降法1. :选取初始点xo ,给定误差2. 计算⼀阶梯度。

若⼀阶梯度⼩于误差,停⽌迭代,输出3. 取()()()k k p f x =?4. 10t ()(), 1.min k k k k k k k k k k t f x t p f x tp x x t p k k +≥+=+=+=+进⾏⼀维搜索,求,使得令转第⼆步例题:求min (x-2)^4+(x-2*y)^2.初始值(0,3)误差为0.1(1)编写⼀个⽬标函数,存为f.mfunction z = f( x,y )z=(x-2.0)^4+(x-2.0*y)^2;end(2)分别关于x 和y 求出⼀阶梯度,分别存为fx.m 和fy.m function z = fx( x,y ) z=2.0*x-4.0*y+4.0*(x-2.0)^3;end和function z = fy( x,y )z=8.0*y-4.0*x;end(3)下⾯是脚本⽂件,⼀维搜索⽤的是黄⾦分割法Tic 计算时间eps=10^(-4);误差err=10;dt=0.01;x0=1.0;初始值y0=1.0;mm=0;while err>eps 黄⾦分割法dfx=-fx(x0,y0);dfy=-fy(x0,y0);tl=0;tr=1;确定⼀维搜索的区间h=3;nn=0;gerr=10;geps=10^(-4);while gerr>gepstll=tl+0.382*abs(tr-tl);trr=tl+0.618*abs(tr-tl);iff(x0+tll*h*dfx,y0+tll*h*dfy)>f(x0+trr*h*dfx,y0+trr*h*dfy) tl=tll;elsetr=trr;endgerr=abs(tl-tr); 区间的长度之差tt=0.5*(tl+tr);nn=nn+1;步数增加if nn>200 迭代终⽌条件breakendendx0=x0+tt*h*dfx; 重新迭代y0=y0+tt*h*dfy;err=sqrt(fx(x0,y0)^2+fy(x0,y0)^2);mm=mm+1;步数增加if mm>700 迭代步数超过700,终⽌breakendendres=[x0,y0];输出最后的x,y。

最速下降法(sd);共轭梯度法

最速下降法(sd);共轭梯度法

最速下降法(sd);共轭梯度法
最速下降法(SD)和共轭梯度法(CG)都是求解非线性优化问题中的常用算法。

最速下降法是基于梯度方向的一种搜索方法,在每一步所需找到函数在当前点的最陡方向,并沿着该方向走一步,直到达到要求的精度为止。

该方法速度快,收敛性好,但容易陷入“zigzag”现象,即由于步长过大或过小,导致序列在搜索方向上反复飞奔而收敛缓慢,同时,最速下降法对函数“弯曲性”敏感,函数梯度变化太快时收敛缓慢。

共轭梯度法是一种基于梯度方向的线性搜索方法,其优势在于快速收敛,准确性高。

其核心思想是,由于函数在一般性条件下不是QUADRATIC FUNCTION,因此,图像往往不是一个明显的"碗状",而是一个复杂的非线性图形。

在这种情况下,最速下降法很容易落入“zigzag”现象,收敛速度慢。

而共轭梯度法可以从不同方向进行极小值点的搜索,进而明显提高收敛速度。

总之,最速下降法适用于方向比较简单的情况,而共轭梯度法适用于方向较为复杂的情况。

根据不同的情况进行选择,可以有效地提高求解的效率和精度。

matlab共轭梯度法求解方程组

matlab共轭梯度法求解方程组

主题:matlab共轭梯度法求解方程组近年来,随着科学技术的不断发展,数学建模和计算机仿真成为科学研究和工程技术领域的重要手段。

在实际应用中,我们常常需要解决线性方程组的求解问题,而共轭梯度法作为一种高效的迭代求解方法,广泛应用于信号处理、图像处理、地球物理勘探和优化问题等领域。

本文将介绍如何利用matlab中的共轭梯度法求解线性方程组的基本原理和实际操作方法。

1. 共轭梯度法的基本原理共轭梯度法是一种迭代法,用于求解对称正定线性方程组Ax=b。

该方法的核心思想是通过一系列的迭代操作,逐步逼近方程组的解,直到满足一定的精度要求。

在每一步迭代中,共轭梯度法利用残差和方向向量的共轭性质,不断寻找最优的步长,从而实现方程组的求解。

2. matlab中共轭梯度法的基本调用方法在matlab中,调用共轭梯度法求解线性方程组非常简单。

需要将方程组的系数矩阵A和右端向量b输入到matlab中,然后利用内置函数conjugateGradient进行求解。

具体的调用方法如下:x = conjugateGradient(A, b, x0, maxIter, tol)其中,A为系数矩阵,b为右端向量,x0为初始解向量,maxIter为最大迭代次数,tol为精度要求。

调用完毕后,matlab将返回方程组的近似解x。

3. 共轭梯度法在实际工程中的应用共轭梯度法作为一种高效的求解方法,在工程技术领域得到了广泛的应用。

以图像处理为例,图像处理中经常需要解决大规模的线性方程组,而共轭梯度法能够高效地求解这类问题,提高了图像处理算法的效率和稳定性。

另外,在地球物理勘探中,共轭梯度法也被广泛应用于三维数据的快速处理和解释。

可以说,共轭梯度法在实际工程中发挥着重要的作用。

4. 共轭梯度法的优缺点分析尽管共轭梯度法具有非常高的效率和稳定性,但是该方法也存在一些缺点。

该方法只适用于对称正定的线性方程组,对于一般的线性方程组并不适用。

共轭梯度法的收敛速度受到方程条件数的影响,对于病态问题,可能收敛速度较慢。

MATLAB共轭梯度发

MATLAB共轭梯度发

最优化方法学院:专业班级:学生姓名:学号:共轭梯度法一.实验目的:(1).熟悉使用共轭梯度法求解无约束非线性规划问题的原理;(2).在掌握原理的基础上熟练运用此方法解决问题;(3).学会利用计算机语言编写程序来辅助解决数学问题;(4).解决问题的同时分析问题,力求达到理论与实践的相统一;(5).编写规范的实验报告.二.实验原理:<算法原理>:共轭梯度法为求解线性方程组而提出。

后来,人们把这种方法用于求解无约束最优化问题,使之成为一种重要的最优化方法。

共轭梯度法的基本思想是把共轭性与最速下降方法相结合,利用已知点处的梯度构造一组共轭方向,并沿这组方向进行搜索,求出目标函数的极小点。

根据共轭方向的基本性质,这种方法具有二次终止性。

在各种优化算法中,共轭梯度法是非常重要的一种。

其优点是所需存储量小,具有步收敛性,稳定性高,而且不需要任何外来参数。

共轭方向无约束最优化方法的核心问题是选择搜索方向.在本次实验中,我们运用基于共轭方向的一种算法—共轭梯度法三.算法流程图:四.实验结果:(1).实验函数f=(3*x1-cos(x2*x3)-1/2)^2+(x1^2-81*(x2+0.1)+sin(x3)+1.06)^2+(exp(-x1*x2)+20*x3+ 1/3*(10*3.14159-3))^2;(2).源程序:syms x1 x2 x3 r;f=(3*x1-cos(x2*x3)-1/2)^2+(x1^2-81*(x2+0.1)+sin(x3)+1.06)^2+(exp(-x1*x2)+20* x3+1/3*(10*3.14159-3))^2;x=[x1,x2,x3];df=jacobian(f,x);df=df.';error=0.000001;x0=[0,0,0]';g1=subs(df,x,x0);k=0;while(norm(g1)>error)if k==0d=-g1;elsebta=g1'*g1/(g0'*g0);d=-g1+bta*d0;endy=subs(f,x,x0+r*d);result=jintuifa(y,r);result2=gold(y,r,result);step=result2;x0=x0+step*d;g0=g1;g1=subs(df,x,x0);d0=d;k=k+1;end;kx0min=subs(f,[x1 x2 x3],x0)(3).实验结果k =26x0 =0.4996-0.0900-0.5259min =1.1496e-018(4).结果分析:由FR法得到的最优解为(0.4996, -0.0900, -0.5259),迭代的次数为32,由于精度的降低(ε= 61-),实际证明对结果也有一定的影响。

用MATLAB实现共轭梯度法求解实例

用MATLAB实现共轭梯度法求解实例

康福 0031一.无约束优化方法无约束优化方法的必要性一般机械优化设计问题,都是在一定的限制条件下追求某一指标为最小,它们都属于约束优化问题。

但是为什么要研究无约束优化问题? (1)有些实际问题,其数学模型本身就是一个无约束优化问题。

(2)通过熟悉它的解法可以为研究约束优化问题打下良好的基础。

(3)约束优化问题的求解可以通过一系列无约束优化方法来达到。

所以无约束优化问题的解法是优化设计方法的基本组成部分,也是优化方法的基础。

(4)对于多维无约束问题来说,古典极值理论中令一阶导数为零,但要求二阶可微,且要判断海赛矩阵为正定才能求得极小点,这种方法有理论意义,但无实用价值。

和一维问题一样,若多元函数F(X)不可微,亦无法求解。

但古典极值理论是无约束优化方法发展的基础。

共轭梯度法目前已研究出很多种无约束优化方法,它们的主要不同点在于构造搜索方向上的差别。

(1)间接法——要使用导数,如梯度法、(阻尼)牛顿法、变尺度法、共轭梯度法等。

(2)直接法——不使用导数信息,如坐标轮换法、鲍威尔法单纯形法等。

用直接法寻找极小点时,不必求函数的导数,只要计算目标函数值。

这类方法较适用于解决变量个数较少的(n ≤20)问题,一般情况下比间接法效率低。

间接法除要计算目标函数值外,还要计算目标函数的梯度,有的还要计算其海赛矩阵。

搜索方向的构成问题乃是无约束优化方法的关键。

1(0,1,2,)k k kk sk α+=+=L xx共轭梯度法是沿着共轭方向进行搜索,属于共轭方向法中的一种,该方法中每一个共轭向量都是依赖于迭代点处的负梯度而构造出来。

共轭梯度法作为一种实用的迭代法,它主要有下面的优点:(1)算法中,系数矩阵A的作用仅仅是用来由已知向量P 产生向量W=AP ,这不仅可充分利用A的稀疏性,而且对某些提供矩阵A较为困难而由已知向量P 产生向量W=AP 又十分方便的应用问题是很有益的。

(2)不需要预先估计任何参数就可以计算,这一点不像SOR 等; (3)每次迭代所需的计算,主要是向量之间的运算,便于并行化。

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

西京学院数学软件实验任务书
课程名称数学软件实验班级数0901
学号0912020107姓名李亚强
实验课题线性方程组的最速下降法与共轭梯度法
实验目的熟悉线性方程组的最速下降法与共轭梯度法
实验要求
运用Matlab/C/C++/Java/Maple/Mathematica等其中
一种语言完成
实验内容线性方程组的最速下降法线性方程组的共轭梯度法
成绩教师
实验五实验报告
一、实验名称:最速下降法与共轭梯度法解线性方程组。

二、实验目的:进一步熟悉理解掌握最速下降法与共轭梯度法解法思路,提高matlab 编程能力。

三、实验要求:已知线性方程矩阵,应用最速下降与共轭梯度法在相关软件编程求解线性方程组的解。

四、实验原理:
1.最速下降法:
从某个初始点)0(X 出发,沿)(X f 在点)0(X 处的负梯度方向
)
0()0()0()(AX b X f r -=-∇=求得)(X f 的极小值点)1(X ,即
)(min )0()0(0
r X f λλ+>然后从)1(X 出发,重复上面的过程得到)2(X 。

如此下去,
得到序列{)(k X }
)
(...)()()()1()0(k X f X f X f >>>可以证明,从任一初始点)0(X 出发,用最速下降法所得
到的序列{)(k X }均收敛于问题使X 最小化)(X f 的解,也就是方程组b AX =的解。

其收敛速度取决于1
1λλλλ+-n n ,其中1λ,n λ分别
为A 的最小,最大特征值。

最速下降法迭代格式:给定初值)0(X ,)(k X 按如下方法决定:
())
()(1)(k )()()()(k )
()(X ,,)(k k k k T k k T k k k k r X Ar r r r AX b X f r λλ+=>
<><=-=-∇=+2.共轭梯度法
其基本步骤是在点)(k X 处选取搜索方向)(k d ,使其与前一次的搜索方向)1(-k d 关于A 共轭,即
(1)()(1),0
k k k d d Ad --<>=然后从点)(k X 出发,沿方向)(k d 求得)(X f 的极小值点)1(+k X ,即
)
(min )()()(0)1(k d X f X f k k λλ+=>+如此下去,得到序列{)(k X }。

不难求得0,)1()(>=<-k k Ad d 的解为
)()1()1()()()
()1(,,k k k k k k k d Ad d d AX b X X ><>-<+=--+注意到)(k d 的选取不唯一,我们可取
)
1(1)()()(--+-∇=k k k k d X f d β由共轭的定义0,)1()(>=<-k k Ad d 可得:
>
<><-=----)1()1()1()(1,,k k k k k Ad d Ad r β共轭梯度法的计算过程如下:第一步:取初始向量)0(X ,计算
⎪⎪⎩
⎪⎪⎨⎧+=><><-=-=-∇==(0)0(0)(1))0()0()0()0(0(0)
(0)(0)(0)d X X ,,X )X (r d λλAd d Ad r A b f 第1+k 步:计算
⎪⎪⎪⎪⎪⎩
⎪⎪⎪⎪⎪⎨⎧+=><><-=+=><><-=-=-∇=+------(k)0(k)1)(k )()()()()
1(1(k))()1()1()1()(1(k)
(k)(k)d X X ,,r ,,X )X (r λλββk k k k k k k k k k k k k Ad d Ad r d d Ad d Ad
r A b f 五、实验内容:
%最速下降法
function [x,k]=fastest(A,b,eps);
x0=zeros(size(b),1);
x=x0;
k=0;
m=1000;
tol=1;
while tol>=eps
r=b-A*x0;
q=dot(r,r)/dot(A*r,r);
x=x0+q*r;
k=k+1;
tol=norm(x-x0);
x0=x;
if k>=m
disp('迭代次数太多,可能不收敛!');
return;
end
end
x
k
%共轭梯度法
function[k,x]=gong_e(A,b)
esp=input('请输入允许误差esp=');
x0=input('请输入初始值x0=');
k=0;
r0=b-A*x0;%求出dangqian梯度while norm(r0)>esp
r0=b-A*x0;
k=k+1;
if k==1
p0=r0;
else
lamda=(r0'*r0)/(p0'*A*p0);
r1=r0-lamda*A*p0;
p0=r0+(r0'*r0)/(r1'*r1)*p0;
x1=x0+lamda*p0;
x0=x1;
r0=r1;
end
end
x=r0;
k;
end
六、实验结果:
A=[520;641;125];
b=[1018-14]';
eps=1.0e-6;
x=
-0.8750
7.1875
-5.5000
k=
60。

相关文档
最新文档