《數(shù)值分析》課程設(shè)計(jì)—作業(yè)

上傳人:仙*** 文檔編號:29758503 上傳時(shí)間:2021-10-08 格式:DOC 頁數(shù):54 大?。?.79MB
收藏 版權(quán)申訴 舉報(bào) 下載
《數(shù)值分析》課程設(shè)計(jì)—作業(yè)_第1頁
第1頁 / 共54頁
《數(shù)值分析》課程設(shè)計(jì)—作業(yè)_第2頁
第2頁 / 共54頁
《數(shù)值分析》課程設(shè)計(jì)—作業(yè)_第3頁
第3頁 / 共54頁

下載文檔到電腦,查找使用更方便

15 積分

下載資源

還剩頁未讀,繼續(xù)閱讀

資源描述:

《《數(shù)值分析》課程設(shè)計(jì)—作業(yè)》由會(huì)員分享,可在線閱讀,更多相關(guān)《《數(shù)值分析》課程設(shè)計(jì)—作業(yè)(54頁珍藏版)》請?jiān)谘b配圖網(wǎng)上搜索。

1、 《數(shù)值分析》課程設(shè)計(jì)—作業(yè) (本文檔內(nèi)的有些運(yùn)行結(jié)果,限于篇幅,使文檔結(jié)構(gòu)更和諧、緊湊,已做相關(guān)的改動(dòng),程序代碼沒變) 實(shí)驗(yàn)一 1.1 水手、猴子和椰子問題:五個(gè)水手帶了一只猴子來到南太平洋的一個(gè)荒島上,發(fā)現(xiàn)那里有一大堆椰子。由于旅途的顛簸,大家都很疲憊,很快就入睡了。第一個(gè)水手醒來后,把椰子平分成五堆,將多余的一只給了猴子,他私藏了一堆后便又去睡了。第二、第三、第四、第五個(gè)水手也陸續(xù)起來,和第一個(gè)水手一樣,把椰子分成五堆,恰多一只猴子,私藏一堆,再去入睡,天亮以后,大家把余下的椰子重新等分成五堆,每人分一堆,正好余一只再給猴子,試問原先共有幾

2、只椰子? 試分析椰子數(shù)目的變化規(guī)律,利用逆向遞推的方法求解這一問題。 解: 1、 問題分析:對于本題,比較簡單,我們只需要判斷原來椰子的個(gè)數(shù)及每個(gè)人私藏了一份之后剩下的是否能被5除余1,直到最后分完。 function fentao(n) a=cat(1,7); for j=n:-1:1 a(1)=j;i=1; while i<7 a(i+1)=4*(a(i)-1)/5; i=i+1; end if a(7)==fix(a(7)) a, end end end 2、 問題

3、求解:通過matlab建立M文件,有如下程序: n=input(input n:); for x=1:n p=5*x+1; for k=1:5 p=5*p/4+1; end if p==fix(p) break; end end disp([x,p]) 或者 input n:2000 1023 15621 >> fenta

4、o(20001) a = 15621 12496 9996 7996 6396 5116 4092 對于第一個(gè)程序,n取2000;對于第二個(gè)程序,n取20001,就能得到我們想要的結(jié)果,即原先一共有15621個(gè)椰子,最終平均每人得4092個(gè)椰子。 1.2 當(dāng)時(shí),選擇穩(wěn)定的算法計(jì)算積分. 解: 一、問題分析:由知: 以及: 得遞推關(guān)系: , 但是通過仔細(xì)觀察就能知道上述遞推公式每一步都將誤差放大十倍,即使初始誤 差很小,但是誤差的傳播會(huì)逐步擴(kuò)大,也就

5、是說用 它構(gòu)造的算法是不穩(wěn)定的,因此我們改進(jìn)上述遞推公式(算法)如下: 通過比較不難得出該誤差是逐步縮小的,即算法是穩(wěn)定的。 2、 問題求解:為了利用上面穩(wěn)定的算法,需要我們估計(jì)初值的值。 因?yàn)? 所以當(dāng)n=100的時(shí),我們有: 9.090909090909091e-0049.645173882220725e-004 于是可取他們的平均值,有=9.368041486564908e-004,利用上述穩(wěn)定算法,可求得相應(yīng)的值和程序如下(限于篇幅,這里只給出了前后各十個(gè)連續(xù)的數(shù)據(jù),詳細(xì)的數(shù)據(jù)會(huì)連同作業(yè)一并上交):

6、 積分計(jì)算對照表 1 n 計(jì)算值 準(zhǔn)確值 誤差 1 0.059182761 0.059182761 0 2 0.040292947 0.040292947 0 3 0.029402893 0.029402893 1.17997E-16 4 0.022708713 0.022708713 1.5278E-16 5 0.018333956 0.018333956 0 6 0.01531285 0.01531285 1.13285E-16 7 0.013125037 0.013125037 0 8 0.0

