2023年数值分析实验报告新编

河北工业大学《数值分析》课程实验实验报告专业:计算机科学与技术班级:计1305班姓名:岳帅鹏学号：133007完毕日期:2023.1.11实验三一、实验规定分别用复化梯形公式和复化辛卜生公式计算f(x)=sin(x)/x的积分，并与准确值比较判断精度。二、实验目的通过实际计算来真实体会复合梯形公式和复合辛普森公式两种数值求积方法的精确度。三、算法描述1、复合梯形公式环节一：拟定函数的区间为【a,b】，并拟定其节点的个数，算出每个社区间的步长h;环节二：根据拟定的h,代入复合梯形公式，算出其相应的结果。从2、复合辛普森公式a环节一:拟定函数的区间为【a,b】，并拟定其节点的个数，算出每个步长h;环节二:根据拟定的h,代入复化辛普生公式,算出相应的积分结果。四、源程序#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("将复化梯形区间划分:");。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\n”,Tn)；oprimf("将复化辛卜生区间划分为:")；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\n",Sn);0dl=(lloat)(fabs(0.9460831-Tn));8d2=(float)(fabs(0.9460831-Sn));printf("复化梯形误差:%f\n".d1);-printf("复化辛卜生误差：%f\n",d2);，if(d1>d2)printf("复化梯形求法精度低于复化辛卜生求法\n”)；。elseif(dl==d2)g。printf("亚化梯形求法精度等于兔化辛卜生求法\n");e1se(d1<d2)g;printf("复化梯形求法精度高于复化辛卜生求法\n")；))五、运营结果|MC:\Progra>Files\IicrosoftVisualStudio\ByProjects\3\Debug\3.exe*置复化梯彭区间划分复化梯形输出：0.945832西复化〒4生区间划分为复化老卜生物出：0.741290复化梯形误群0.000251复位辛卜生谟差：0.204793复化稔形求法精度高于复化辛卜生求法Pressanykeytocontinue.六、对算法的理解与感悟€实验四一、实验规定用改善欧拉方法解初值问题y，=x+y；y(0)=1。0<x<L取步长h=0.1计算，并与准确值y=・x-l-2e'相比较。二、实验目的1)熟悉求解常微分方程初值问题的有关方法和理论，重要是改善欧拉方法。2)会编制上述方法的计算程序。3)针对实习题编制程序，并上机计算其所需要的结果。4)通过对各种求解方法的计算实习，体会各种解法的功能、优缺陷及合用场合,会选取适当的求解方法三、算法描述解一阶常微分方程初值问题y'=f(x,y)a<x<b.y(xo)=Yo将区间[a,b]作n等分,取步长h="卫n欧拉公式为Yi+i=yi+hf(Xi，yJ梯形公式为yi+i=yi+^[f(x,,yl)+f(xi+1,yi+1)]改进欧拉法采用公式yi+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);}五、运营结果药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搜狗拼音半:六、对算法的理解与感悟实验五一、实验规定分别用下列方法求f(x)=x3・3x・l=0在x0=2附近的根。根的准确值为x*=1.87938524…，规定准确到四位有效数字，并对比各种算法的计算量。(1)二分法；(2)简朴迭代法；(3)牛顿迭代法二、实验目的1、通过对二分法和牛顿迭代法作编程练习和上机运算,进一步体会它们在方程求根中的不同特点；2、比较两者的计算速度和计算精度。三、算法描述1、二分法环节一准备计算f(x)在有根区间[a,b]端点处的值f(a),f(b).环节二二分计算f(a)在区间中点(a+b)/2处的值f[(a+b)/2]。环节三判断若f[(a+b)/2]=0,则(a+b)/2为根，计算结果结束，否则检查；若f[(a+b)/2]<0,则以(a+b)/2代替b,否则替代a。反复执行环节二和环节三，直到区间[a,b]的长度小于允许误差，此时的中点即为所求。2、简朴迭代法环节一准备提供一个初始值x0；环节二迭代计算迭代值xl二a。)；环节三控制检查|xl-x0|;若|xl-xO|>£(£为用户输入的精度)，则以xl替换x0继续迭代；当Ixl-xOl时终止计算，取xl作为所求的结果。3、牛顿迭代法给定初值班，，为根的允许误差，〃为|/(幻|的允许误差,N为迭代次数的允许值。(1)假如/(x)=0或迭代次数大于N,则算法失败,结束;否则执行(2)。(2)计算，下施-/(X)/f(幻(3)若卜「不卜&或|/(九)|<〃，则输出尤，程序结束；否则执行(4)。

