1、數(shù)值分析上機實驗1、已知a與b_ 12.38412215237-1.0610742.11523719.141823-3.125432-1.061074-3.12543215.5679141.112336-1.012345323848a =-0135842.1897362.0314540.7187191.5638491.8367421.742382-0.784165-1.0567813.0678131.1123480.336993-2.0317433.123124-1.010103b = 2.187436933.992318-25.1734171.11233641135840.718719-1.

2、0123452.1897361.5638493.1238482.0314541.83674227.1084374.101011-3.7418564.10101119.8979180.431637-3.7418560.4316379.7893652.1010233111223-0.103458-0.718282.121314-1.103456-0.0375851.7843170.2384170.846716951.784317-86.6123431.7423823.067813-2.031743-0.7841651.1123483.123124-1.0567810.336993-l0101032

3、.101023-0.71828-0.037585-3.1112232.1213141.784317-0.103458-1.1034560.23841714.71384653.123789-2.2134743.12378930.7193344.446782-2.2134744.44678240.00001 _1.11012304.71934556784392(2)用超松弛法求解bx二b (取松弛因子<o=1.4, x(0>=0,迭代9次)。(3)用列主元素消去法求解bx=bo解：(2)用超松弛法求解bx=b：(一) 、超松弛法的理論依據(jù)為：對于j方法或gs方法有時收斂速度較慢，松弛法

4、針對加 快收斂速度提出來的。其基本思想是在gs方法已求111 x(m),x(m-1)的基礎(chǔ)上，經(jīng)過重新組合 得到新序列，而此新序列收斂速度加快。組合時川到松弛因了 co,當(dāng)co>l時稱為超松弛法,當(dāng)3>1時，可以加快收斂速度。若a為對稱正定陣，則當(dāng)松弛因子3滿足032時，松 弛法收斂。本題中(0=1.4。超松弛法為：xi仲)=(1(0)治伸")+co( < 1?必嚴(yán))+ £ xj+gi)y=l；=/+1(二) 、計算程序(c語言)#include"stdio.h”#includcnmath.h"main()intdouble t=0,s

5、=0,w=1.4,g9,b99,x9=0,0,0,0,0,0,0,0,0;給出初始條件和松弛因子double b9=2.1874369,33.992318,-25.173417,0.84671695,1.784317,-86.612343,1.1101230,4.719345,- 5.6784392,a910= 12.38412 ,2.115237 廠 1.061704,1.112336 ,0.113584,0.718719 ,1.742382 ,3.067813,2031743, 2.115237 ,19.141823,-3.125432,-1.012345,2.189736 ,1.5638

6、49 ,-0.784165 ,1.112348 ,3.1231241.061074,3125432,15.567914,3.123848 ,2.031454 ,1.836742 ,-1.056781 ,0.336993 ,-1.010103 ，112336 ,-1.012345,3.123848 ,27.108437,4.101011 ,-3.741856,2j01023 ,-0.71828 ,-0.037585-0.113584,2.189736 ,2.031454 ,4.010111 ,19.897918,0.431637 ,-3.111223 ,2.121314 ,1.784317 0.

7、718719 ,1.563849 ,1.836742 ,-3.741856,0.431637 ,9.789365 ,-0.103458 ,-1.103456,0.238417 1.742382 ,-0.784165,-1.056781,2.101023 ,-3.111223,-0.10345& 14.7138465,3.123789 ,-2.21347 4,3.067813 ,112348 ,0.336993 ,-0.71828 ,2.121314 ,-1.103456,323789 ,30.719334,4.446782 ，2.031743,3.123124 ,-1.010103,-

8、0.037585,1.784317 ,0.238417 ,2213474 ,4.446782 ,40.00001 ；for(i=0;i<=8;i+)gi=bi/aii;歹bi/如，5=0*/bii=o;fbr(i=0;i<=8;i+)for(j=0;j<=8;j+) if(i!=j) bij=-1.0*aij/aii;/*算出心臨 ih薩/for(m=l;m<=20;m+)/*外層發(fā)循環(huán)控制迭代次數(shù)，本題要求迭代勺/*次數(shù)為9次，但迭代9次精度不高，所以就多迭代了幾次*/fdr(i=0;i<=8;i+)for(j=0;j<=i-l;j4-+)一1t=t+bi

9、j*xuj；/*算出工%xj伽 */冃for(j=汁 l;jv=8;j+)s=s+bij*xj；嚴(yán)算出 工“xt” */；=/+!xi=(l -w)*xi+w*(t+s+gi);/*算出每次迭代后的值可t=s=o.o;嚴(yán)給s, t清零,若不清零會影響后面的計算*/printf("the rusult is:n");for(i=0;i<=8;i4-+) printf(u%lfnn,xi);/*打印出計算結(jié)果 */getch();(三) 、計算結(jié)果打?。簁l0=0.489997 x13=2.544899 x23=0.362510x3=-1.420572x43=0.4195

10、06 k5j=-9.957563x63=0.464643 xt7=-0.429914 x8=-0.193821(四) 、問題討論：松弛方法的矩陣形式為：x(m)=(e- 3 l) 1 (1 3 )e4- 3 r)x(m- 1)+(e- 3 l) 1 3 g若記l3=(e3l)1(13)e+3r),sor收斂的充要條件是s(l3)<1.且若a為對稱正定陣, 則當(dāng)松弛因了滿足0<3<2時，sor方法收斂。此題屮矩陣b是對稱正定陣，且是三對角的，所以選擇合適的松馳因子3,收斂速度是很 快的。其最佳松弛因子可用下式計算：3 二2/1 +j1(£()彳在迭代9次的時候得出的結(jié)

11、果與準(zhǔn)確值有一定的誤差，但當(dāng)?shù)螖?shù)增加到一定次數(shù)，其結(jié) 果的準(zhǔn)確度便相當(dāng)?shù)母摺?3)、用列主元素消去法求解bx=b(一) 、理論依據(jù)：其基木思想是選取絕對值盡量大的元素作為主元素，進(jìn)行行與列的交換，再進(jìn)行回代,求出方程的解。將方陣a和向量b寫成c= (ab)。將c的第1列屮第1行的元索與其下面的此列的元索逐 一進(jìn)行比較，找到最大的元素勺將第j行的元素少第1行的元素進(jìn)行交換，然后通過行 變換，將笫1列中笫2到笫n個元素都消成0。將變換后的矩陣c的笫二列中第二行的元 素與其下而的此列的元素逐一進(jìn)行比較，找到最大的元素c；?,將第k行的元素與第2行的 元素進(jìn)行交換，然后通過行變換，將第2列中第3到

