2023年數(shù)值分析實驗報告新編

河北工業(yè)大學(xué)《數(shù)值分析》課程實驗實驗報告專業(yè):計算機科學(xué)與技術(shù)班級:計1305班姓名:岳帥鵬學(xué)號：133007完畢日期:2023.1.11實驗三一、實驗規(guī)定分別用復(fù)化梯形公式和復(fù)化辛卜生公式計算f(x)=sin(x)/x的積分，并與準確值比較判斷精度。二、實驗?zāi)康耐ㄟ^實際計算來真實體會復(fù)合梯形公式和復(fù)合辛普森公式兩種數(shù)值求積方法的精確度。三、算法描述1、復(fù)合梯形公式環(huán)節(jié)一：擬定函數(shù)的區(qū)間為【a,b】，并擬定其節(jié)點的個數(shù)，算出每個社區(qū)間的步長h;環(huán)節(jié)二：根據(jù)擬定的h,代入復(fù)合梯形公式，算出其相應(yīng)的結(jié)果。從2、復(fù)合辛普森公式a環(huán)節(jié)一:擬定函數(shù)的區(qū)間為【a,b】，并擬定其節(jié)點的個數(shù)，算出每個步長h;環(huán)節(jié)二:根據(jù)擬定的h,代入復(fù)化辛普生公式,算出相應(yīng)的積分結(jié)果。四、源程序#include<stdio.h>#inc1ude<math.h>voidmain(){?intm,n,k,i;f1oatTn,Sn,dl,d2,a=O,b=0,c=0;用oatxll0OOJ,y[I000],z[1(X)0J,w[l000];Primf("將復(fù)化梯形區(qū)間劃分:");。scanf("％d"，&n);?for(i=1；i<n;i++)(。?x[i]=f1oat(1.0/n*i)；gy[i]=(f1oat)(sin(x[i])/x[ij);a+y[i]；。)*Tn=(float)((1,0+sin(1.0)/1.0+2*a)/2.0/n);printf("復(fù)化梯形輸出：%f\n”,Tn)；oprimf("將復(fù)化辛卜生區(qū)間劃分為:")；scanf("%d",&m);)。for(k=1;k<m；k++)(z[k]=f1oat(1.0/m*k);??z[k+l]=float(1.0/m*(k+1));m[k+1/2]=float(z[k]+z[k+l])/2;8w[k]=(f1oat)(sin(z[k])/z[k]);ow[k+1/2]=(float)(sin(z[k+l/2])/((z[k]+z[k+1J)/2));。b=b+w[k];c=c+w[k+1/2];Sn=(float)((l.0+2*b+4*c+sin(I.0)/l.0)/6.0/m);叩rintf("復(fù)化辛卜生輸出：%f\n",Sn);0dl=(lloat)(fabs(0.9460831-Tn));8d2=(float)(fabs(0.9460831-Sn));printf("復(fù)化梯形誤差:%f\n".d1);-printf("復(fù)化辛卜生誤差：%f\n",d2);，if(d1>d2)printf("復(fù)化梯形求法精度低于復(fù)化辛卜生求法\n”)；。elseif(dl==d2)g。printf("亞化梯形求法精度等于兔化辛卜生求法\n");e1se(d1<d2)g;printf("復(fù)化梯形求法精度高于復(fù)化辛卜生求法\n")；))五、運營結(jié)果|MC:\Progra>Files\IicrosoftVisualStudio\ByProjects\3\Debug\3.exe*置復(fù)化梯彭區(qū)間劃分復(fù)化梯形輸出：0.945832西復(fù)化〒4生區(qū)間劃分為復(fù)化老卜生物出：0.741290復(fù)化梯形誤群0.000251復(fù)位辛卜生謨差：0.204793復(fù)化稔形求法精度高于復(fù)化辛卜生求法Pressanykeytocontinue.六、對算法的理解與感悟€實驗四一、實驗規(guī)定用改善歐拉方法解初值問題y，=x+y；y(0)=1。0<x<L取步長h=0.1計算，并與準確值y=?x-l-2e'相比較。二、實驗?zāi)康?)熟悉求解常微分方程初值問題的有關(guān)方法和理論，重要是改善歐拉方法。2)會編制上述方法的計算程序。3)針對實習(xí)題編制程序，并上機計算其所需要的結(jié)果。4)通過對各種求解方法的計算實習(xí)，體會各種解法的功能、優(yōu)缺陷及合用場合,會選取適當?shù)那蠼夥椒ㄈ?、算法描述解一階常微分方程初值問題y'=f(x,y)a<x<b.y(xo)=Yo將區(qū)間[a,b]作n等分,取步長h="衛(wèi)n歐拉公式為Yi+i=yi+hf(Xi，yJ梯形公式為yi+i=yi+^[f(x,,yl)+f(xi+1,yi+1)]改進歐拉法采用公式y(tǒng)i+1=yi+hf(xi,yi)Yi>i=yi+^(nxi,yi)+f(xM,yi.l)]或者yp=Yi+hf(XQi)yc=y,+hf(x+li,yp)=5仇+九)四、源程序#inc1ude<stdio.h>#includc<math.h>#defineM_E2.doubleyp(doub1ex,doubIey,doubleh)；doubleyc(douhlex,doubley,doubleh);doub1eyx(doub1ex)：intmain(){doub1exn=0：doubleyn=l,yxn;doub1eh=0.1;inti;printf("/**********實驗四:\n用改善歐拉方法解初值問題y'=x+y；y(0)=l<.0<x<1,取步長h=0.1計算,\n并與準確值y=-x-l-2ex相比較\n\n");printf("\t\txn\t\tyn\t\ty(xn)\n°)；for(i=();i<IO；i++){?yn=0.5*(yp(xn,yn,h)+yc(xn4-h,yn,h));?xn=xn+h;?yxn=yx(xn);printf("%20f%20f%2Of\n",xn,yn,yxn);)rcturn0;)doubleyp(doublexjoub1ey,doub1eh){returny+h*(x+y);doubleyc(doublex.doubley.doub1eh)returny+h*(x+h+y);)doubleyx(doublex)|return-x-1.0-2.0*exp(x);}五、運營結(jié)果藥WDOCUIEJnSAIDSETTINGS\AD?INISTRATOR\桌面,新建文件夾,實瞼4\Debug\實???P「/XXXXMXXXMX頭四：用改進歐拉方法解初值問題y，=x+y；y<0>=lo0<x〈l,取步長h=0.1計算,棄與程確值y=-xT-2ex相比較xnyny<xn>0.1000001.110000-3.3103420.2000001.241000-3.6428060.3000001.395100-3.9997180.4000001.574610-4.3836490.5000001.782071-4.7974430.6000002.020278-5.2442380.7000002.292306-5.7275050.8000002.601537-6.2510820.9000002.951690-6.8192061.0000003.346859-7.436564Pressanykeytocontinue搜狗拼音半:六、對算法的理解與感悟?qū)嶒炍逡?、實驗?guī)定分別用下列方法求f(x)=x3?3x?l=0在x0=2附近的根。根的準確值為x*=1.87938524…，規(guī)定準確到四位有效數(shù)字，并對比各種算法的計算量。(1)二分法；(2)簡樸迭代法；(3)牛頓迭代法二、實驗?zāi)康?、通過對二分法和牛頓迭代法作編程練習(xí)和上機運算,進一步體會它們在方程求根中的不同特點；2、比較兩者的計算速度和計算精度。三、算法描述1、二分法環(huán)節(jié)一準備計算f(x)在有根區(qū)間[a,b]端點處的值f(a),f(b).環(huán)節(jié)二二分計算f(a)在區(qū)間中點(a+b)/2處的值f[(a+b)/2]。環(huán)節(jié)三判斷若f[(a+b)/2]=0,則(a+b)/2為根，計算結(jié)果結(jié)束，否則檢查；若f[(a+b)/2]<0,則以(a+b)/2代替b,否則替代a。反復(fù)執(zhí)行環(huán)節(jié)二和環(huán)節(jié)三，直到區(qū)間[a,b]的長度小于允許誤差，此時的中點即為所求。2、簡樸迭代法環(huán)節(jié)一準備提供一個初始值x0；環(huán)節(jié)二迭代計算迭代值xl二a。)；環(huán)節(jié)三控制檢查|xl-x0|;若|xl-xO|>￡(￡為用戶輸入的精度)，則以xl替換x0繼續(xù)迭代；當Ixl-xOl時終止計算，取xl作為所求的結(jié)果。3、牛頓迭代法給定初值班，，為根的允許誤差，〃為|/(幻|的允許誤差,N為迭代次數(shù)的允許值。(1)假如/(x)=0或迭代次數(shù)大于N,則算法失敗,結(jié)束;否則執(zhí)行(2)。(2)計算，下施-/(X)/f(幻(3)若卜「不卜&或|/(九)|<〃，則輸出尤，程序結(jié)束；否則執(zhí)行(4)。