7、114765 0.0114765 0 9 0.010193063 0.010193063 0 91 0.00099999 0.00099999 2.68937E-11 92 0.00098911 0.00098911 2.71896E-10 93 0.000978464 0.000978464 2.74854E-09 94 0.000968045 0.000968045 2.77813E-08 95 0.000957845 0.000957846 2.80771E-07 96 0.000947862 0.00094785

8、9 2.83729E-06 97 0.000938051 0.000938078 2.86687E-05 98 0.000928766 0.000928497 0.000289646 99 0.000916421 0.00091911 0.002926039 100 0.000936804 0.000909911 0.029556214 I=cat(1,100); J=cat(1,10

9、0); K=cat(1,100); I(100)=9.368041486564908e-004; format long; %求近似值 for n=99:-1:1 I(n)=((1-exp(-n))/n-I(n+1))/10; end %求精確值 for n=1:100 syms x; k=n; J(k)=int(exp(-k*x)/(exp(-x)+10),x,0,1); end %求誤差 for n=1:100 K(n)=abs((I(n)-J(n))/J(n)); end n=1:100; A=[n;I;J;K

10、];B=A xlswrite(1_2.xls,B) >> format long >> min=(1-exp(-100))/(11*100), max=(1-exp(-100))/(100*(exp(-1)+10)) min = 9.090909090909091e-004 max = 9.645173882220725e-004 >> I100=(min+max)/2 I100 = 9.368041486564908e-004 1.3 繪制靜態(tài)和動(dòng)態(tài)的

11、Koch分形曲線 問題描述:從一條直線段開始,將線段中間的三分之一部分用一個(gè)等邊三角形的另兩條邊代替,形成具有5個(gè)結(jié)點(diǎn)的新的圖形;在新的圖形中,又將圖中每一直線段中間的三分之一部分都用一個(gè)等邊三角形的另兩條邊代替,再次形成新的圖形,這時(shí),圖形中共有17個(gè)結(jié)點(diǎn)。這種迭代繼續(xù)進(jìn)行下去可以形成Koch分形曲線。在迭代過程中,圖形中的結(jié)點(diǎn)將越來越多,而曲線最終顯示細(xì)節(jié)的多少取決于所進(jìn)行的迭代次數(shù)和顯示系統(tǒng)的分辨率。Koch分形曲線的繪制與算法設(shè)計(jì)和計(jì)算機(jī)實(shí)現(xiàn)相關(guān)。 圖1.1 Koch曲線的形成過程 解: 一、算法分析: 考慮由直線段(2個(gè)點(diǎn))產(chǎn)生第一個(gè)圖形(5個(gè)點(diǎn))的過程。上

12、圖中,設(shè)和分別為原始直線段的兩個(gè)端點(diǎn),現(xiàn)需要在直線段的中間依次插入三個(gè)點(diǎn),,。顯然位于線段三分之一處,位于線段三分之二處,點(diǎn)的位置可看成是由點(diǎn)以點(diǎn)為軸心,逆時(shí)針旋轉(zhuǎn)600而得。旋轉(zhuǎn)由正交矩陣 實(shí)現(xiàn)。 算法根據(jù)初始數(shù)據(jù)(和點(diǎn)的坐標(biāo)),產(chǎn)生圖1中5個(gè)結(jié)點(diǎn)的坐標(biāo)。結(jié)點(diǎn)的坐標(biāo)數(shù)組形成一個(gè)矩陣,矩陣的第一行為的坐標(biāo),第二行為的坐標(biāo)……,第五行為的坐標(biāo)。矩陣的第一列元素分別為5個(gè)結(jié)點(diǎn)的坐標(biāo),第二列元素分別為5個(gè)結(jié)點(diǎn)的坐標(biāo)。 進(jìn)一步考慮該曲線形成過程中結(jié)點(diǎn)數(shù)目的變化規(guī)律。設(shè)第次迭代產(chǎn)生的結(jié)點(diǎn)數(shù)為,第次迭代產(chǎn)生的結(jié)點(diǎn)數(shù)為,則和中間的遞推關(guān)系為。 2、 問題求解:

13、 ⑴、程序代碼 p=[0 0;10 0]; %P為初始兩個(gè)點(diǎn)的坐標(biāo),第一列為x坐標(biāo),第二列為y坐標(biāo) n=2; %n為結(jié)點(diǎn)數(shù) A=[cos(pi/3) -sin(pi/3);sin(pi/3) cos(pi/3)]; %旋轉(zhuǎn)矩陣 for k=1:6 d=diff(p)/3; %diff計(jì)算相鄰兩個(gè)點(diǎn)的坐標(biāo)之差,得到相鄰兩點(diǎn)確定的向量 %則d就計(jì)算出每個(gè)向量長度的三分之一,與題中將線段三等分對應(yīng) m=4*n-3; %迭代公式 q=p(1:n-1,:); %以原點(diǎn)為起點(diǎn),前n-1個(gè)點(diǎn)的坐標(biāo)為終

14、點(diǎn)形成向量 p(5:4:m,:)=p(2:n,:); %迭代后處于4k+1位置上的點(diǎn)的坐標(biāo)為迭代前的相應(yīng)坐標(biāo) p(2:4:m,:)=q+d; %用向量方法計(jì)算迭代后處于4k+2位置上的點(diǎn)的坐標(biāo) p(3:4:m,:)=q+d+d*A; %用向量方法計(jì)算迭代后處于4k+3位置上的點(diǎn)的坐標(biāo) p(4:4:m,:)=q+2*d; %用向量方法計(jì)算迭代后處于4k位置上的點(diǎn)的坐標(biāo) n=m; subplot(2,3,k) plot(p(:,1),p(:,2),r) %繪出每相鄰兩個(gè)點(diǎn)的連線 axis([0,3,0

15、,3]) %迭代后新的結(jié)點(diǎn)數(shù)目 %axis off end ⑵、運(yùn)行結(jié)果: 實(shí)驗(yàn)二 2.1 小行星軌道問題:一天文學(xué)家要確定一顆小行星繞太陽運(yùn)行的軌道,他在軌道平面內(nèi)建立以太陽為原點(diǎn)的直角坐標(biāo)系,在五個(gè)不同的對小行星作了五次觀察,測得軌道上五個(gè)點(diǎn)的坐標(biāo)數(shù)據(jù)(單位:萬公里)如下表所示: P1 P2 P3 P4 P5 X坐標(biāo) 5

16、3605 58460 62859 66662 68894 Y坐標(biāo) 6026 11179 16954 23492 68894 由開普勒第一定律知,小行星軌道為一橢圓,橢圓的一般方程可表示為: 現(xiàn)需要建立橢圓的方程以供研究。 (1) 分別將五個(gè)點(diǎn)的數(shù)據(jù)代入橢圓一般方程中,寫出五個(gè)待定系數(shù)滿足的等式,整理后寫出線性方程組 以及方程組的系數(shù)矩陣和右端項(xiàng)b; (2) 用MARLAB求低階方程的指令A(yù)\b求出待定系數(shù); (3) 分別用直接法、Jacobi迭代法、Gauss-Seidel迭代法求出待定系數(shù). 解: (1)、①、程序代碼: X=[53605 584

17、60 62859 66662 68894]; Y=[6026 11197 16954 23492 68894]; B=[X.*X;2*X.*Y;Y.*Y;2*X;2*Y]; %B是系數(shù)矩陣的轉(zhuǎn)置 A=B;%A是系數(shù)矩陣 A=vpa(A,16) syms a1 a2 a3 a4 a5; a=[a1;a2;a3;a4;a5]; b=[-1;-1;-1;-1;-1]; b=sym(A*a) ②、運(yùn)行結(jié)果: A = [ 2873496025.0, 646047460.0, 36312676.0, 107

18、210.0, 12052.0] [ 3417571600.0, 1309153240.0, 125372809.0, 116920.0, 22394.0] [ 3951253881.0, 2131422972.0, 287438116.0, 125718.0, 33908.0] [ 4443822244.0, 3132047408.0, 551874064.0, 133324.0, 46984.0] [ 4746383236.0, 9492766472.0, 4746383236.0, 137788.0, 137788.0] b

19、 = 2873496025.0*a1 + 646047460.0*a2 + 36312676.0*a3 + 107210.0*a4 + 12052.0*a5 3417571600.0*a1 + 1309153240.0*a2 + 125372809.0*a3 + 116920.0*a4 + 22394.0*a5 3951253881.0*a1 + 2131422972.0*a2 + 287438116.0*a3 + 125718.0*a4 + 33908.0*a5 4443822244.0*a1 + 3132047408.0*a2 + 551874064.0

20、*a3 + 133324.0*a4 + 46984.0*a5 4746383236.0*a1 + 9492766472.0*a2 + 4746383236.0*a3 + 137788.0*a4 + 137788.0*a5 即是方程組的形式如下: -1=2873496025.0*a1 + 646047460.0*a2 + 36312676.0*a3 + 107210.0*a4 + 12052.0*a5 -1=3417571600.0*a1 + 1309153240.0*a2 + 125372809.0*a3 + 116920.0*a4 + 22394.0*a