12、第n個元素都消成0。以此方法將矩陣 的左下部分全都消成0。(二) 、計算程序：#includc<stdlib.h>#include<stdio.h>main()int i,j,k;float abc,a99=12.38412,2.115237, 1.061074,l.112336, 0.113584,0.718719,1.742382,3.067813,2031743,2.115237,19.141823,3125432, 1.012345,2.189736,1.563849,-0.784165,1.112348,3.123124,j.061074,-3.125432,

13、15.567914,3.123848,2.031454,1.836742, 1.056781,0.336993, 1.010103,1.112336,-1.012345,3.123848,27.108437,4.101011,-3.741856,2.101023,-0.71828,-0.037585,-0.113584,2.189736,2.031454,4.101011,19.897918,0.431637,-3.111223,2.121314,1.784317, 0.71819,1.563849,1.836742, 3.741856,0.431637,9.789365, 0.103458,

14、l.103456,0.238417, 1.742382,-0.784165,-1.056781,2.101023,-3.111223,-0.103458,14.7138465,3 123789,-2.2123474 ，3.067813,1.112348,0.336993, 0.71828,2.121314,1.103456,3.123789,30.719334,4.446782,-2.031743,3.123124,-1.010103,-0.037585,1.784317,0.238417,-2.213474,4.446782,40.00001 ,b9=2.1874369,33.99231&a