（4）令乂=④，轉(zhuǎn)向（1）。四、源程序1、二分法#inc1ude<iostream>usingnamespacestd;f1oatjisuan(f1oatx);intmain()3C0Ut<<”************3C0Ut<<”************3C0Ut<<”**************3C0Ut<<”**************使用二分法求解**************”<<endl;^cout?"請輸入f(x)的有根區(qū)間［a,b］的下限和上限:floata,b,c;cin?a?b;cout<V”請輸入允許誤差：”；f1oatw;笛in?w;floatfl,fr,fz；ofl=jisuan(a);fr=jisuan(b);for(inti=l;(b-a)>=w；i++)°{fz=jisuan((a+b)/2);00ocout?H二分法求解函數(shù)的近似根為："<<(a+b)/2?endl;break;)else(bcoutv<"第”《i?”個分點為:"；ooprintf(n%f\nu,(a+b)/2)；oooC=(a+b)/2;。if(fz*fl<0)oooob=(a+b)/2；。eIse。a=(a+b)/2;00}}。cout<<"二分法求解函數(shù)的近似根為printf(H%f\n",c);cout<V”迭代次數(shù)為：H?i-l<<endl;return0；floatjisuan(f1oatx)f1oatfx;fx=x*x*x-3*x-1;。returnfx;)2、簡樸迭代法include<iostream>inc1ude<math.h>usingnamespacestd;intmain()(cout<〈”請輸入允許誤差:";f1oatw;cin>>w;afloatx0,xl;笛OUt?”請輸入初值:"；cin>>xO;=x0;intk=0;實驗實驗規(guī)定10000用兩種不同的順序計算￡ri"。1-^834,分析其誤差的變化。"=1二、實驗?zāi)康?、通過上機編程，復(fù)習(xí)鞏固以前所學(xué)程序設(shè)計語言；2、通過上機計算，了解舍入誤差所引起的數(shù)值不穩(wěn)定性。3、通過上機計算，了解運算順序?qū)τ嬎憬Y(jié)果的影響，從而盡量避免大數(shù)吃小數(shù)的現(xiàn)象。三、算法描述一方面用float精度正序相加得到10000項的和；)接著用float精度逆序相加得到10000項的和；3)然后用double精度正序相加得到10000項的和;小)再用double精度逆序相加得到10000項的和;④5)根據(jù)結(jié)果,進行比較分析問題。四、源程序#inc1ude<stdio.h>#include<math.h>voidmain()(°inti;?floatsl=0?s2=0,dl,d2;for(i=l；i<=l0000;i++)sl=sH-1.0f/(i*i);ofor(i=l0000;i>=l;i—)do{0=x1;1=pow((3*xO+1),1/3.0)；obk++;o}whi1e(fabs(x1—x0)>w);，coutV<“簡樸迭代法求解近似根為「<VxlVVendl<<”迭代次數(shù)為:"<vk?endl;return0;)3、牛頓迭代法#inc1ude<iostream>#include<math.h>usingnamespacestd;intmain()(ocoutvv”請輸入初始近似值:"；floatx0;?cin?x0;eCOUt<V”請輸入允許誤差:“；f1oatw；?cin?w；floatx1,f0,f;of0=x0*x0*x0-3*x0-1;f=3*x0*x0-3;oxl=x0-f0/f;intk；ofor(k=l;fabs(x1-x())>w；k++)exO=X1；0fo=x0火xO*x0-3*x0-1；bf=3*xO*xO-3;xl=x0—fO/f;)ecoutV<"牛頓迭代法求解近似根為:"<<xl<<end1"”迭代次數(shù)為:“<<kVvendl;return0;五、運營結(jié)果*C:\Progra>Files\licrosoftVisualStudioMyProjects\123\Debug\123.exe***************使用二分法求解**************情掩入￡<x）的有根區(qū)間L的下限和上限：01常播入研凝:0.0005使1個分點為：。.500000笫2個分點為：0.750000第3個分點為：0.8750001-4個分點為：0.9375005個分點為：0.9687506個分點為：0.9843757個分點為992188用8個分點為：0.996094使9個分點為：。.998047塞。個分點為：0.999023弟11.個分點為：0.999512二分法求解函數(shù)的近似根為：。.999512其代次數(shù)為PressanykeytocontinueI'cT*C;\Progra?Files\licrosoftVisualStudio\lyProjects\123\Debug\123.exe-HEIB請輸入允誑誤差:。.。。即耀全颼求解近似根為工87945選代次數(shù)為：6Pressanykeytocontinue.六、對算法的理解與感悟?qū)嶒灹鶎嶒炓?guī)定分別用高斯列主元消去法和直接三角分解法(LU分解)求方程組的解系數(shù)矩陣：10787常向量：10TOC\o"1-5"\h\z75658861096759107精確解為：(-60,102,?27,16)二、實驗?zāi)康?)熟悉求解線性方程組的有關(guān)理論和方法；2)會編制列主元消去法、LU分解法、雅可比及高斯-塞德爾迭代法的程序;3)通過實際計算，進一步了解各種方法的優(yōu)缺陷，選擇合適的數(shù)值方法。三、算法描述1、列主元高斯消去法算法將方程用增廣矩陣[A|b]=(ajj)ny包表達1)消元過程對k=l,2,n-1①選主元，找ikw{k,k+l,…,n}使得max..“卜k②假如小卜=0,則矩陣A奇異，程序結(jié)束；否則執(zhí)行③。③假如ikwk,則互換第k行與第ik行相應(yīng)元素位置，a.<->a.a.<->a.19ikja.<->a.19a.<->a.19ikjj=k,—,n+l④消元,對i=k+1,—,n計算4=。認/akk對j=l+l,…，n+1計算%=%一-2)回代過程①若明“=0,則矩陣A奇異，程序結(jié)束；否則執(zhí)行②。②怎；對i=n-1,—,2,1,計算(n\七二%〃+1-Xauxj/%Ij=i+l72、矩陣直接三角分解法將方程組Ax=b中的A分解為A=LU,其中L為單位下三角矩陣，U為上三角矩陣,則方程組Ax=b化解為2個方程組Ly=b,Ux=y,具體算法如下(1)對j=1,2,3,…,n計算％=a|j,對i=2,3,…,n計算。=au/aH(2)對k=1,2,3,:a.Mj=k,k+1,??,nilWu^=a&q=ik-l3種=1<+1,1<+2,?一,11計算口二匹一工1%)/%q=1k-l⑶=3,對k=2,3,???,[]計算丫卜=bkq=i(4)xn=yn/unn,^k=n-l,n-2,---,2,li+Mxk=(yk-1kqXq)/也q=k+i注：由于計翔的公式與計購的公式形式上一樣，故可直接對增廣矩陣aila12…ainai,n+l[A|b]=a?'%…a2na,.向_^nl^n2…^nn^n,n+1_施行算法(2)(3),此時U的第n+1列元素即為y。四、源程序#include<stdio.h>main()(doublea|4|[4|={10,7,8,7,7,5,6,5,8,6,10,9,7,5,9,10);doubleb[5]={10,8,6,7},x[4]={0},m,1[4][4],u[4][4],t=0,y[4];inti,j,k,choice;printf("矩陣為:\n");for(i=0;i<4；i++){for(j=0;j<4;j++)printf[i][j])；printf(n\nu);}printf("選擇方法(l:高斯列主元消去法2：LU分解法3:結(jié)束)：\n"”？");scanf(n%dn,&choice);while(choice!=3)switch(choice)case1:(for(i=0;i<4;i++)for(j=i+l;j<4；j++)(m=a|j][i|/a[i][i];for(k=i;k<4；k++)(a[j][k]-=a[i][k]*m;}b[j]—=b[i]*m;)m=();for(i=3；i>=0;i--){for(j=i+l;j<4;j++){m+=aLi])x[i]=(b[i]-m)/a[i][i];m=();}printf("\n方程解為：”)；printf(H(H);for(i=0;i<4;i++)printf(n%.2Ifn,x[i]);printf(")H)；}break；case2:(for(i=0;i<4;i++){u[0][i]=a[0][i];l[i][0]=a[i][0]/u[0][0];)for(k=l;k<4;k++)(for(i=k;i<4；i++)(for(j=O;j<=k-1;j++)t=t+lLk]U]*uU][i];u[k][i]=a[k][i]-t;t=0;)for(i=k+l;i<4;i++)for(j=0;j<=k—1;j++)t=t+l[i][j]*u[j][k];l[i][k]=(a[i][k]-t)/u[k][k];t=0;})y[0]=b[0]；for(i=1；i<4；i++){for(j=0；j<=i-l；j++)t+=l：i][J]^y[J];y[i]=b[i]-t;t=0;}x[3]=y[3]/u[3][3]；for(i=2;i>=0;i--){for(j=i+l;j<=3;j++)t=t+u[i][j]*x[j];x[i]=(y[i]-t)/u[i][i];t=0;)printf(”\n方程解為：")；