21、5 -1=3951253881.0*a1 + 2131422972.0*a2 + 287438116.0*a3 + 125718.0*a4 + 33908.0*a5 -1=4443822244.0*a1 + 3132047408.0*a2 + 551874064.0*a3 + 133324.0*a4 + 46984.0*a5 -1=4746383236.0*a1 + 9492766472.0*a2 + 4746383236.0*a3 + 137788.0*a4 + 137788.0*a5 (2) 、用指令A(yù)\b求出待定系數(shù)的值如下: >> b=[-1;-1;-

22、1;-1;-1];A;a=A\b a = 0.000000000080151337137572070705034479362305 -0.000000000099919136582879577308894350837601 -0.00000000017639794583181340912897589334276 -0.000012556065177690492529375516352892 0.000015497775048603393791640753306649 (3) 、1)、直接法

23、 ①、程序代碼: function [RA,RB,n,X]=bilizy(A,b) B=[A,b];n=length(b); RA=rank(A);RB=rank(B); rank_cha=RB-RA; if rank_cha>0, disp(因?yàn)镽A!=RB,所以此方程組無解.) return end if RA==RB if RA==n disp(因?yàn)镽A=RB=n,所以此方程組有唯一解.) X=zeros(n,1);C=zeros(1,n+1); T=0;D=B;s=max(abs(D(1:

24、n,:))); S=s; for p=1:n-1 [Y,j]=max((abs(B(p:n,p)))./S(p:n)); J=min(j); C=B(p,:);B(p,:)=B(J+p-1,:);B(J+p-1,:)=C; T=S(p);S(p)=S(J+p-1);S(J+p-1)=T; for k=p+1:n

25、 m=B(k,p)/B(p,p); B(k,p:n+1)=B(k,p:n+1)-m*B(p,p:n+1); end end b=B(1:n,n+1); A=B(1:n,1:n);X(n)=b(n)/A(n,n); for q=n-1:-1:1 X(q)=(b(q)-sum(A(q,q+1:n)*X(q+1:n)))/A(q,q); end else disp(請注意:

26、因?yàn)镽A=RB>b=[-1;-1;-1;-1;-1];[RA,RB,n,X]=bilizy(A,b) 因?yàn)镽A=RB=n,所以此方程組有唯一解. RA = 5 RB = 5 n = 5 X = 1.0e-004 * 0.0000 -0.0000 -0.0000 -0.1256 0.1550 >> X=vpa(X,1

27、6) X = 0.00000000008015133713756955 -0.0000000000999191365828799 -0.0000000001763979458318154 -0.00001255606517769044 0.00001549777504860351 2)、Jacobi迭代法 ①、程序代碼: %P:范數(shù)的名稱,P=1,2,inf %error:近似解a的誤差 %maxl:迭代的最大次數(shù) function [a]=Jac

28、obi(A,b,a0,P,error,maxl) [n n]=size(A); a=zeros(n,1); for k=1:maxl for j=1:n a(j)=(b(j)-A(j,[1:j-1,j+1:n]) a0([1:j-1,j+1:n]))/A(j,j); end erra=norm(a-a0,P); a0=a; if(erra

29、=error) disp(請注意:Jacobi迭代次數(shù)已經(jīng)超過最大迭代次數(shù)maxl。) end ②、運(yùn)行結(jié)果: >>X=[53605 58460 62859 66662 68894]; Y=[6026 11197 16954 23492 68894]; B=[X.*X;2*X.*Y;Y.*Y;2*X;2*Y]; %B是系數(shù)矩陣的轉(zhuǎn)置 A=B;%A是系數(shù)矩陣 b=[-1;-1;-1;-1;-1]; a0=[0;0;0;0;0]; >> [a]=Jacobi

30、(A,b,a0,1,0.0001,100) 近似解a分別是: a = -0.0000000003480081375786834 -0.0000000007638525188999265 -0.000000003479009721870011 -0.000007500525036752573 -0.000007257526054518536 3)、Gauss-Seidel迭代法 ①、程序代碼: %P:范數(shù)的名稱,P=1,2,inf;error:近似解x的誤差;max

31、l:迭代的最大次數(shù) function [a]=GS(A,b,a0,P,error,maxl) [n n]=size(A); a=zeros(n,1); for k=1:maxl for j=1:n aa=0; for i=1:n if ij aa=aa+A(j,i)*a0(i); end en

32、d a(j)=(b(i)-aa)/A(j,j); end erra=norm(a-a0,P); a0=a; if (erra=error) disp(請注意:Gauss-Seidel迭代次數(shù)已經(jīng)超過最大迭代次數(shù)maxl.) end

