文档库 最新最全的文档下载
当前位置:文档库 › Matlab积分

Matlab积分

Matlab积分
Matlab积分

MATLAB数值积分

1 定积分的MATLAB符号计算 (2)

2 变限积分的MATLAB符号计算 (2)

3 矩形公式的MATLAB程序 (2)

4 梯形公式的MATLAB程序 (3)

5 辛普森(SIMPSON)数值积分的MATLAB程序 (5)

6 牛顿-科茨(NEWTON-COTES)数值积分和误差分析的MATLAB程序 (6)

7 利用三次样条求表格型数值积分的MATLAB方法 (11)

8 利用拉格朗日插值等方法求表格型数值积分的MATLAB方法 (12)

9 龙贝格(ROMBERG)公式及其MATLAB程序 (15)

10 自适应积分及其MATLAB程序 (17)

11 高斯(GAUSS)型积分公式及其MATLAB程序 (17)

11.1在[-1,1]上的高斯-勒让德积分及其MATLAB程序 (17)

11.2在[A,B]上的高斯-勒让德积分公式及MATLAB程序 (20)

12 拉道积分公式和洛巴托积分公式及其MATLAB程序 (21)

13 无穷积分的符号计算及其MATLAB程序 (28)

14 无穷积分的近似计算及其MATLAB程序 (30)

14.1累积求和无穷积分法及其MATLAB程序 (30)

14.2无穷区间的截断法及其MATLAB程序 (31)

14.3高斯-拉盖尔求积公式及其MATLAB程序 (33)

14.3高斯—埃尔米特求积公式及其MATLAB程序 (38)

15 无界函数反常积分的符号计算及其MATLAB程序 (41)

16 无界函数反常积分的近似计算及其MATLAB程序 (42)

16.1无界函数的反常积分的“挖去”法及其MATLAB程序 (42)

16.2高斯—切比雪夫求积公式及其MATLAB程序 (43)

17 多重积分的计算及其MATLAB程序 (45)

17.1二重积分的符号计算及其MATLAB程序 (45)

17.2二重积分的梯形公式及其MATLAB程序 (46)

17.3矩形域上的辛普森公式及其MATLAB程序 (49)

17.4一般域上二重积分的数值计算及其MATLAB程序 (51)

17.5三重积分的计算及其MATLAB程序 (51)

1 定积分的MATLAB 符号计算

例 由x y sin =,2

3,21,cos =-==x x x y 所围成的平面区域D .求平面区域D 的面积S .

解 输入作函数图形的程序

>> x=-1:0.001:2; F1= sin(x); F2=cos(x);

plot(x ,F1,'b-',x ,F2,'g-'),

axis([-1,pi/4+1,-1.3,1.3]),

xlabel('x'), ylabel('y'),