15、mp; -25.173417,0.84671695,1.784317,-86.612343,1.1101230,4.719345,-5.6784392;for(k=0;k<9;k+)abc=akk;for(j=0;j<9;j+)aku=akj/abc;bk=bk/abc;for(i=0;i<9;i+)abc=aik;for0=0;j<9;j+)if(i!=k) aiu=aiu-abc*aku;if(i!=k) bi=bi-abc*bk;for(i=0;i<9;i+)printf(”x%d=%fn”,i,bi);(三) 、計算結(jié)果打印:x0=0.491863xtl3

16、=2.546492x 21=0.361278 x3j=-l.421235xl4j=0.419984xt5=-9.958960xl6j=0.462853xt7j=-0.430627xt81=-0.193904(四) 、問題討論:a. 由于選主元，使|lij|最小，這樣去乘方程的每一系數(shù)時，系數(shù)的舍入誤差不至擴人， 并防止溢出與停機。b. 由于選主元，回代時作除數(shù)主元aii(i)的絕對值也是最大,這樣擴大的誤差也是最小.c. 若detaho貝j選主元后的主元aii(i)工0,這是因為若aii(i)=0則必冇aji(i)=o (j>i)這樣按 行列式的laplace展開式就有deta=o而矛庇

17、因此用主元消去法中進(jìn)行下去，不止中斷停機. 同吋也由于選主元,aii(i)接近零的概率減少.運行結(jié)果基本與準(zhǔn)確值無界，因為這種算法的無誤差的算法。上機實驗2用三次樣條插值求函數(shù)值一、解題條件和理論依據(jù)及算法1、解題條件己知十點兩數(shù)值及起點和終點的導(dǎo)數(shù)值，要求用三次樣條插值求f (4. 563)及f' (4. 563) 的近似值.x12i5f（x）00.693147181.09861231.38629441.6094378x678910f（x）1.79175951.94591012. 0794452.19722462.3025851f（x）f（x） = lf（x）二0. 12、理論依據(jù)及

18、算法為使問題一般化，本題的三次樣條求解使用了對于任意分化的三彎矩插值法。若記各節(jié) 點間距叫一嚴(yán)xj-xj_ 21,2,m,根據(jù)理論分析，三次樣條插值歯數(shù)為：）3疋 i （xf -x）x -xs -節(jié)（十心w于）（寸）其中兀丘(兀1宀)。而式中的諸可由如卜-的矩陣方程來確定:_ 21mq2右m£ ai2右 ai2_dn_這里6 （兒-6 （uh + h j h.-yn _ -（yn hzhnhi-2 = 1 “ =7y /» +化d。）!3 71 ） = 6/（%._!，兀 j,xj+i）一 y.v-1 ）1（j = 1,2,，n -1）對于上述矩陣方程，由于其系數(shù)矩陣為一三

19、對角陣，故我們用追趕法（thomas算法） 來求解。-般地，對矩陣方程5xb2兀2d.an-bn-lanbn _dn_有如下算法：= x，厲=d】,/小=予丄，0小=乞+| -厶+心< c+i =心+】-厶+0, i =/? - 1% =邑-, x. = s, - c,'v/ + 1 , i = n -”0”， 隊由此可求ii1. “j ,故'(x)中各參數(shù)均已求出故可求出/(4.563) = s(4.563),廣(4.563) = s'(4.563)o二、程序清單：#include<stdio.h>#include<math.h>#def

20、ine n 9void main()double xn+l=l,2,3,4,5,6,7,8,9,10,yn+l=0,0.69314718,1.0986123,1.3862944,1.6094378,1.7917595,1.9459101,2.079445,2.1 972246,2.3025851,hn+l,dn+l,an+l,cn+l,bn+l=2,2,2,2,2,2,2,2,2,2,sn+l,tn+l,ln+l,mn+l, f,fl ；int i;for(i=l;i<=n;i卄)/*使用對任意分化的三彎矩插值法*/hi-l=xi-xi-l;d0=6/h0*(yl-y0)/h0-l);d