33、 ②、運(yùn)行結(jié)果: >> X=[53605 58460 62859 66662 68894];Y=[6026 11197 16954 23492 68894]; B=[X.*X;2*X.*Y;Y.*Y;2*X;2*Y]; A=B;b=[-1;-1;-1;-1;-1];a0=[0;0;0;0;0]; >> [a]=GS(A,b,a0,inf,0.00001,100) 近似解a分別是: a = 0.0000000001678142444990526 0.000000002024408102712696 0.000000000369101766105352

34、 -0.000009750490974607461 -0.0001554717527254895 2.3 設(shè) (1) 將矩陣A進(jìn)行LU分解A=LU,其中U是上三角矩陣,L是主對角線上的元素都是1的下三角矩陣。 解: >> A=[20 2 3;1 8 1;2 -3 15]; >> [L,U]=lu(A) L = 1.0000 0 0 0.0500 1.0000 0 0.1000 -0.4051 1.0000

35、 U = 20.0000 2.0000 3.0000 0 7.9000 0.8500 0 0 15.0443 (2) 利用上述分解分別求解方程組 并由此求出逆矩陣。 >> b1=[1;0;0];b2=[0;1;0];b3=[0;0;1]; >> X1=U\(L\b1), X2=U\(L\b2),X3=U\(L\b3),inv(A) X1 = 0.0517 -0.0055 -0.0080 解: inv(A) =

36、 0.0517 -0.0164 -0.0093 -0.0055 0.1237 -0.0072 -0.0080 0.0269 0.0665 X3 = -0.0093 -0.0072 0.0665 X2 = -0.0164 0.1237 0.0269 (3) 用LU分解求下列線性方程組的解 (方程組的精確解是) 解: >> A=[4 2 -3 -1 2 1 0 0 0 0; 8 6 -5 -3 6 5 0 1 0 0; 4

37、 2 -2 -1 3 2 -1 0 3 1; 0 -2 1 5 -1 3 -1 1 9 4; -4 2 6 -1 6 7 -3 3 2 3; 8 6 -8 5 7 17 2 6 -3 5; 0 2 -1 3 -4 2 5 3 0 1; 16 10 -11 -9 17 34 2 -1 2 2; 4 6 2 -7 13 9 2 0 12 4; 0 0 -1 8 -3 -24 -8 6 3 -1]; >> b=[5;12;3;2;3;46;13;38;19;-21]; >> [L U]=lu(A)

38、; >> x=U\(L\b) x = 1.0000 -1.0000 0.0000 1.0000 2.0000 0.0000 3.0000 1.0000 -1.0000 2.0000 結(jié)果是: 2.4 (1) 用追趕法求解方程組 (a) , (b) , 解: (a)、程序代碼: %追趕法 function x=zhuiganfa(A,b) n=rank(A); for i=1:n if A(i,i

39、)==0 return; end end d=ones(n,1);a=ones(n-1,1);c=ones(n-1); for i=1:n-1 a(i,1)=A(i+1,i);c(i,1)=A(i,i+1);d(i,1)=A(i,i); end d(n,1)=A(n,n); %解Ly=b的解,解保存在b中 for i=2:n d(i,1)=d(i,1)-(a(i-1,1)/d(i-1,1))*c(i-1,1); b(i,1)=b(i,1)-(a(i-1,1)/d(i-1,1))*b(

40、i-1,1); end %求解Uy=x的解 x(n,1)=b(n,1)/d(n,1); for i=n-1:-1:1 x(i,1)=(b(i,1)-c(i,1)*x(i+1,1))/d(i,1); end 0.1667 0.1667 0.1667 0.1667 0.1667 0.1667 0.1666 0.1668 0.1662 0.1683 0.1607 0.1891

41、 0.0829 0.4793 x = 0.4793 0.0829 0.1891 0.1607 0.1683 0.1662 0.1668 0.1666 0.1667 0.1667 0.1667 0.1667 0.1667 0.1667 0.1667 運(yùn)行結(jié)果: >> clear all >> A=zeros(30,30); A(30,30)=4; b=ones(30,1); b(1,1)=2;b(30,1)=2

42、; for i=1:29 A(i,i+1)=1;A(i+1,i)=1; A(i,i)=4; end >> x=zhuiganfa(A,b) >> A=zeros(20,20); A(19,19)=5;A(20,19)=-2; A(19,20)=-2;A(20,20)=5; b=zeros(20,1); b(1,1)=1; for i=1:18 A(i,i+2)=1;A(i+2,i)=1; A(i,i+1)=-2;A(i+1,i)=-2; A(i,

43、i)=5; end >> [L,U]=lu(A); >> x=vpa(U\(L\b),6) 0.0000504376 0.000106309 0.0000292324 -0.0000143431 -0.0000126677 -0.00000146436 0.00000249126 0.00000131198 -9.33542*10^(-8) -2.99738*10^(-7) x = 0.242121 0.0943914

