目录。上机作业12上机作业26上机作业310上机作业411上机作业5:
上机作业1:
x1、分别用不动点迭代与newton法求解方程xe20的正根与负根。
2、用newton法求解方程xsinx0的根,再用steffensen’s method加速其收敛。要求:1)抄题;2)迭代公式(初值)或简单原理;3)程序,结果;(打印)4)结果分析。
解:1、方法一:不动点迭代法:1)迭代公式:pn+1=g(pn),n=0,1,……2)程序:
function[p,k]=mystablepoint(n,p0,tol);if(nargin==2)tol=1.0e-4;endk=1;while kdisp(k);disp(p)function f=g(x):if x>=-2&x<=0f=exp(x)-2;elseif x>=0&x<=2f=log(x+2);end
3)结果:[p,k]=mystablepoint(1000,-1);
p,k]=mystablepoint (1000,1);8
4)结果分析:不动点迭代法的原理比较简单,易于程序实现。方法二:newton法:1)迭代公式:pn1pn2)程序:
function xn = mynewton(f,a,b,tol)if(nargin==3)tol=
endf1=subs(sym(f),findsym(sym(f)),a);f2=subs(sym(f),findsym(sym(f)),b);if(f1==0)xn=a;
endif(f2==0)xn=b;
endif(f1*f2>0)
disp('两端点函数值大于0');return;
elseeps=1;
fun=diff(sym(f));
fa=subs(sym(f),findsym(sym(f)),a);fb=subs(sym(f),findsym(sym(f)),b);
f(pn)f'(pn)
n=0,1,……
dfa=subs(sym(fun),findsym(sym(fun)),a);dfb=subs(sym(fun),findsym(sym(fun)),b);if(dfa>dfb)xn=a-fa/dfa;else
xn=b-fb/dfb;end
while(eps>tol)x1=xn;
fx=subs(sym(f),findsym(sym(f)),x1);dfx=subs(sym(fun),findsym(sym(fun)),x1);xn=x1-fx/dfx;eps=abs(xn-x1);endend
3)结果:xn = mynewton('x+2-exp(x)',2,-1)
xn =-1.8414
xn=mynewton('x+2-exp(x)',1,1.5)xn =
4)结果分析:牛顿迭代法收敛于根的速度最快) steffensen加速原理:pk1)
p0k1)
p1p2k1)k1)p0
k1)k1)
2p1p0k1)
newton法结果:[xn,n]=newtonm('x-sin(x)',1,1)
xn =1.4990e-006n =32
2)程序:function [xs,n]=mysteffensen(f,a,b,eps)
if(nargin==3)
tol=end
f1=subs(sym(f),findsym(sym(f)),a);f2=subs(sym(f),findsym(sym(f)),b);if(f1==0)
xs=a;endif(f2==0)
xs=b;endif(f1*f2>0)
return;
elseeps=1;
fa=subs(sym(f),findsym(sym(f)),a);fb=subs(sym(f),findsym(sym(f)),b);faa=subs(sym(f),findsym(sym(f)),fa+a);xs=a-fa*fa/(faa-fa);n=0;while(eps>tol)n=n+1;x1=xs;
fx=subs(sym(f),findsym(sym(f)),x1);v=fx+x1;
fxx=subs(sym(f),findsym(sym(f)),v);xs=x1-fx*fx/(fxx-fx);eps=abs(xs-x1);end
end3)结果:[xs,n]=mysteffensen('x-sin(x)',1,1)
xs =1.7141e-006n =
3)结果分析:steffensen加速极大地改善了原迭代的收敛性质,它加速了牛顿法的收敛速度,但计算量有所增加。
上机作业2:
分别用jacobi、seidel、sor(w=1.1,1.2,1.3,1.4,1.5)方法求解方程组ax=b,这里:
abe=10x000’
要求:抄题,公式,程序、计算结果(终止迭代步数k、近似解x(k)),结果分析(三种迭代列表)。解:
方法一:jacobi迭代法1)公式:xk1)i
1aiii1nijbi
aj1xk)j
aji1ijxj)
k)i1,2,,n,k0,1,2,x(k+1)=(i-d-1a)x(k)+ d-1b
2)程序:b=diag([-12,-12,-12,-12,-12,-12,-12,-12,-12,-12]);
d=ones(10);a=b+d;
b=[1;1;1;1;1;1;1;1;1;1];x0=[0;0;0;0;0;0;0;0;0;0];x=x0;for k=1:100
for i=1:10
x(i)=(b(i)-a(i,:)x0)/a(i,i)+x0(i);end
if norm(x-x0)<0.00001disp(x);disp(k)return;endx0=x;end
3)计算结果:迭代步数k=53;近似解x(53)=[0.5000; -0.
5000; -0.5000; -0.5000;-0.
5000;-0.5000; -0.5000; -0.
5000; -0.5000; -0.5000]方法二:
gauss-seidel迭代法1)公式:x
k1)i1aiii1nij
biaj1x
k1)jaji1ij
xj)k)
i1,2,,n,k0,1,2,x(k+1)=[l+d)-1u]x(k)+ l+d)-1b
2)程序:b=diag([-12,-12,-12,-12,-12,-12,-12,-12,-12,-12]);
d=ones(10);a=b+d;
b=[1;1;1;1;1;1;1;1;1;1];x0=[0;0;0;0;0;0;0;0;0;0];x=x0;for k=1:100
for i=1:10
x(i)=(b(i)-a(i,:)x)/a(i,i)+x(i);end
if norm(x-x0)<0.00001disp(x);disp(k)return;endx0=x;end
3)计算结果:迭代步数k=29;近似解x(53)=[0.5000; -0.
5000; -0.5000; -0.5000;-0.
5000;-0.5000; -0.5000; -0.
5000; -0.5000; -0.5000]方法三:
sor迭代法1)公式:x
k1)i1)xk)i
aiii1nijbi
aj1xk1)j
aji1ijxj)
k)x(k+1)=(d+ωl)-1[(1-ω)d-ωu]x(k)+ωd+ωl)-1b2)程序:○
1编写m函数文件 y=f(w)
b=diag([-12,-12,-12,-12,-12,-12,-12,-12,-12,-12]);d=ones(10);a=b+d;
b=[1;1;1;1;1;1;1;1;1;1];x0=[0;0;0;0;0;0;0;0;0;0];x=x0;for k=1:100
for i=1:10
x(i)=w*(b(i)-a(i,:)x)/a(i,i)+x(i);end
if norm(x-x0)<0.00001disp(x);disp(k)return;endx0=x;end
2调用分别求出迭代结果f(1.1)f(1.2)f(1.3)f(1.4)
f(1.5)
3)计算结果:○1当w=1.1时,终止迭代步数k=24;近似解x(17)=[0.5000;
2当w=1.2时,○终止迭代步数k=19;近似解x(17)=[0.5000;-0.
5000; -0.5000; -0.5000; -0.
5000;-0.5000; -0.5000; -0.
5000;-0.5000; -0.5000]
3当w=1.3时,终止迭代步数k=14;近似解x○(13)=[0.5000;-0.
5000; -0.5000; -0.5000; -0.
5000;-0.5000; -0.5000; -0.
5000;-0.5000; -0.5000]
4当w=1.4时,终止迭代步数k=13;近似解x○(12)=[0.5000;-0.
5000; -0.5000; -0.5000; -0.
5000;-0.5000; -0.5000; -0.
5000;-0.5000; -0.5000]
5当w=1.5时,终止迭代步数k=17;近似解x○(16)=[0.5000;-0.
5000; -0.5000; -0.5000; -0.
5000;-0.5000; -0.5000; -0.
5000;-0.5000; -0.5000]
表1:三种迭代方法对比表。
迭代方法jacobi迭代法gauss-seidel迭代法。
w=1.1w=1.2
sor迭。w=1.3
代法。w=1.4w=1.5
131714迭代次数。
x=[-0.5000;-0.5000;-0.
5000;-0.5000;-0.5000;-0.
5000; -0.5000; -0.5000; -0.
5000; -0.5000]
近似解。结果分析:由以上计算过程可以看出,这三个迭代法都是收敛的,但不同的是jacobi迭代法的收敛速度是最慢的,sor迭代法收的敛速度是最快,且随着松弛因子w的变化而不同,从表中可以看出,当w=1.
4的时候收敛速度最快。
方程的近似解为x=[-0.5000;-0.5000;-0.
5000;-0.5000;-0.5000;-0.
5000; -0.5000;-0.5000; -0.
5000; -0.5000]。
上机作业3:
用newton法求方程组。
22x10xy802xyx10y80
在(0,0)附近的根。
要求:1)抄题;2)迭代公式(初值)或简单原理;3)程序,结果;(打印)4)结果分析。解:
1)迭代公式:x(k+1)=x(k)-[f’(x(k)) 1f(x(k))2)程序:x0=[0;0];
for k=1:10
z=-[2*x0(1)-10,2*x0(2);x0(2)^2+1,2*x0(1)*x0(2)-10]\[x0(1)^2-10*x0(1)+x0(2)^2+8;x0(1)*x0(2)^2+x0(1)-10*x0(2)+8];x=x0+z;
if norm(z)<0.000001disp(x);disp(k)return;endx0=xend
3)计算结果见下表:
表2:迭代结果表。
kx(k)y(k)
4)结果分析:由上表可以看出,newton法迭代的收敛速度是比较快的,在第三迭代就收敛于方程组的解x*=(1,1)t。
上机作业4:
目的:观察lagrange插值的runge现象。对于函数f(x)1125x
进行lagrange
插值。取不同等分n(n=5,n=10),将区间[-1,1]n等分,取等距节点。把插值多项式的曲线画在同一张图上进行比较。解:
n1)公式:ln(x)=
k0ykln,k(x)
2)程序:function f = language(x,y)
syms t;
if(length(x) =length(y))
n = length(x);
elsedisp('x和y的维数不等!')return;
endf = 0.0;for(i = 1:n)
l = y(i);for(j = 1:i-1)
l = l*(t-x(j))/x(i)-x(j));end;for(j = i+1:n)
l = l*(t-x(j))/x(i)-x(j));end;
f = f + l;simplify(f);
end数据:
n=5:x`f(x)n=10:xf(x)
3)结果:x1= -1:0.4:1;
> y1= [0.0385,0.1,0.5,0.5,0.1,0.0385];>f1 = language(x1,y1)f1 =
125/384*(-77/800*t-231/4000)*(t+1/5)*(t-1/5)*(t-3/5)*(t-1)+625/384*(1/4*t+1/4)*(t+1/5)*(t-1/5)*(t-3/5)*(t-1)-625/96*(5/8*t+5/8)*(t+3/5)*(t-1/5)*(t-3/5)*(t-1)+625/64*(5/12*t+5/12)*(t+3/5)*(t+1/5)*(t-3/5)*(t-1)-625/96*(1/16*t+1/16)*(t+3/5)*(t+1/5)*(t-1/5)*(t-1)+625/384*(77/4000*t+77/4000)*(t+3/5)*(t+1/5)*(t-1/5)*(t-3/5)>>collect(f1)ans =
5317/3072*t^2+7385/6144*t^4+145231/256000>> x2= -1:0.2:1;
> y2= [0.038,0.059,0.
1,0.111,0.5,1,0.
5,0.111,0.1,0.
059,0.038];>f2 = language(x2,y2)f2 =
78125/145152*(-19/100*t-19/125)*(t+3/5)*(t+2/5)*(t+1/5)*t*(t-1/5)*(t-2/5)*(t-3/5)*
t-4/5)*(t-1)-390625/72576*(59/200*t+59/200)*(t+3/5)*(t+2/5)*(t+1/5)*t*(t-1/5)*(t-2/5)*(t-3/5)*(t-4/5)*(t-1)+390625/8064*(1/4*t+1/4)*(t+4/5)*(t+2/5)*(t+1/5)*t*(t-1/5)*(t-2/5)*(t-3/5)*(t-4/5)*(t-1)-390625/2016*(37/200*t+37/200)*(t+4/5)*(t+3/5)*(t+1/5)*t*(t-1/5)*(t-2/5)*(t-3/5)*(t-4/5)*(t-1)+390625/864*(5/8*t+5/8)*(t+4/5)*(t+3/5)*(t+2/5)*t*(t-1/5)*(t-2/5)*(t-3/5)*(t-4/5)*(t-1)-390625/576*(t+1)*(t+4/5)*(t+3/5)*(t+2/5)*(t+1/5)*(t-1/5)*(t-2/5)*(t-3/5)*(t-4/5)*(t-1)+390625/576*(5/12*t+5/12)*(t+4/5)*(t+3/5)*(t+2/5)*(t+1/5)*t*(t-2/5)*(t-3/5)*(t-4/5)*(t-1)-390625/864*(111/1400*t+111/1400)*(t+4/5)*(t+3/5)*(t+2/5)*(t+1/5)*t*(t-1/5)*(t-3/5)*(t-4/5)*(t-1)+390625/2016*(1/16*t+1/16)*(t+4/5)*(t+3/5)*(t+2/5)*(t+1/5)*t*(t-1/5)*(t-2/5)*(t-4/5)*(t-1)-390625/8064*(59/1800*t+59/1800)*(t+4/5)*(t+3/5)*(t+2/5)*(t+1/5)*t*(t-1/5)*(t-2/5)*(t-3/5)*(t-1)+390625/72576*(19/1000*t+19/1000)*(t+4/5)*(t+3/5)*(t+2/5)*(t+1/5)*t*(t-1/5)*(t-2/5)*(t-3/5)*(t-4/5)>>collect(f2)ans =
插值多项式的曲线图。
n=5插值曲线n=10插值曲线。
上机作业5:
用romberg求积公式计算积分:解:
ecosxdx
x使精度达到10.
t1(h)t(h)hm
1)公式:tm()4tm(h)
2t(h),m1,2,m1m14
2)程序:function [r,k,t]=romberg(fun,a,b,tol)k=0;n=1;h=b-a;
t=h/2*(fun(a)+fun(b));err=1;
while err>=tolk=k+1;h=h/2;tmp=0;for i=1:n
tmp=tmp+fun(a+(2*i-1)*h);end
t(k+1,1)=t(k)/2+h*tmp;for j=1:k
t(k+1,j+1)=t(k+1,j)+(t(k+1,j)-t(k,j))/4^j-1);endn=n*2;
err=abs(t(k+1,k+1)-t(k,k));endr=t(k+1,4);实现过程:fun=@(x)exp(x)*cos(x)
r,k,t]=romberg(fun,1,3,1e-6)
3)计算结果如下表:步长h210.50.250.0625
0.03125-10.4102-10.
4031-10.4031-10.4031-10.
4031-10.4031所以其积分的近似值为-10.4031。
matlab作业
2011029170002王柳。a 一个问题的病态性如何,与求解它的算法有关系。错 b 无论问题是否病态,好的算法都会得到它好的近似解。错 c 计算中使用更高的精度,可以改善问题的病态性。错 d 用一个稳定的算法计算一个良态问题,一定会得到它好的近似解。对 e 浮点数在整个数轴上是均匀分布。错 f ...
matlab作业
matlab语言 第3次作业 字符串,单元数组和结构体 专业 海洋技术 海洋测绘方向 姓名 张体强学号 1026222 1 如何将一个char 数据类型的向量转化为相应的double 型数据类型的数据向量。从式1 到8,判断这些语句是否正确。如果它们正确,那么将产生什么结果?这题不要在电脑中做。1....
matlab作业
电子与通信工程学院。通信系统 实验报告。2013 2014 学年第1学期。调频 fm 系统调制解调 专业 通信工程。班级 通信111 班。学号 姓名 指导教师姓名陈多瑜。2013年 11 月日。1.频率调制或调频 fm 1 设调制信号为m t 调频信号的数学表达式为。例如 m t 的时域波形为。m ...