title('y=sinx , y=cosx 和x=-0.5及x=1.5所围成的平面区

域的图形')

运行后屏幕显示图形.

求平面区域D 的面积S .输入计算面积S 的程序

>> syms x

f1= cos(x)-sin(x); f2=-f1; S1 =int(f1,x,-0.5,pi/4);

S2=int(f2,x, pi/4,1.5); S=S1+S2,Sj= double (S)

运行后屏幕显示计算面积的值 S 及其近似值 Sj 如下

S =

2*2^(1/2)+sin(1/2)-cos(1/2)-sin(3/2)-cos(3/2)

Sj =

1.36203791318826

2 变限积分的MATLAB 符号计算

例 已知?=2

)(x x x F e )2sin(3t t +d t ,求)('x F .

解 输入程序:

>> syms x t

F1=int(exp(t)*sin(2+sqrt(t^3)),x,0);

F2=int(exp(t)*sin(2+sqrt(t^3)),0,x^2);

Fi= F1+ F2;

dF=diff(Fi)

运行后屏幕显示计算变限积分)(x F 的导数.

3 矩形公式的MATLAB 程序

(一) 函数sum 的调用格式

调用格式一:sum(X)

调用格式二:sum (X,DIM)

例 用MATLAB 和矩形公式计算?

π20e x sin d x ,并与精确值比较.

解 将]2/,0[π分成20等份,步长为40/π,输入程序

>> h=pi/40; x=0:h:pi/2; y=exp(sin(x));

z1=sum(y(1:20))*h, z2=sum(y(2:21))*h,

运行后屏幕显示矩形公式计算结果分别如下

z1 = z2 =

3.0364 3.1713

求定积分的精确值,输入程序

>> syms x

F=int(exp(sin(x)),x,0, pi/2), Fs= double (F),

wz1=abs( Fs-z1), wz2= abs( Fs-z2)

运行后屏幕显示定积分的精确值Fs 和与用矩形公式计算结果的绝对误差wz 1、wz 2.

(二) 函数cumsum 的调用格式

调用格式一:cumsum(X)

调用格式二:cumsum (X,DIM)

例 用MATLAB 的函数sum 和 cumsum 及矩形公式计算?π

2

0e x x sin -d x ,并与精确值比较.

解 将]2/,0[π分成20等份,步长为40/π,输入程序如下(注意sum 和 cumsum 的用法)

>> h=pi/40; x=0:h:pi/2; y=exp(-x).*sin(x);

z1=sum(y(1:20))*h,z2=sum(y(2:21))*h,

z=cumsum(y); z11=z(20)*h, z12=(z(21)-z(1))*h,

运行后屏幕显示计算结果分别如下

z1 = z2 = z11 = z12 =

0.3873 0.4036 0.3873

0.4036

求定积分的精确值,输入程序

>> syms x

F=int(exp(-x)*sin(x),x,0, pi/2)

Fs= double (F) ,wz1=abs( Fs-z1), wz2= abs( Fs-z2)

运行后屏幕显示定积分的精确值Fs 和用矩形公式计算结果的绝对误差wz 1、wz 2分别如下

F = Fs =

1/2*(-1+exp(pi)^(1/2))/exp(pi)^(1/2) 0.3961

wz1 = wz2 =

0.0088

0.0075

4 梯形公式的MATLAB 程序

(一) 根据梯形公式和估计误差公式自己编写MATLAB 程序计算定积分

例 分别取80/,800/,0008/πππ=h ,用梯形公式计算定积分

?π=2

/0I e x sin d x ,并与精确值比较.然后观察h 对计算结果的有效数字和绝对

误差的影响.

解 编写并输入如下程序

>>h=pi/8000;a=0;b=pi/2;x=a:h:b;n=length(x),

y=exp(sin(x));

z1=(y(1)+y(n))*h/2; z2=sum(y(2:n-1))*h;

z8000=z1+z2,

syms t

f=exp(sin(t)); intf=int(f,t,a,b), Fs=double(intf), Juewucha8000=abs(z8000-Fs)

运行后屏幕显示取8000/π=h 时,积分区间]2/,0[π上等距节点的个数n ,用梯形公式计算定积分I 的值z 8000和精确值intf 的近似值Fs 及其绝对误差Juewucha 8000.

(二) 用函数trapz 计算定积分

调用格式一:Z =trapz(Y)

调用格式二:Z =trapz(X,Y)

调用格式三:Z = trapz (X,Y,DIM) 或 trapz (Y,DIM)

(三) 用函数cumtrapz 计算定积分

调用格式一:Z =cumtrapz (Y)

调用格式二:Z =cumtrapz (X,Y)

调用格式三:Z = cumtrapz (X,Y,DIM) 或 cumtrapz (Y,DIM)

例 用MATLAB 的函数trapz 和cumtrapz 分别计算?π2

/0e x x sin -d x ,精

确到410-,并与矩形公式比较.

解 将]2/,0[π分成20等份,步长为40/π,输入程序如下(注意trapz(y)是单位步长, trapz(y)*h=trapz(x,y)):

>> h=pi/40; x=0:h:pi/2; y=exp(-x).*sin(x);

z1=sum(y(1:20))*h, z2=sum(y(2:21))*h, z=(z1+z2)/2 z3=trapz(y)*h, z3h=trapz(x,y), z3c=cumtrapz(y)*h, 运行后屏幕显示用矩形公式(9.3),(9.4)计算结果z 1、z 2和二者的平均数z 、函数trapz 和cumtrapz 分别计算结果z 3、z 3c .

(四)梯形数值积分的MATLAB 主程序

梯形数值积分的MATLAB 主程序

function T=rctrap(fun,a,b,m)

n=1;h=b-a; T=zeros(1,m+1); x=a;

T(1)=h*(feval(fun,a)+feval(fun,b))/2;

for i=1:m

h=h/2; n=2*n; s=0;

for k=1:n/2

x=a+h*(2*k-1); s=s+feval(fun,x);

end

T(i+1)=T(i)/2+h*s;

end

T=T(1:m);

例 用rctrap 计算?ππ=2

021I e 22x -d x ,递归14次,并将计算结果与精

确值比较.

解 输入程序

>>T=rctrap(@fun,0,pi/2,14), syms t

fi=int(exp((-t^2)/2)/(sqrt(2*pi)),t,0, pi/2);

Fs= double(fi), wT= double(abs(fi-T))

运行后屏幕显示I 精确值F s ,用rctrap 计算I 的递归值T 和T 与精确值F s 的绝对误差w T

5 辛普森(Simpson )数值积分的MATLAB 程序

调用格式一:quad(‘fun’,a,b)

调用格式二:quad(‘fun’,a,b,tol)

调用格式三:[Q,FCNT] = quad (...)

调用格式四:quad(‘fun’,a,b, tol,TRACE)

调用格式五:quad(‘fun’,a,b, tol,TRACE,P1,P2, …)

复合辛普森(Simpson )数值积分的MATLAB 主程序

function y=comsimpson(fun,a,b,n)

z1=feval (fun,a)+ feval (fun,b);m=n/2;

h=(b-a)/(2*m); x=a;

z2=0; z3=0; x2=0; x3=0;

for k=2:2:2*m

x2=x+k*h; z2= z2+2*feval (fun,x2);

end

for k=3:2:2*m

x3=x+k*h; z3= z3+4*feval (fun,x3);

end

y=(z1+z2+z3)*h/3;

例 用comsimpson.m 和quad.m 分别计算定积分?π=1021I e 22x -d x ,

取精度为410-,并与精确值比较.

解 输入程序

>> [Q1,FCNT14] = quad(@fun,0,1,1.e-4,3),

Q2 =comsimpson (@fun,0,1,10000)

syms x

fi=int(exp( (-x.^2)./2)./(sqrt(2*pi)),x,0, 1); Fs= double (fi)

wQ1= double (abs(fi-Q1) ), wQ2= double (abs(fi-Q2) ) 运行后屏幕显示I 的精确值F s ,用comsimpson.m 和quad.m 分别计算I 的近似值Q 2、Q 1和迭代次数FCNT14,取精度分别为410-,Q 2、Q 1分别与精确值F s 的绝对误差wQ 2, wQ 1如下

9 0.0000000000 2.71580000e-001

0.1070275100

11 0.2715800000 4.56840000e-001 0.1597942242

13 0.7284200000 2.71580000e-001

jacobi G-S,超松弛迭代法MATLAB程序

function iteration A=[10,1,2,3,4; 1,9,-1,2,-3; 2,-1,7,3,-5; 3,2,3,12,-1; 4,-3,-5,-1,15]; b=[12,-27,14,-17,12]'; x0=[0,0,0,0,0]'; tol=1e-12; disp('jacobi迭代法的结果和次数如下:') [x,k]=Fjacobi(A,b,x0,tol) disp('G-S迭代法的结果和次数如下:':') [x,k]=Fgseid(A,b,x0,tol) disp('超松弛的结果和次数如下:':') [x,k]=Fsor(A,b,x0,1.2,tol) disp('共轭梯度法的结果和次数如下:':') [x,k]=Fcg(A,b,x0,tol) %jacobi迭代法 function [x,k]=Fjacobi(A,b,x0,tol) max=300; D=diag(diag(A)); L=-tril(A,-1); U=-triu(A,1); B=D\(L+U); f=D\b; x=B*x0+f; k=1; while norm(x-x0)>=tol x0=x; x=B*x0+f; k=k+1; if(k>=max) disp('μü′ú3?1y300′?£?·?3ì×é?é?ü2?ê?á2'); return; end end %G-S迭代法 function [x,k]=Fgseid(A,b,x0,tol) max=300; D=diag(diag(A)); L=-tril(A,-1); U=-triu(A,1); G=(D-L)\U; f=(D-L)\b; x=G*x0+f; k=1; while norm(x-x0)>=tol x0=x; x=G*x0+f; k=k+1; if(k>=max) disp('μü′ú3?1y300′?£?·?3ì×é?é?ü2?ê?á2'); return; end

MATLAB代码 解线性方程组的迭代法

解线性方程组的迭代法 1.rs里查森迭代法求线性方程组Ax=b的解 function[x,n]=rs(A,b,x0,eps,M) if(nargin==3) eps=1.0e-6;%eps表示迭代精度 M=10000;%M表示迭代步数的限制值elseif(nargin==4) M=10000; end I=eye(size(A)); n=0; x=x0; tol=1; %迭代过程 while(tol>eps) x=(I-A)*x0+b; n=n+1;%n为最终求出解时的迭代步数tol=norm(x-x0); x0=x; if(n>=M) disp('Warning:迭代次数太多,可能不收敛!'); return; end end 2.crs里查森参数迭代法求线性方程组Ax=b的解 function[x,n]=crs(A,b,x0,w,eps,M) if(nargin==4) eps=1.0e-6;%eps表示迭代精度 M=10000;%M表示迭代步数的限制值 elseif(nargin==5) M=10000; end I=eye(size(A)); n=0; x=x0; tol=1; %迭代过程 while(tol>eps) x=(I-w*A)*x0+w*b; n=n+1;%n为最终求出解时的迭代步数tol=norm(x-x0); x0=x;

if(n>=M) disp('Warning:迭代次数太多,可能不收敛!'); return; end end 3.grs里查森迭代法求线性方程组Ax=b的解 function[x,n]=grs(A,b,x0,W,eps,M) if(nargin==4) eps=1.0e-6;%eps表示迭代精度 M=10000;%M表示迭代步数的限制值 elseif(nargin==5) M=10000; end I=eye(size(A)); n=0; x=x0; tol=1;%前后两次迭代结果误差 %迭代过程 while(tol>eps) x=(I-W*A)*x0+W*b;%迭代公式 n=n+1;%n为最终求出解时的迭代步数tol=norm(x-x0); x0=x; if(n>=M) disp('Warning:迭代次数太多,可能不收敛!'); return; end end 4.jacobi雅可比迭代法求线性方程组Ax=b的解 function[x,n]=jacobi(A,b,x0,eps,varargin) if nargin==3 eps=1.0e-6; M=200; elseif nargin<3 error return elseif nargin==5 M=varargin{1}; end D=diag(diag(A));%求A的对角矩阵 L=-tril(A,-1);%求A的下三角阵

MATLAB样例之雅克比迭代法

要求: 下面分别使用雅克比迭代法和高斯-赛德尔迭代法求一个方程组的近似解用的线性方程组是按实验要求给的: 7*x1+x2+2*x3=10 x1+8*x2+2*x3=8 2*x1+2*x2+9*x3=6 雅克比迭代法的matlab代码:(老师写的) A=[7,1,2;1,8,2;2,2,9]; b=[10;8;6]; if(any(diag(A))==0) error('error,pause') end eps=input('误差限eps='); N=input('迭代次数N='); D=diag(diag(A)); B=inv(D)*(D-A); f=inv(D)*b; K=0; x0=zeros(size(b)); while 1 x1=B*x0+f K=K+1; fprintf('第-次迭代的近似解为',K) disp(x1'); if norm(x1-x0,inf)N fprintf('迭代超限') end x0=x1; end 高斯-赛德尔迭代法matlab代码:(自己改的)

A=[7,1,2;1,8,2;2,2,9]; b=[10;8;6]; if(all(diag(A))==0) error('error,pause') end eps=input('误差限eps='); N=input('迭代次数N='); D=diag(diag(A)); B=inv(D)*(D-A); f=inv(D)*b; K=0; x0=zeros(size(b)); x00=x0; while 1 x11=B*x0+f; x00(1,1)=x11(1,1); x12=B*x00+f; x00(2,1)=x12(2,1); x13=B*x00+f; x00(3,1)=x13(3,1); x1=x00 K=K+1; fprintf('第-次迭代的近似解为',K) disp(x1'); if norm(x1-x0,inf)N fprintf('迭代超限') end x0=x1; end

MATLAB实现迭代法最佳松弛因子的选取

迭代法最佳松弛因子的选取 一、问题提出: 针对矩阵430341014A ?? ??=-?? ??-?? ,b=[24;30;-24],用SOR 迭代求解。并选出最佳松弛 因子。理论分析 1.24ω==≈。做出()L ωρ关于ω函数 的图像。 二、理论基础 选取分裂矩阵M 为带参数的下三角矩阵)(1 wL D w M -=, 其中w>0为可选择的松弛因子. 于是,由 ?????+=+f Bx x x k k ) ()1()0() (初始向量 (k=0,1,…,)可构造一个迭代法,其迭代矩阵为A wL D w I L w 1)(---≡ =).)1(()(1wU D w wL D +--- 从而得到解Ax=b 的主次逐次超松弛迭代法. 解Ax=b 的SOR 方法为 ?????+=+f Bx x x k k ) ()1()0() (初始向量 (k=0,1,…,) (1) 其中 w L =).)1(()(1wU D w wL D +---(2) b wL D w f 1)(--= 下面给出解Ax=b 的SOR 迭代法的分量计算公式.记 ,),...,,...,() () () (1)(T k n k i k k x x x x = 由(1)式可得 ,))1(()()()1(wb x wU D w x wL D k k ++-==-+ ).()()()1()()1(k k k k k Dx Ux Lx b w Dx Dx -+++=++ (3) 由此,得到解Ax=b 的SOR 方法的计算公式

?????????==--+==∑∑-==++.),1,0;,...,2,1(/)(,),...,(11) (1)()1()0()0(1)0(为松弛因子 w k n i a x a x a b w x x x x x ii i j n i j k j ij k j ij i k i k i T n (4) 或 ?? ?? ? ??????==--=??+==∑∑-==++.,...),1,0;,...,2,1()/(,,),...,(.11)()1() () 1()0()0(1)0(为松弛因子w k n i a x a x a b w x x x x x x x i j n i j ii k j ij k j ij i i i k i k i T n (5) ※ 若要求选取出最佳松弛因子,则有两种方法: ⑴、 给出w 的最佳范围,当取不同的w 值时,会求出不同的谱半径R 的值, 然后判断出值最小的谱半径。那么这个最小的谱半径所对应的w ,即为所求最佳松弛因子。 ⑵、 给出w 的最佳范围,当取不同的w 值时,由(2)式进行迭代,看它们在 相同精度范围内的迭代次数,找出迭代次数最低的那一个,其所应用的w 即为最佳松弛因子。 三、实验内容: 从表格中可以看出,迭代次数随着松弛因子的增长而呈现先减后增的趋势,当谱半径最小时,其迭代次数最小。则表示出谱半径最小时,其松弛因子为最佳松弛因子。

matlab迭代法代码

matlab 迭代法代码 1、%用不动点迭代法求方程x-e A x+4=0的正根与负根,误差限是 10A-6% disp(' 不动点迭代法 '); n0=100; p0=-5; for i=1:n0 p=exp(p0)-4; if abs(p-p0)<=10(6) if p<0 disp('|p-p0|=') disp(abs(p-p0)) disp(' 不动点迭代法求得方程的负根为 :') disp(p); break; else disp(' 不动点迭代法无法求出方程的负根 .') end else p0=p; end end

if i==n0 disp(n0) disp(' 次不动点迭代后无法求出方程的负根') end p1=1.7; for i=1:n0 pp=exp(p1)-4; if abs(pp-p1)<=10(6) if pp>0 disp('|p-p1|=') disp(abs(pp-p1)) disp(' 用不动点迭代法求得方程的正根为 ') disp(pp); else disp(' 用不动点迭代法无法求出方程的正根 '); end break; else p1=pp; end end if i==n0

disp(n0) disp(' 次不动点迭代后无法求出方程的正根 ') end 2、%用牛顿法求方程x-e A x+4=0的正根与负根,误差限是disp(' 牛顿法') n0=80; p0=1; for i=1:n0 p=p0-(p0-exp(p0)+4)/(1-exp(p0)); if abs(p-p0)<=10(6) disp('|p-p0|=') disp(abs(p-p0)) disp(' 用牛顿法求得方程的正根为 ') disp(p); break; else p0=p; end end if i==n0 disp(n0) disp(' 次牛顿迭代后无法求出方程的解 p1=-3; for i=1:n0 p=p1-(p1-exp(p1)+4)/(1-exp(p1)); 10A-6 ') end

【良心出品】不动点迭代法matlab程序

实验四 姓名:木拉丁。尼则木丁班级:信计08-2 学号:20080803405 实验地点:新大机房 实验目的:通过本实验学习利用MATLAB不动点迭代法,抛物线法,斯特芬森迭代法解非线性方程组,及其编程实现,培养编程与上机调试能力。 实验要求:①上机前充分准备,复习有关内容,写出计算步骤,查对程序; ②完成实验后写出完整的实验报告,内容应该包括:所用的算法语言, 算法步骤陈述,变量说明,程序清单,输出计算结果,结果分析等等; ③用编好的程序在Matlab环境中执行。 迭代法 MATLAB程序: function pwxff(f,x0,x1,x2,d,n) f=inline(f); x(1)=x0; x(2)=x1; x(3)=x2; w1=(f(x(2))-f(x(3)))/(x(2)-x(3)); t1=(f(x(1))-f(x(3)))/(x(1)-x(3)); t2=(f(x(1))-f(x(2)))/(x(1)-x(2)); w2=1/(x(1)-x(2))*(t1-t2); w=w1+w2*(x(3)-x(2));

for k=3:n x(k+1)=x(k)-2*f(x(k))/(w+sqrt(w^2-4*f(x(k))*w2)); if abs(x(k+1)-x(k))

线性方程组的迭代解法(Matlab)

第六章线性方程组的迭代解法 2015年12月27日17:12 迭代法是目前求解大规模稀疏线性方程组的主要方法之一。包括定常迭代法和不定常迭代法,定常迭代法的迭代矩阵通常保持不变,包括有雅可比迭代法(Jacobi)、高斯-塞德尔迭代法(Gauss-Seidel)、超松弛迭代法(SOR) 1.雅可比迭代法(Jacobi) A表示线性方程组的系数矩阵,D表示A的主对角部分,L表示下三角部分,U表示上三角部分。 A=D+L+U 要解的方程变为Dx+Lx+Ux=b x=D^(-1)(b-(L+U)x) 所以Jocabi方法如下: Matlab程序 function [x,iter] =jacobi(A,b,tol) D=diag(diag(A)); L=D-tril(A); U=D-triu(A); x=zeros(size(b)); for iter=1:500 x=D\(b+L*x+U*x); error=norm(b-A*x)/norm(b); if(error

二分法、简单迭代法的matlab代码实现教学文案

实验一非线性方程的数值解法(一)信息与计算科学金融崔振威201002034031 一、实验目的: 熟悉二分法和简单迭代法的算法实现。 二、实验内容: 教材P40 2.1.5 三、实验要求 1 根据实验内容编写二分法和简单迭代法的算法实现 2 简单比较分析两种算法的误差 3 试构造不同的迭代格式,分析比较其收敛性 (一)、二分法程序: function ef=bisect(fx,xa,xb,n,delta) % fx是由方程转化的关于x的函数,有fx=0。 % xa 解区间上限 % xb 解区间下限 % n 最多循环步数,防止死循环。 %delta 为允许误差 x=xa;fa=eval(fx); x=xb;fb=eval(fx); disp(' [ n xa xb xc fc ]'); for i=1:n xc=(xa+xb)/2;x=xc;fc=eval(fx); X=[i,xa,xb,xc,fc]; disp(X), if fc*fa<0 xb=xc; else xa=xc; end if (xb-xa)

k=0; while abs(x-x0)>eps & k> fplot('[x^5-3*x^3-2*x^2+2]',[-3,3]);grid 得下图: 由上图可得知:方程在[-3,3]区间有根。 (2)、二分法输出结果 >> f='x^5-3*x^3-2*x^2+2' f = x^5-3*x^3-2*x^2+2 >> bisect(f,-3,3,20,10^(-12)) 2.0000 - 3.0000 0 -1.5000 0.0313

二分法、简单迭代法的matlab代码实现教学文案

二分法、简单迭代法的m a t l a b代码实现

实验一非线性方程的数值解法(一)信息与计算科学金融崔振威 201002034031 一、实验目的: 熟悉二分法和简单迭代法的算法实现。 二、实验内容: 教材P40 2.1.5 三、实验要求 1 根据实验内容编写二分法和简单迭代法的算法实现 2 简单比较分析两种算法的误差 3 试构造不同的迭代格式,分析比较其收敛性 (一)、二分法程序: function ef=bisect(fx,xa,xb,n,delta) % fx是由方程转化的关于x的函数,有fx=0。 % xa 解区间上限 % xb 解区间下限 % n 最多循环步数,防止死循环。 %delta 为允许误差 x=xa;fa=eval(fx); x=xb;fb=eval(fx); disp(' [ n xa xb xc fc ]'); for i=1:n xc=(xa+xb)/2;x=xc;fc=eval(fx); X=[i,xa,xb,xc,fc]; disp(X), if fc*fa<0 xb=xc; else xa=xc;

end if (xb-xa)eps & k> fplot('[x^5-3*x^3-2*x^2+2]',[-3,3]);grid 得下图:

迭代解法的matlab实现

解线性方程组b AX =的迭代法是从初始解出发,根据设计好的步骤用逐次求出的近似解逼近精确解.在第三章中介绍的解线性方程组的直接方法一般适合于A 为低阶稠密矩阵(指n 不大且元多为非零)的情况,而在工程技术和科学计算中常会遇到大型稀疏矩阵(指n 很大且零元较多)的方程组,迭代法在计算和存贮两方面都适合后一种情况.由于迭代法是通过逐次迭代来逼近方程组的解,所以收敛性和收敛速度是构造迭代法时应该注意的问题.另外,因为不同的系数矩阵具有不同的性态,所以大多数迭代方法都具有一定的适用范围.有时,某种方法对于一类方程组迭代收敛,而对另一类方程组迭代时就发散.因此,我们应该学会针对具有不同性质的线性方程组构造不同的迭代. 4.1 迭代法和敛散性及其MATLAB 程序 4.1.2 迭代法敛散性的判别及其MATLAB 程序 根据定理4.1和谱半径定义,现提供一个名为pddpb.m 的M 文件,用于判别迭代公H=eig(B);mH=norm(H,inf); if mH>=1 disp('请注意:因为谱半径不小于1,所以迭代序列发散,谱半径mH 和B 的所 有的特征值H 如下:') else disp('请注意:因为谱半径小于1,所以迭代序列收敛,谱半径mH 和B 的所有 的特征值H 如下:') end mH 4.1.3 与迭代法有关的MATLAB 命令 (一) 提取(产生)对角矩阵和特征值 可以用表4–1的MATLAB 命令提取对角矩阵和特征值. (二) 提取(产生)上(下)三角形矩阵

可以用表4–2的MATLAB命令提取矩阵的上三角形矩阵和下三角形矩阵. (三)稀疏矩阵的处理 对稀疏矩阵在存贮和运算上的特殊处理,是MA TLAB进行大规模科学计算时的特点和优势之一.可以用表4–3的MATLAB命令,输入稀疏矩阵的非零元(零元不必输入),即可进行运算. 4.2 雅可比(Jacobi)迭代及其MATLAB程序 4.2.2 雅可比迭代的收敛性及其MATLAB程序 [n m]=size(A); for j=1:m a(j)=sum(abs(A(:,j)))-2*(abs(A(j,j))); end for i=1:n if a(i)>=0 disp('请注意:系数矩阵A不是严格对角占优的,此雅可比迭代不一定收敛') return end end if a(i)<0 disp('请注意:系数矩阵A是严格对角占优的,此方程组有唯一解,且雅可比迭代收敛') end 例4.2.2 用判别雅可比迭代收敛性的MATLAB主程序,判别由下列方程组的雅可比迭

matlab Jacobi迭代法Gauss-seidel和SOR迭代

1.Jacobi迭代法 例1 用jacobi迭代法求解代数线性代数方程组,保留四位有效数字(err=1e-4) 其中A=[8 -1 1;2 10 -1;1 1 -5];b=[1 ;4; 3]。 解:编写jacobi迭代法的函数文件,保存为jacobi.m function [x,k]=jacobi(A,b,x0,eps,N) % 求解Ax=b;x0为初始列向量;eps为误差容限;N为最大迭代次数 % 输出x为近似解;k为迭代次数 n=length(A); x=zeros(n,1); for k=1:N for i=1:n ――――――― end if norm(x-x0,inf)

end x0=x; end 编写主程序如下 format long clear A=[8 -1 1;2 10 -1;1 1 -5]; b=[1 ;4; 3]; x0=[0.125; 0.4 ;-0.6 ]; % x0为初始列向量N为最大迭代次数err=1e-4; % err为误差容限 N=25; % N为最大迭代次数 [x,k]=jacobi(A,b,x0,err,N) 得到结果如下 x = 0.22492315625000 0.30561995000000 -0.49388680000000

k = 6 2.Gauss-seidel迭代法 例2 用Gauss-seidel迭代法求解代数线性代数方程组,保留四位有效数字(err=1e-4) 其中A=[8 -1 1;2 10 -1;1 1 -5];b=[1 ;4; 3]。 解:编写Gauss-seidel迭代法的函数文件,保存为gaus.m function [x,k]=gaus(A,b,x0,eps,N) % 求解Ax=b;x0为初始列向量;eps为误差容限;N为最大迭代次数% 输出x为近似解;k为迭代次数 n=length(A); x=zeros(n,1); for k=1:N for i=1:n ―――――― end if norm(x-x0,inf)

线性方程组迭代法matlab

第三章 大作业 考虑线性方程组 12345615310002213131002201310010013101130013122150001322x x x x x x ????- ? ? ? ??? ? ?-- ? ? ? ? ? ? ?-- ? ?= ? ? ?-- ? ? ? ? ? ?-- ? ? ? ??? ? ? ? ?-???? (1)分别用Jacobi 方法和G-S 方法求解上述方程组,并比较它们的收敛快慢; (2)用SOR 方法求解,试分析松弛因子ω的选取对方法收敛的影响,并找到最佳松弛因子 J 法: function [X]=jf(A,B,X) m=rank(A); X=reshape(X,length(X),1); B=reshape(B,length(B),1); k=1; % if m~=size(A,1) %求行数 disp('A 不可逆'); return end % Bj=zeros(m,m); % while k

% k=1; while k<=m Bj(k,1:k-1)=-A(k,1:k-1)/A(k,k); Bj(k,1+k:m)=-A(k,k+1:m)/A(k,k); f(k)=B(k)/A(k,k); k=k+1; end f=f'; % k=0; t=norm((X-f),inf); while t>0.00000001 %回代求解 X1=Bj*X+f; t=norm((X1-X),inf); X=X1; k=k+1; end; X=X'; disp('迭代次数') k

最新matlab迭代法代码

1、%用不动点迭代法求方程 x-e^x+4=0的正根与负根,误差限是10^-6% disp('不动点迭代法'); n0=100; p0=-5; for i=1:n0 p=exp(p0)-4; if abs(p-p0)<=10^(-6) if p<0 disp('|p-p0|=') disp(abs(p-p0)) disp('不动点迭代法求得方程的负根为:') disp(p); break; else disp('不动点迭代法无法求出方程的负根.') end else p0=p; end end if i==n0 disp(n0) disp('次不动点迭代后无法求出方程的负根') end p1=1.7; for i=1:n0 pp=exp(p1)-4; if abs(pp-p1)<=10^(-6) if pp>0 disp('|p-p1|=') disp(abs(pp-p1)) disp('用不动点迭代法求得方程的正根为') disp(pp); else disp('用不动点迭代法无法求出方程的正根'); end break; else p1=pp; end end if i==n0 disp(n0)

disp('次不动点迭代后无法求出方程的正根') end 2、%用牛顿法求方程 x-e^x+4=0的正根与负根,误差限是10^-6 disp('牛顿法') n0=80; p0=1; for i=1:n0 p=p0-(p0-exp(p0)+4)/(1-exp(p0)); if abs(p-p0)<=10^(-6) disp('|p-p0|=') disp(abs(p-p0)) disp('用牛顿法求得方程的正根为') disp(p); break; else p0=p; end end if i==n0 disp(n0) disp('次牛顿迭代后无法求出方程的解') end p1=-3; for i=1:n0 p=p1-(p1-exp(p1)+4)/(1-exp(p1)); if abs(p-p1)<=10^(-6) disp('|p-p1|=') disp(abs(p-p1)) disp('用牛顿法求得方程的负根为') disp(p); break; else p1=p; end end if i==n0 disp(n0) disp('次牛顿迭代后无法求出方程的解') end

雅可比迭代法的MATLAB程序

雅可比迭代法的MA TLAB程序: Function[x,k,index]=Jacobi(A,b,ep,it-max) % 求线性方程组的雅可比法; % A为方程组的系数矩阵; % b为方程组的右端项; % x为方程组的解; % ep为精度要求,缺省值为le-5; % it_max为最大迭代次数,缺省值为100; % k为迭代次数; % index 为指标变量,index=0表示计算失败,index=1表示计算成功; if nargin<4 it_max=100; end if nargin<3 ep=le-5; end n=length(A);k=0;x=zeros(n,1);y=zeros(n,1);index=1; while k<=it_max for i=1:n if abs (A(i,i))

matlab(迭代法_牛顿插值)

实验报告容: 一:不动点迭代法解方程 二:牛顿插值法的MA TLAB实现 完成日期:2012年6月21日星期四 数学实验报告一 日期:2012-6-21

所以,确定初值为x0=1 二:不断迭代 算法: 第一步:将f(x0)赋值给x1 第二步:确定x1-x0的绝对值大小,若小于给定的误差值,则将x1当做方程的解,否则回到第一步 编写计算机程序: clear f=inline('0.5*sin(x)+0.4'); x0=1; x1=f(x0); k=1; while abs(x1-x0)>=1.0e-6 x0=x1; x1=f(x0); k=k+1; fprintf('k=%.0f,x0=%.9f,x1=%.9f\n',k,x0,x1) end 显示结果如下: k=2,x0=0.820735492,x1=0.765823700 k=3,x0=0.765823700,x1=0.746565483 k=4,x0=0.746565483,x1=0.739560873 k=5,x0=0.739560873,x1=0.736981783

k=6,x0=0.736981783,x1=0.736027993 k=7,x0=0.736027993,x1=0.735674699 k=8,x0=0.735674699,x1=0.735543758 k=9,x0=0.735543758,x1=0.735495216 k=10,x0=0.735495216,x1=0.735477220 k=11,x0=0.735477220,x1=0.735470548 k=12,x0=0.735470548,x1=0.735468074 k=13,x0=0.735468074,x1=0.735467157 >>。。。 以下是程序运行截图:

matlab(迭代法-牛顿插值)

实验报告内容: 一:不动点迭代法解方程 二:牛顿插值法的MATLAB实现 完成日期:2012年6月21日星期四 数学实验报告一 日期:2012-6-21

所以,确定初值为x0=1 二:不断迭代 算法: 第一步:将f(x0)赋值给x1 第二步:确定x1-x0的绝对值大小,若小于给定的误差值,则将x1当做方程的解,否则回到第一步 编写计算机程序: clear f=inline('0.5*sin(x)+0.4'); x0=1; x1=f(x0); k=1; while abs(x1-x0)>=1.0e-6 x0=x1; x1=f(x0); k=k+1; fprintf('k=%.0f,x0=%.9f,x1=%.9f\n',k,x0,x1) end 显示结果如下: k=2,x0=0.820735492,x1=0.765823700 k=3,x0=0.765823700,x1=0.746565483 k=4,x0=0.746565483,x1=0.739560873 k=5,x0=0.739560873,x1=0.736981783

k=6,x0=0.736981783,x1=0.736027993 k=7,x0=0.736027993,x1=0.735674699 k=8,x0=0.735674699,x1=0.735543758 k=9,x0=0.735543758,x1=0.735495216 k=10,x0=0.735495216,x1=0.735477220 k=11,x0=0.735477220,x1=0.735470548 k=12,x0=0.735470548,x1=0.735468074 k=13,x0=0.735468074,x1=0.735467157 >>。。。 以下是程序运行截图:

matlab迭代法代码

matlab迭代法代码 1、%用不动点迭代法求方程 x-e^x+4=0的正根与负根,误差限是10^-6% disp('不动点迭代法'); n0=100; p0=-5; for i=1:n0 p=exp(p0)-4; if abs(p-p0)<=10^(-6) if p<0 disp('|p-p0|=') disp(abs(p-p0)) disp('不动点迭代法求得方程的负根为:') disp(p); break; else disp('不动点迭代法无法求出方程的负根.') end else p0=p; end end if i==n0 disp(n0)

disp('次不动点迭代后无法求出方程的负根') end p1=1.7; for i=1:n0 pp=exp(p1)-4; if abs(pp-p1)<=10^(-6) if pp>0 disp('|p-p1|=') disp(abs(pp-p1)) disp('用不动点迭代法求得方程的正根为') disp(pp); else disp('用不动点迭代法无法求出方程的正根'); end break; else p1=pp; end end if i==n0 disp(n0) disp('次不动点迭代后无法求出方程的正根') end 2、%用牛顿法求方程 x-e^x+4=0的正根与负根,误差限是10^-6 disp('牛顿法')

n0=80; p0=1; for i=1:n0 p=p0-(p0-exp(p0)+4)/(1-exp(p0)); if abs(p-p0)<=10^(-6) disp('|p-p0|=') disp(abs(p-p0)) disp('用牛顿法求得方程的正根为') disp(p); break; else p0=p; end end if i==n0 disp(n0) disp('次牛顿迭代后无法求出方程的解') end p1=-3; for i=1:n0 p=p1-(p1-exp(p1)+4)/(1-exp(p1)); if abs(p-p1)<=10^(-6) disp('|p-p1|=') disp(abs(p-p1)) disp('用牛顿法求得方程的负根为')

不动点迭代法matlab程序

不动点迭代法 实验目的:通过本实验学习利用MATLAB不动点迭代法,抛物线法,斯特芬森迭代法解非线性方程组,及其编程实现,培养编程与上机调试能力。 实验要求:①上机前充分准备,复习有关内容,写出计算步骤,查对程序; ②完成实验后写出完整的实验报告,内容应该包括:所用的算法语言, 算法步骤陈述,变量说明,程序清单,输出计算结果,结果分析等等; ③用编好的程序在Matlab环境中执行。 迭代法 MATLAB程序: function pwxff(f,x0,x1,x2,d,n) f=inline(f); x(1)=x0; x(2)=x1; x(3)=x2; w1=(f(x(2))-f(x(3)))/(x(2)-x(3)); t1=(f(x(1))-f(x(3)))/(x(1)-x(3)); t2=(f(x(1))-f(x(2)))/(x(1)-x(2)); w2=1/(x(1)-x(2))*(t1-t2); w=w1+w2*(x(3)-x(2)); for k=3:n x(k+1)=x(k)-2*f(x(k))/(w+sqrt(w^2-4*f(x(k))*w2));

if abs(x(k+1)-x(k))

相关文档
相关文档 最新文档