?d1=(float)(fabs(l.644834-s1));。d2=(float)(fabs(1.644834-s2));叩rinif("正向求和結(jié)果為％f\n,誤差為％f\n\n",sl,dl);printf("反向求和結(jié)果為％f\n,誤差為％f\n\n”,s2,d2);df(d1<d2)。printf("正向求和誤差小于負向求和誤差\n”)；elseif(dl=d2)叩rimf("正向求和誤差等于負向求和誤差\n")；elseprintf("正向求和誤差大于負向求和誤差\n")；I五、運營結(jié)果*C:\Progra>Files\IicrosoftVisualStudioMyProjects\vucha\Debug\vuch...罡茍曩嘉蕊為i.644725，誤愛為e.000109反向求和結(jié)果為1.644834,誤差為0.000000正向求和誤差等于負向求和誤差Pressanykeytocontinue六、對算法的理解與感悟printf(H(n);for(i=0;i<4;i++)printf("％.21fn,x[i]);printf(H)");}break；)printf(H\n\n選擇方法(1:高斯列主元消去法2:LU分解法3：結(jié)束):\n"”？");scanf(n%dn,&choice);)五、運營結(jié)果c：-C:\Progra?FilesMicrosoftVisualStudio\lyProjects\123\Debug\123.exe^花陣為:10.0000?70500*8.00P.005但擇方法8.00?.0000?50010.009.00.009.0010.00Q：高斯列主元消去法2：LU分解法3：結(jié)束）：方程解為：（-60.00102.00-27.0016.00j他擇方法Q：高斯列主元消去法2：川分解法3：結(jié)束）：,2而■程解為：（-60.00102.00-27.0016.00）卜擇方法（「高斯列主元消去法2：LU分解法3：結(jié)束）：ressanykeytocontinue六、對算法的理解與感悟?qū)嶒灦?、實驗?guī)定1、拉格朗日插值按下列數(shù)據(jù)Xi?3.0-1.01.02.03.0y.-1.01.52.02.01.0作二次插值，并求X|=-2,X2=0,X3=2.75時的函數(shù)近似值2牛頓插值按下列數(shù)據(jù)Xj0.300.420.500.580.660.72y.-1.044031.084621.118031.156031.198171.23223作五次插值，并求X[=0.46,x2=0.55,x3=0.60時的函數(shù)近似值.二、實驗?zāi)康?、熟悉拉格朗日插值多項式和牛頓插值多項式，注意其不同特點;三、算法描述拉格朗日插值(1)輸入x.,y,.(i=0,1,2,n),令L(x”)=0;(2)對二0,1,2,...,n計算、-A-x-x：l(x)=n——-i=o評L.<-L*+1j(x)y,牛頓插值多項式(1)輸入n,x,y(i=0,1,2...,n);對k=l,2,3n,i=1,2,...?k計算各階差商f(x.,x1…，xA);計算函數(shù)值Nn(x)=f(x0)+f[x(),xj(x—x(l)+...+f[xpx1,...,xn](x-xfl)(x-X])...(X-XM_()四、源程序.拉格朗日#include<iostream>#inc1ude<iomanip>#include<math.h>usingnamespaces(d;floatjisuanmokuai(int,float,float[],float[])；voidmain(){41oata,xl50],y[50],1;。intm.n;ocout?"^目:按下列數(shù)據(jù)"v<endl;ocout<<"x：-3.0-LO1.02.03.0"?end1;8cout?"y:l.01.52.02.01.5"?end1;ocoutvV"作二次插值，求x=-2,x=0,x=2.75時的函數(shù)近似值."vvend1;?coutVv"輸入插值次數(shù):"vVend1；?cin>>n；◎coutvv"插入計算次數(shù):"《end1;ocin>>m;?for(inti=0；i<=m；i++)?for(intj=O;j<=n；j++)。{。?inth99=i+l,b99=j+l；。ocoutvv”輸入第"?h99vv”個數(shù)的第”vvb99<v"個節(jié)點及函數(shù)值：H?end1；°。cin>>x[j]>>y[j]；00|2.牛頓#inc1ude<stdio.h>#defineM6floatw(f1oatX,intn,f1oatx[])；floatquoticnt(intk,inti,floatf[][M],fIoatx[],floaty[]);floatnewton(f1oatX,intn,f1oatfoatx(l,floaty[]);voidmain(){afloatx[M],y[M];ox[0]=0.30f；y[0]=l.04403f；x[l]=0.42f；y[l]=1.08462f；x[2]=0.50f;y[2]=l.11803f;ox[3]=0.58f；y[3]=1.15603f；?x[4]=0.66fM4]=1.19817f;x[5]=0.72f；y[51=1.23223f;oprintf("x:0.300.420.500.580.660.72\n");叩rintf("y:1.044031.084621.118031.156031.198171,23223\n")；printf("做五次插值，并求x=0.46,x=0.55,x=0.60時的函數(shù)近似值.\n\n");for(inti=O；i<M；i++)floatN[3]；N[0]=ne