《廣東工業(yè)大學(xué)數(shù)值計算試卷及答案.doc》由會員分享,可在線閱讀,更多相關(guān)《廣東工業(yè)大學(xué)數(shù)值計算試卷及答案.doc(21頁珍藏版)》請在裝配圖網(wǎng)上搜索。
1、
數(shù)值計算引論
學(xué) 院 機電工程學(xué)院
專 業(yè) 機械設(shè)計制造及其自動化
年級班別 2014級(6)班
學(xué) 號 3114000271
學(xué)生姓名 劉就杰
2016年 11 月
一 編寫雅可比迭代法求解線性方程組的程
2、序,要求附有算例(20分)。
(可能的算例包括基本的驗證性算例、方程系數(shù)隨機生成的一般算例、用于算法對比的比較性算例等,對各算例的結(jié)果進行分析。)
雅可比迭代法的matlab程序如下
function x=Jacobi(A,b,x0,tol)
%雅可比迭代法解線性方程組
%A為系數(shù)矩陣,b為右端項,x0為初始向量,tol為誤差精度
sprintf(USAGE:Jacobi(A,b,x0,tol))
D=diag(diag(A));%diag(x) 返回由向量x的元素構(gòu)成的對角矩陣
U=triu(A,1);%triu(A)提取矩陣A的上三角部分生成上三角矩陣
L=tril(
3、A,-1);%tril(A)提取矩陣A的下三角部分生成下三角矩陣
B=-D\(L+U);%B為迭代矩陣
dl=D\b;
x=B*x0+dl;
n=1;
while norm(x-x0)>=tol
x0=x;
x=B*x0+dl;
n=n+1;
end
n %n為迭代次數(shù)
高斯-賽德爾迭代法的matlab程序如下:
function x=Guass_seidel(A,b,x0,tol)
%高斯-賽德爾迭代法解線性方程組
%A為系數(shù)矩陣,b為右端項,x0為初始向量,tol為誤差精度
sprintf(USAGE:Guass_seidel
4、(A,b,x0,tol))
D=diag(diag(A));%diag(x) 返回由向量x的元素構(gòu)成的對角矩陣
U=triu(A,1);%triu(A)提取矩陣A的上三角部分生成上三角矩陣
L=tril(A,-1);%tril(A)提取矩陣A的下三角部分生成下三角矩陣
G=-(D+L)\U;%G為迭代矩陣
dl=(D+L)\b;
x=G*x0+dl;
n=1;
while norm(x-x0)>=tol
x0=x;
x=G*x0+dl;
n=n+1;
end
n %n為迭代次數(shù)
調(diào)用編好的程序求解方程組:
A=[5 -1 -1 -1
5、 ;-1 10 -1 -1;-1 -1 5 -1;-1 -1 -1 10];
b=[-4;12;8;34];
x0=[0;0;0;0];
tol=1e-6;
x=Jacobi(A,b,x0,tol)
x=Guass_seidel(A,b,x0,tol)
實驗結(jié)果如下:
ans =
USAGE:Jacobi(A,b,x0,tol)
n =
20
x =
1.0000
2.0000
3.0000
4.0000
ans =
USAGE:Guass_seidel(A,b,x0,t)
n =
12
x =
6、 1.0000
2.0000
3.0000
4.0000
取相同的初始值達到同樣的精度10-6,雅可比迭代需要迭代20次,而高斯-賽德爾迭代法只需12次。
實驗總結(jié):通過這次實驗,對雅可比迭代法以及高斯-賽德爾迭代法求解線性方程組的基本原理有了進一步的理解,同時了解了雅可比和高斯-賽德爾迭代法的優(yōu)點,即雅可比和高斯-賽德爾在求解線性方程組的過程中具有更快的收斂速度,而高斯-賽德爾比雅可比的收斂速度更快(即取相同的初始值,達到同樣精度所需的迭代次數(shù)較少)。
二 編寫分段二次拉格朗日插值的程序,要求附有算例(20分)。
(對 在節(jié)點0,0.2
7、,0.4,0.6,0.8,1.0上進行插值,求x=0.7處的值,繪出被插值函數(shù)與插值函數(shù)的圖形,予以對比。)
建立如下拉格朗日插值函數(shù):
function y=lagrange(x0,y0,x);
n=length(x0);
m=length(x);
for i=1:m
z=x(i);
s=0.0;
for k=1:n
p=1.0;
for j=1:n
if j~=k
p=p*(z-x0(j))/(x0(k)-x0(j));
en
8、d
end
s=p*y0(k)+s;
end
y(i)=s;
end
在matlab中用拉格朗日插值求0.7處的值
>> exp(0.7)
ans =
2.013752707470477
>> lagrange(x,y,0.7)
ans =
2.013751960394443
繪出被插值函數(shù)與插值函數(shù)的圖形
x=[0 0.2 0.4 0.6 0.8 1.0];
y=exp(x);
x0=[-5:0.001:5];
y0=lagrange(x,y,x0);
9、y1=exp(x0);
plot(x0,y0,r)
hold on
plot(x0,y1,g)
紅線為插值函數(shù),綠線是被插值函數(shù),由圖像可以知道,在區(qū)間(-2,2)是較好擬合的,當(dāng)超出這個范圍后就會偏差越來越大。
三 編寫復(fù)化辛普森積分的程序,要求附有算例(20分)。
(對定積分,計算精度達到)
復(fù)化辛普森積分的程序
function S=bianfuhuasimpson(fx,a,b,eps,M)
% 變步長復(fù)合simpson求積公式
% fx -- 求積函數(shù)(函數(shù)文件)
% a, b -- 求積區(qū)間
% eps -- 計算
10、精度
% M--最大允許輸出劃分數(shù)
n=1;
h=(b-a)/n;
T1=h*(feval(fx,a)-feval(fx,b))/2;
Hn=h*feval(fx,(a+b)/2);
S1=(T1+2*Hn)/3;
n=2*n;
% 最好與倒數(shù)第三行保持一致(變步長)
while n<=M
T2=(T1+Hn)/2;
Hn=0;
h=(b-a)/n;
for j=1:n
x(j)=a+(j-1/2)*h;
y(j)=feval(fx,x(j));
Hn=Hn+y(j);
end
Hn=h*Hn;
11、
S2=(T2+2*Hn)/3;
fprintf( n=%2d S2=%-12.9f S2-S1=%-12.9f\n,n,S2,abs(S2-S1));
if abs(S2-S1)
12、h=0.1, 0.01, 0.001,0.0001),運行兩種算法的程序,并將結(jié)果繪制成圖形,進行比較、分析。若要解的精度達到,應(yīng)采取什么措施?)
程序:
%Euler法
F=x^2+x-y;
a=0;
b=0.5;
h=0.1;
n=(b-a)/h;
X=a:h:b;
Y=zeros(1,n+1);
Y(1)=0;
for i=2:n+1
x=X(i-1);
y=Y(i-1);
Y(i)=Y(i-1)+eval(F)*h;
end
%隱式Euler法
Y1=zeros
13、(1,n+1);
Y1(1)=0;
for i=2:n+1
x=X(i);
y=Y1(i-1);
Y1(i)=(y+x*h^2+h*x)/(h+1);
end
%準確解
temp=[];
f=dsolve(Dy=x^2+x-y,y(0)=0,x);
df=zeros(1,n+1);
for i=1:n+1
temp=subs(f,x,X(i));
df(i)=double(vpa(temp));
end
disp( 步長 Euler法 隱式Euler法 準確值);
disp([X,Y,Y1,df])
14、;
%畫圖觀察效果
figure;
plot(X,df,k-,X,Y,--r,X,Y1,.-b);
grid on;
title(Euler法和隱式Euler法解常微分方程);
legend(準確值,Euler法,隱式Euler法);
結(jié)果為:
1.當(dāng)h=0.1時
Euler法和,隱式Euler法與實際值相差較大,誤差隨x增大而增大。由圖可以看出h=0.1時,不能得到精確結(jié)果。Euler法精度為,隱式Euler法精度為。
2.當(dāng)h=0.01時
Euler法和隱式Euler法
15、與實際值相差不大,誤差隨x增大而增大,但還是相對準確。由圖可以看出h=0.01時,能得到精確結(jié)果.Euler法精度為,隱式Euler法精度為。
3.當(dāng)h=0.001時,
Euler法和,隱式Euler法與實際值十分相近,誤差很小。由圖可以看出,三條函數(shù)曲線重疊。h=0.001時,能得到精確結(jié)果,,.Euler法精度為,隱式Euler法精度為。
4.當(dāng)h=0.0001時
Euler法和,隱式Euler法與實際值十分相近,誤差很小。由圖可以看出,三條函數(shù)曲線重疊。h=0.0001時,能得到精確結(jié)果,,Eul
16、er法精度為,Euler預(yù)測-校正法精度為。
總體來看,隱式Euler法和Euler法精度沒有大的差別,步數(shù)每縮小10倍,精度提高10倍。若要解的精度達到,可以取h=0.000001,隱式Euler法和Euler法能達到要求。
五 編寫簡單迭代法求解非線性方程的程序,要求附有算例(20分)。
(對方程,求x=0.5附近的根,寫出5種迭代格式,至少兩種格式收斂;對每種格式的計算結(jié)果進行分析、比較,并盡可能把不收斂格式轉(zhuǎn)化為收斂形式)
format long
f=inline(此處用用下面構(gòu)造高的5種式子替代)
disp(x=);
17、x=feval(f,0.5);
disp(x);
Eps=1E-5;
i=1;
while 1
x0=x;
i=i+1;
x=feval(f,x);
disp(x);
if ~isreal(x)%不是實數(shù)不進行迭代
disp(出現(xiàn)復(fù)數(shù))%提示出現(xiàn)復(fù)數(shù)
break;
end
if x>1E10;%發(fā)散不進行迭代
disp(發(fā)散)%提示發(fā)散
break;
end;
if abs(x-x0)
18、> Untitled3
f =
內(nèi)聯(lián)函數(shù):
f(x) = (3*x-x^2-1)^(1/3)
x=
0.629960524947437
0.789995893719114
0.906899308403962
0.964856599317808
0.987723757473192
0.995840405544922
0.998605758099246
19、
0.999534387968673
0.999844699607740
0.999948222482311
0.999982739635881
0.999994246412883
0.999998082122915
i =
13
x =
0.999998082122915
由程序運行結(jié)果知道該方法收斂,用了13次就達到求解。
(2)xk+1=(3*xk-xk^3-1)^(1/2)
>> Untitled3
f =
內(nèi)聯(lián)函數(shù):
20、 f(x) = (3*x-x^3-1)^(1/2)
x=
0.612372435695794
0.779408521701848
0.929920594922635
0.992779330733286
0.999921978095163
0.999999990869111
1.000000000000000
i =
7
x =
1.000000000000000
由程序運行結(jié)果知道,該方法收斂,用了7次就達到求解。相比第一種方式更快,更精確。
21、
(3)xk+1=(xk^3+xk^2+1)/3
>> Untitled3
f =
內(nèi)聯(lián)函數(shù):
f(x) = ((x^3+x^2+1)/3)
x=
0.458333333333333
0.435450424382716
0.424061968869700
0.418695667957701
0.416235317082047
0.415121791135062
0.414620807126294
0.414396016061541
0.414295274559236
22、 0.414250151156419
0.414229944730161
0.414220897204816
i =
12
x =
0.414220897204816
由程序運行結(jié)果知道該方法收斂,用了12次就達到求解。
(4)xk+1=(-1)/(xk^2+xk-3)
>> Untitled3
f =
內(nèi)聯(lián)函數(shù):
f(x) = ((-1)/(x^2+x-3))
x=
0.444444444444444
0.42408376
23、9633508
0.417350219079978
0.415201597459088
0.414523917231877
0.414310962753779
0.414244121582387
0.414223149438888
0.414216569954722
i =
9
x =
0.414216569954722
由程序運行結(jié)果知道該方法收斂,用了9次就達到求解。比方式3快且精確。
(5)xk+1=(xk-(xk^3+xk^2-3*xk
24、+1)/(3*xk^2+2*xk-3))
>> Untitled3
f =
內(nèi)聯(lián)函數(shù):
f(x) = (x-(x^3+x^2-3*x+1)/(3*x^2+2*x-3))
x=
0.400000000000000
0.413953488372093
0.414213470906413
0.414213562373084
i =
4
x =
0.414213562373084
由程序運行結(jié)果知道該方法收斂,用了次就達到求解。其速度與效果優(yōu)于第3,第4種方式。
廣東工業(yè)大學(xué)試卷用紙,共21頁,第21頁