44、 -0.0218242 -0.0313623 -0.00694256 0.00488693 0.00358612 0.000214823 -0.000784527 -0.000357862 (b)、 (2) 設(shè)計(jì)實(shí)驗(yàn)驗(yàn)證Hilbert矩陣的病態(tài)性,其中 format rat; for n=1:10 A=hilb(n); tiaojianshuA=norm(A,1)*norm(inv(A),1) en

45、d 解: tiaojianshuA = 1 tiaojianshuA = 27 tiaojianshuA = 748 tiaojianshuA = 28375 tiaojianshuA = 943656 tiaojianshuA = 29070279 tiaojianshuA = 985194890 tiaojianshuA = 33872791987 tiaojians

46、huA = 1099651961846 tiaojianshuA = 35353690673203 從上面的結(jié)果可以看出,Hilbert矩陣的條件數(shù)的增長速度驚人,才十階,矩陣的條件數(shù)就為35353690673203,這是一個(gè),可怕的數(shù)字,由它組成的方程組肯定是病態(tài)的。 實(shí)驗(yàn)三 3.1 用Taylor級數(shù)的第項(xiàng)多項(xiàng)式,分別在區(qū)間和上逼近正弦函數(shù),并用計(jì)算機(jī)繪出上面31個(gè)函數(shù)的圖形。 解: 一、程序代碼: syms x; y=sin(x); ezplot(y,[0,pi]); grid on;

47、 axis([0,pi,0,1.5]); hold on; for m=1:2:30 p=taylor(y,x,m); ezplot(p,[0,pi]); axis([0,pi,0,1.5]); end 二、函數(shù)圖形: 3.2 已知四個(gè)點(diǎn)A,B,C,D的具體位置A(0, 0),B(0, 3), C(8, 1), D(10, 5),求兩個(gè)點(diǎn)H1,H2的具體位置,使AH1+BH1+H1H2+H2C+H2D為最短。 解:function [ d ] =ju_li(a,b) d=sqrt(

48、sum((a-b).^2)); return end syms x1 y1 x2 y2; A=[0 0];B=[0 3];C=[8 1];D=[10 5];H1=[x1 y1];H2=[x2 y2]; sum_d=ju_li(A,H1)+ju_li(B,H1)+ju_li(H1,H2)+ju_li(C,H2)+ju_li(D,H2) [x1 y1 x2 y2]= solve(diff(sum_d,x1),diff(sum_d,y1),diff(sum_d,x2),diff(sum_d,y2),x1,y1,x2,y2) su

49、m_d = 3*x1^2 - 2*x1*x2 + 3*x2^2 - 36*x2 + 3*y1^2 - 2*y1*y2 - 6*y1 + 3*y2^2 - 12*y2 + 199 x1 = 9/4 y1 = 27/4 x2 = 15/8 y2 = 21/8 實(shí)驗(yàn)四 4.1 方程的根實(shí)際上是兩個(gè)函數(shù)的交點(diǎn)的橫坐標(biāo),用計(jì)算機(jī)繪出兩個(gè)函數(shù)在區(qū)間的圖形,觀察圖形,分析它們的交點(diǎn)分布規(guī)律及特點(diǎn),試寫出方程的全部實(shí)根所在的隔根區(qū)間;并用二分法求出每一個(gè)根的近似值。 解: 一、畫函數(shù)圖像如下: >> syms

50、x y1 y2; >> y1=sin(x);y2=1/x; >> ezplot(y1,[-6,6]);grid on; hold on; ezplot(y2,[-6,6]) >> gtext(y1=sin(x)) >> gtext(y2=1/x) 2、 用二分法求近似解: 通過上面的函數(shù)圖像可以知道,y1=sin(x)與y2=1/x在區(qū)間[-6,6]上共有四個(gè)交點(diǎn),即是函數(shù)y=xsin(x)在區(qū)間[-6,6]上共有四個(gè)零點(diǎn),

51、調(diào)用matlab中的函數(shù)【>> ginput(4) 】 可以大概得到四個(gè)根的大概區(qū)間為:[-4,-2]、[-2,-0.5]、[0.5,2]、[2,4],通過二分法求四個(gè)點(diǎn)的近似值的過程如下: %用二分法求非線性方程f(x)=0的根;fun為函數(shù)的表達(dá)式 %a,b為左右端點(diǎn),eps為精度,x為近似根,k為二分次數(shù) function [x,k]=eff(fun,a,b,eps) if nargin<4 eps=1e-5;%如果輸入自變量數(shù)目<4,默認(rèn)eps=1e-5 end fa=feval(fun,a);fb=feval(fun,b); if f