21、n=6/hn-l*(0.1 -(yn-yn-l)/hn-l);fbr(i=l;i<=n-l;i+)di=6/(hil+hi)*(yi+lyi)/hi(yiyil)/hil); ai=hi-l/(hi-l+hi);ci=l-ai;fbr(i=o;i<=n-l;i+)li+ l=ai+l/si; si+l=bi+l-li+l*ci; ti+l=di+l-li+l*ti;mn=tn/sn;for(i=n-l ;i>=o;i)mi=(ti-ci*mi+l)/si;fm3*pow(x4-4.563),3)/6/h3/*求算 4.563 這點的函數(shù)值*/+m4*pow(4.563-x3)

22、,3)/6/h3+(y3-m 3*h3*h3/6)*(x4-4.563)/h3+(y4m4*h *h3/6)*(4.563-x )/h ;fl=-3*m3*pow(x4-4.563),2)/6/h3/*求算 4.563 這點的一階導(dǎo)數(shù)值*/+3*m4*pow(4.563x 3),2)/6/h3(y3-m3*h3*h3/6)/h3 +(y4-m4*h3*h3/6)/h3;printf(uf(4.563)=%lf f(4.563)=%lfn”,f,fl);三、運行結(jié)果:四、問題討論:樣條插值效果比lagrange插值好，沒有runge現(xiàn)象.計算比高次插值函數(shù)簡單，效果比 高次插值函數(shù)好。且光滑性比

23、低次分段插值函數(shù)好。它集小了高次插值函數(shù)，低次插值函 數(shù)，低次分段插值函數(shù)的優(yōu)點。對各個插值點(1, 2, 3. 10)分別運行以上程序，得到的結(jié)果與所給的插值點函數(shù)值 除了 1, 10兩點的值有較小的差異外，其余點的值都相當(dāng)接近。這也說明了樣條插值算法 的穩(wěn)定性。上機實驗3 newton法求方程4、用newton法求方程：x7-28x4+14=0在(0. 1, 1. 9)中的近似根(初始近似值取為區(qū)間端點, 迭代6次或誤差小于0. 00001.)一、解題條件和理論依據(jù)及算法1、解題條件用牛頓法求方程£ - 28x4 +14 = 0在(0. 1 19)中的近似根(初始近似值取為區(qū)間端

24、點)，迭代6次或誤差小于0. 00001.2、解題依據(jù)及算法設(shè)函數(shù)在有限區(qū)間a,b上二階導(dǎo)數(shù)存在，且滿足條件：i f(a)f(b)<0ii f：(x)在區(qū)間a,b_上不變號.iii f (x) hoiv |f (c) |/b-aw|f (c) | 其中 c 是 a,b 中使 min|f' (a), f' (b)達(dá)到的一個, 則對任意吋近似值x0a, b, newton迭代過程為xk+i二 (xj 二x-f (xj /f' (xk), k二 1,2, 3 令= x7-28x4 +14,/(0.1) > 0 j(1.9) < 0 廣= 7x6-1 12疋=

25、lxx3-16)<0 fx) = 42x5 -336x2 = 42x2(x3 -8)<0 廠(1.9)/(1.9)>0故以19為起點二、計算程序(c語言)#include "stdio.h"#include "math.h”main()double x 1=1.9,eps=0.00001 ,x0;/* 輸入初值 xl*/dox0=xl;x 1 =x0-(pow(x0,7)-28*pow(x0,4)+14)/(7*pow(x0,6)-l 12*pow(x0,3)；/*利用算法公式迭代*/*結(jié)束條件表述*/while(fabs(x 1 -xo)>

26、;=eps);printff'x l=%5fnh,x 1);三、運行結(jié)果ca “e:down譴籌學(xué)譴籌學(xué)debugcpp3.exe"<1 =0.845497press any key to continue四、問題討論：1、誤差分析：牛頓迭代法是平方收斂的，該題經(jīng)過六次迭代誤差就小于0.00001,和雙點弦割法、單點弦割法比較，牛頓法較精確。.2、在編程的過程中，我們可以用所給的迭代數(shù)中斷程序（如上面程序）， 也可根據(jù)所給誤差小0.00001來控制迭代次數(shù)，即|xk+1-xk|<0.00001從而中斷程序。上機實驗4用romberg算法計算積分值5、用 rombe