（4）令乂=④，转向（1）。四、源程序1、二分法#inc1ude<iostream>usingnamespacestd;f1oatjisuan(f1oatx);intmain()3C0Ut<<”************3C0Ut<<”************3C0Ut<<”**************3C0Ut<<”**************使用二分法求解**************”<<endl;^cout«"请输入f(x)的有根区间［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二分法求解函数的近似根为："<<(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<<"二分法求解函数的近似根为printf(H%f\n",c);cout<V”迭代次数为：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;实验实验规定10000用两种不同的顺序计算£ri"。1-^834,分析其误差的变化。"=1二、实验目的1、通过上机编程，复习巩固以前所学程序设计语言；2、通过上机计算，了解舍入误差所引起的数值不稳定性。3、通过上机计算，了解运算顺序对计算结果的影响，从而尽量避免大数吃小数的现象。三、算法描述一方面用float精度正序相加得到10000项的和；)接着用float精度逆序相加得到10000项的和；3)然后用double精度正序相加得到10000项的和;小)再用double精度逆序相加得到10000项的和;④5)根据结果,进行比较分析问题。四、源程序#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<<”迭代次数为:"<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"”迭代次数为:“<<kVvendl;return0;五、运营结果*C:\Progra>Files\licrosoftVisualStudioMyProjects\123\Debug\123.exe***************使用二分法求解**************情掩入£<x）的有根区间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二分法求解函数的近似根为：。.999512其代次数为PressanykeytocontinueI'cT*C;\Progra»Files\licrosoftVisualStudio\lyProjects\123\Debug\123.exe-HEIB请输入允诳误差:。.。。即耀全飕求解近似根为工87945选代次数为：6Pressanykeytocontinue.六、对算法的理解与感悟实验六实验规定分别用高斯列主元消去法和直接三角分解法(LU分解)求方程组的解系数矩阵：10787常向量：10TOC\o"1-5"\h\z75658861096759107精确解为：(-60,102,・27,16)二、实验目的1)熟悉求解线性方程组的有关理论和方法；2)会编制列主元消去法、LU分解法、雅可比及高斯-塞德尔迭代法的程序;3)通过实际计算，进一步了解各种方法的优缺陷，选择合适的数值方法。三、算法描述1、列主元高斯消去法算法将方程用增广矩阵[A|b]=(ajj)ny包表达1)消元过程对k=l,2,n-1①选主元，找ikw{k,k+l,…,n}使得max..“卜k②假如小卜=0,则矩阵A奇异，程序结束；否则执行③。③假如ikwk,则互换第k行与第ik行相应元素位置，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奇异，程序结束；否则执行②。②怎；对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:结束)：\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("正向求和结果为％f\n,误差为％f\n\n",sl,dl);printf("反向求和结果为％f\n,误差为％f\n\n”,s2,d2);df(d1<d2)。printf("正向求和误差小于负向求和误差\n”)；elseif(dl=d2)叩rimf("正向求和误差等于负向求和误差\n")；elseprintf("正向求和误差大于负向求和误差\n")；I五、运营结果*C:\Progra>Files\IicrosoftVisualStudioMyProjects\vucha\Debug\vuch...罡苟曩嘉蕊为i.644725，误爱为e.000109反向求和结果为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：结束):\n"”？");scanf(n%dn,&choice);)五、运营结果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：结束）：方程解为：（-60.00102.00-27.0016.00j他择方法Q：高斯列主元消去法2：川分解法3：结束）：,2而■程解为：（-60.00102.00-27.0016.00）卜择方法（「高斯列主元消去法2：LU分解法3：结束）：ressanykeytocontinue六、对算法的理解与感悟实验二、实验规定1、拉格朗日插值按下列数据Xi・3.0-1.01.02.03.0y.-1.01.52.02.01.0作二次插值，并求X|=-2,X2=0,X3=2.75时的函数近似值2牛顿插值按下列数据Xj0.300.420.500.580.660.72y.-1.044031.084621.118031.156031.198171.23223作五次插值，并求X[=0.46,x2=0.55,x3=0.60时的函数近似值.二、实验目的1、熟悉拉格朗日插值多项式和牛顿插值多项式，注意其不同特点;三、算法描述拉格朗日插值(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);计算函数值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«"^目:按下列数据"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时的函数近似值."vvend1;«coutVv"输入插值次数:"vVend1；»cin>>n；◎coutvv"插入计算次数:"《end1;ocin>>m;«for(inti=0；i<=m；i++)«for(intj=O;j<=n；j++)。{。«inth99=i+l,b99=j+l；。ocoutvv”输入第"«h99vv”个数的第”vvb99<v"个节点及函数值：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时的函数近似值.\n\n");for(inti=O；i<M；i++)floatN[3]；N[0]=ne