52、a*fb>0 disp([a,b]不是有根區(qū)間,請重新調(diào)整); return; end k=0; while abs(b-a)/2>eps x=(a+b)/2;fx=feval(fun,x); if fx*fa<0 b=x;fb=fx; else a=x;fa=fx; end k=k+1; end x=(a+b)/2; ①、程序代碼: ②、運(yùn)行結(jié)果: >

53、> fun=@(x)(x*sin(x)-1) ;eps=1e-5; >> [x1,k1]=eff(fun,-4,-2,eps),[x2,k2]=eff(fun,-2,-0.5,eps), [x3,k3]=eff(fun,0.5,2,eps),[x4,k4]=eff(fun,2,4,eps) x1 = -2.7726 k1 = 17 x4 = 2.7726 k4 = 17 x3 = 1.1142 k3 = 17 x2 =

54、-1.1142 k2 = 17 4.2 設(shè)方程有3個(gè)實(shí)根 盡可能多地采用下面六種不同計(jì)算格式,求的根或或的近似值,并觀察相應(yīng)格式的收斂性。 (1) (2) (3) (4) (5) (6) 解: >> syms x fun; >> fun=x^3-3*x-1; >> ezplot(fun,[-2.5,2.5]),grid 1、 函數(shù)圖像: 得到該函數(shù)的零點(diǎn)大概在:-1.5、-0.3、1.8的附近。 2、 迭代格式: 程序代

55、碼: %迭代法 %用迭代法求非線性方程f(x)=0的根,fun為函數(shù)ψ(x)的表達(dá)式 %x0為初值,eps為精度(默認(rèn)1e-5),N為最大迭代次數(shù)(默認(rèn)500),x為近似根 function [x,k]=ddf(fun,x0,eps,Nmax) if nargin<4 Nmax=500; end if nargin<3 eps=1e-5; end x=x0;x0=x+2*eps;k=0; while abs(x0-x)>eps&k

56、 warning(已迭代次數(shù)上限); end (1) >> fun=inline((3*x+1)/x^2); >> [x1,k1]=ddf(fun,-1.5,1e-5),[x2,k2]=ddf(fun,-0.3,1e-5),[x3,k3]=ddf(fun,1.8,1e-5), -1.5321 k1 = 28 x2 = NaN k2 = 31 x3 = NaN k3 = 4

57、8 x2 = NaN k2 = 31 x3 = NaN k3 = 48 x1 = -1.5321 k1 = 28 (2) >> fun=inline((x^3-1)/3); >> [x1,k1]=ddf(fun,-1.5,1e-5),[x2,k2]=ddf(fun,-0.3,1e-5),[x3,k3]=ddf(fun,1.8,1e-5), x2 = -0.3473 k2 = 5 x3

58、= -0.3473 k3 = 9 x1 = -0.3473 k1 = 12 (3) >> fun=inline((3*x+1)^(1/3)); >> [x1,k1]=ddf(fun,-1.5,1e-5),[x2,k2]=ddf(fun,-0.3,1e-5),[x3,k3]=ddf(fun,1.8,1e-5), x3 = 1.8794 k3 = 8 x2 = 1.8794 k2 = 12

59、 x1 = 1.8794 + 0.0000i k1 = 12 >> fun=inline(1/(x^2-3)); >> [x1,k1]=ddf(fun,-1.5,1e-5),[x2,k2]=ddf(fun,-0.3,1e-5),[x3,k3]=ddf(fun,1.8,1e-5), (4) x2 = -0.3473 k2 = 5 x3 = -0.3473 k3 = 7 x1 = -0.3473 k1 = 8

60、 (5) >> fun=inline((3+1/x)^(1/3)); >> [x1,k1]=ddf(fun,-1.5,1e-5),[x2,k2]=ddf(fun,-0.3,1e-5),[x3,k3]=ddf(fun,1.8,1e-5), x1 = 1.5396 k1 = 6 x2 = 1.5396 + 0.0000i k2 = 7 x3 = 1.5396 k3 = 5 (6) >> fun=inline(x

61、-(1/3)*(x^3-3*x-1)/(x^2-1)); >> [x1,k1]=ddf(fun,-1.5,1e-5),[x2,k2]=ddf(fun,-0.3,1e-5),[x3,k3]=ddf(fun,1.8,1e-5), x1 = -1.5321 k1 = 3 x2 = -0.3473 k2 = 3 x3 = 1.8794 k3 = 4 4.3 一個(gè)10次多項(xiàng)式 的系數(shù)為 [1, -55, 1320, –18150

62、, 157773, –902055, 3416930, –8409500, 12753576, –10628640, 6328800] 解: (1) 用多項(xiàng)式的求根指令roots求出方程的10個(gè)根; >> clear all >> p=[1 -55 1320 -18150 157773 -902055 3416930 -8409500 12753576 -10628640 6328800]; >> roots(p) ans = 10.6051 + 1.0127i 10.6051 - 1.0127i 8.5850 + 2.7898i