27、rg 算法求 j 3' xa（5x + 7）sin x2dx （允許誤差 £ =0. 000 01）o 解：（一）、理論依據(jù)：romberg算法的計算步驟如一 f:（1）先求出按梯形公式所得的積分值：刁=.f（q）+ .f（b）（2 ）把區(qū)間2等分，求出兩個小梯形面積z和，記為7；，即人二乎/+/的+ 2/（晉）這樣由外推法可得可°）,羅）=4葉1）葉0）4-1（3）把區(qū)間再等分（即f等分），得復(fù)化梯形公式7；，由7；與7；外推可得4?。?）_卩（1）4?t _卩（°）砂=十，酹=十產(chǎn)，如此，若已算出2等分的復(fù)化梯形公式畛,則由richardson外推法，

28、構(gòu)造新序列m二 1,2,，1,k二1, 2,1 -m+l,最后求得昭）。（4） 7；彩昭）或17；-即 "就停止計算，否則回到（3）,計算7；（/+,）, 一般可用如下算法：竿）=¥/+/的2/-1«0=i 冋+軒 £ m+-1） #t需二,加= 1,2, ，/, £ = 1,2,，/一加+ 1（2）7；一7 5）於葉2）/當(dāng)?shù)V）0+門a就停機(二) 、計算程序：#include <math.h>float f(float x)float f=0.0;f=pow(3.0,x)*povv(x,l.4)*(5*x+7)*sin(x*x)

29、;return (f);main()int i=l ,j,k,n=12;float t 12,a= 1.0,b=3.0,s=0.0;t0=0.5*(b-a)*(f(a)+f(b);foro=l;j<n-l;j+)for(k= 1 ;k<=pow(2,j 1 );k+) s+=fi(a+(2*k-l)*(b-a)/pow(2,j);tj=0.5*(tj-1 +(b-a)*s/pow(2,j-1);s=0.0; t 11=(4*t1 -t0)/(float)3;fbr(;fabs(t 11 -t0)>0.00001 ;i+)t0=t11;for(j=l ;j<n-1 i；j

30、+) tj =(pow(4,i)*tj+l -tj)/(pow(4,i)-l);tll=(pow(4,i+ 1)*t1 -t0)/(pow(4,i+1)-1)；(三) 、計算結(jié)果打印：the result is:440.536011（四）、問題討論：程序中由count計算循環(huán)次數(shù)，經(jīng)過5次循壞得到計算結(jié)果。romberg算法的優(yōu)點是：1）把積分化為代數(shù)運算，而實際上只需求t1 （1）,以后用遞推可得。2）算法簡單且收斂速度快，一般4或5次即能達(dá)到要求。3）節(jié)省存儲量，算出的幾可存入7；。romberg算法的缺點是:1）對函數(shù)的光滑性要求較高。2）計算新分點的值時，這些數(shù)值的個數(shù)成倍增加6、用定

31、步長四階runge-kutta求解dyx / dt = dy2 !dt =兒dyjdt = 1000 1000)，2 100兒(0)= 0山0) = 0兒(0) = 0h=0.0005,打印 yi(0025),yi(0045),yi(0085),yi(0l),(i=l,2,3)解：(一)、理論依據(jù)：四階runge-kutta方法：h兒+i =兒 + 三 & + % + 2k 3 + kjok = fn< 心二于億+爭，兒+*他)k3 =/(zh + 2 九兒 + -加2)心二/億+力，兒+風(fēng)3)(二) 、計算程序(c語言)：#include <stdio.h>/*建立三個函數(shù)的關(guān)系/*賦初值*/#include <math.h>f(double y4,double k4)*/k1=1;k2=y ；k3 =1000-1000*y2-100*y3;main()int