63、 8.5850 - 2.7898i 5.5000 + 3.5058i 5.5000 - 3.5058i 2.4150 + 2.7898i 2.4150 - 2.7898i 0.3949 + 1.0127i 0.3949 - 1.0127i (3) 將方程左邊多項(xiàng)式中的9次項(xiàng)的系數(shù)–55改為–56,得到一個(gè)新的10次方程;求解新的方程,觀察根的變化是否很顯著,說明了什么? >> p=[1 -56 1320 -18150 157773 -902055 3416930 -8409

64、500 12753576 -10628640 6328800]; >> roots(p) ans = 21.7335 7.3501 + 7.7973i 7.3501 - 7.7973i 4.3589 + 3.3285i 4.3589 - 3.3285i 5.1831 2.4378 + 2.7974i 2.4378 - 2.7974i 0.3949 + 1.0127i 0.3949 - 1.0127i

65、 實(shí)驗(yàn)五 5.1 給出的函數(shù)表如下: 解: ⑴、分別用Lagrange插值和Newton插值求的近似值; ①、Lagrange插值: 程序代碼: %Lagrange插值法 function f=Lagrange(x,y,x0) syms t; if length(x)==length(y) n=length(x); else disp(x和y的維數(shù)不相等!); return; end

66、 %驗(yàn)錯(cuò) f=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)); %計(jì)算Lagrange基函數(shù) end f=f+l; %計(jì)算Lagrange插值函數(shù) simplify(f); %化簡 end if nargin==3 f=vpa(subs(f,t,x0),10); %計(jì)算插值點(diǎn)的函數(shù)值 else f=collect

展開閱讀全文
溫馨提示:
1: 本站所有資源如無特殊說明,都需要本地電腦安裝OFFICE2007和PDF閱讀器。圖紙軟件為CAD,CAXA,PROE,UG,SolidWorks等.壓縮文件請下載最新的WinRAR軟件解壓。
2: 本站的文檔不包含任何第三方提供的附件圖紙等,如果需要附件,請聯(lián)系上傳者。文件的所有權(quán)益歸上傳用戶所有。
3.本站RAR壓縮包中若帶圖紙,網(wǎng)頁內(nèi)容里面會(huì)有圖紙預(yù)覽,若沒有圖紙預(yù)覽就沒有圖紙。
4. 未經(jīng)權(quán)益所有人同意不得將文件中的內(nèi)容挪作商業(yè)或盈利用途。
5. 裝配圖網(wǎng)僅提供信息存儲(chǔ)空間,僅對用戶上傳內(nèi)容的表現(xiàn)方式做保護(hù)處理,對用戶上傳分享的文檔內(nèi)容本身不做任何修改或編輯,并不能對任何下載內(nèi)容負(fù)責(zé)。
6. 下載文件中如有侵權(quán)或不適當(dāng)內(nèi)容,請與我們聯(lián)系,我們立即糾正。
7. 本站不保證下載資源的準(zhǔn)確性、安全性和完整性, 同時(shí)也不承擔(dān)用戶因使用這些下載資源對自己和他人造成任何形式的傷害或損失。

相關(guān)資源

更多
正為您匹配相似的精品文檔
關(guān)于我們 - 網(wǎng)站聲明 - 網(wǎng)站地圖 - 資源地圖 - 友情鏈接 - 網(wǎng)站客服 - 聯(lián)系我們

copyright@ 2023-2025  zhuangpeitu.com 裝配圖網(wǎng)版權(quán)所有   聯(lián)系電話:18123376007

備案號:ICP2024067431-1 川公網(wǎng)安備51140202000466號


本站為文檔C2C交易模式,即用戶上傳的文檔直接被用戶下載,本站只是中間服務(wù)平臺,本站所有文檔下載所得的收益歸上傳人(含作者)所有。裝配圖網(wǎng)僅提供信息存儲(chǔ)空間,僅對用戶上傳內(nèi)容的表現(xiàn)方式做保護(hù)處理,對上載內(nèi)容本身不做任何修改或編輯。若文檔所含內(nèi)容侵犯了您的版權(quán)或隱私,請立即通知裝配圖網(wǎng),我們立即給予刪除!

五月丁香婷婷狠狠色,亚洲日韩欧美精品久久久不卡,欧美日韩国产黄片三级,手机在线观看成人国产亚洲