东北大学秦皇岛分校数值分析实验报告

数值分析实验报告学 号：[2133326]姓 名：［孙丽］提交日期：[2015-6-25]成 绩：东北大学秦皇岛分校【实验名称】 多项式插值【实验目的】掌握拉格朗日插值和牛顿插值；掌握用mathematica编写多项式插值的应用程序；【实验内容】给定函数f(x)=、云-x2/2,xe［—3,3］，将区间［-3,3］剖分为12等份，取节点x.--3+(i—1)x0.5,7=1,2,...,13，分别构造拉格朗日插值和牛顿插值多项式(x)和p(x)，比较两个多项式的计算时间和精度，最后绘出p(x)和f(x)-p(x)的图形。【实验结果】f[x_]:=Sqrt[2Pi]Exp[-xA2/2]g=Plot[f[x],{x,-3,3},PlotStyle->{Red,Thick},BaseStyle->FontSize->15,AxesLabel->{x,y},AxesStyle->Thick]软件工程实验报告I软件工程i,n}];s]lagrange[x_,xp_,yp_]:=Module[{n,s=0},n=Length[xp];Do[s=s+Product[(x-xp[[j]])/(xp[[i]]-xp[[j]]),{j,i-1}]Product[(x-xp[[j]])/(xp[[i]]-xp[[j]]),{j,i+1,n}]yp[[i]],i,n}];s]newton[x_,xp_,yp_]:=MOdule[{n,p,yp0=yp},n=Length[xp];Do[Do[ypO[[jT]=(yp0[[j]]-yp0[[j-1]])/(xp[[j]]-xp[[j-i]]),{j,n,i+1,-1}],{i,n-1}];p=yp0[[n]];Do[p=p(x-xp[[i]])+yp0[[i]],{i,n-1,1,-1}];p]xd=Range[-3・,3・,0・5];yd=Map[f,xd];dat=Transpose[{xd,yd}];n=Length[dat]l[x_]=lagrange[x,xd,yd];Clear[p]p[x_]=newton[x,xd,yd];Plot[l[x],{x,xd[[1]],xd[[n]]}]//Timing东北大学秦皇岛分校计算机与通信工程学院Plot[p[x],{x,xd[[1]],xd[[n]]}]//Timing4・h[x_]=InterpolatingPolynomial[dat,x];PlOt[h[x],{x,xd[[1]],xd[[n]]}]//Timing软件工程实验报告g=Plot[l[x]-p[x],{x,xd[[1]],xd[[n]]},PlotRang»All,PlotStyle^Blue,BaseStyle^{Bold,FontSize^14},AxesStyle^Thick,AxesLabel^{x,error}]Plot[h[x]-p[x],{x,xd[[1]],xd[[n]]},PlotRange>All]

Clear[p1]p1[x_]=p[x]Clear[p]g=Plot[p1[x],{x,xd[[1]],xd[[n]]},PlotStyle>{Red,Thick},BaseStyle^FontSize^15,AxesStyle^Thick,AxesLabelT{x,p[x]},EpilogT{PointSize[Large],Blue,Point/2@dat}]谿软件工程实验报告II13问冈谿软件工程实验报告II13问冈Epiloj(PointSize[Large].rBlue,Point/@(iat}]Out[99]=Untitled-1*In[98]!=Clear[pl]pl[x]=P[x]0.0278461+(3.+x)(0.164575+(2.5+x)(0.293629+(2.+x)(0.131509+(1.5+x)(-0.074708+(1.+x)(-0.0323582+(0.5+x)(0.0194839+(0.000229414+(-0.00235005+(0.000673994+(-0.0000220602+(-0.0000369849+0.0000123283(-2.5+x))(-2.+x))(-1.5+x))(-1.+幻)(-0-S+x))(0.In[98]:=Clear[p]In[98]:=g=Plot[pl[x],{x,^[[l]],xd[[n]]},PlotStylet{Red,Thick),BaseStyle->Fontsize->15,AxesStyle->Thick,AxesLabel->{略p[x]},errorPlot[f[x]-plEx],{X,Md[[l]],[[II]]},Plotnaw->|100%*9・error=・g=Plot[f[x]-p1[x],{x,xd[[1]],xd[[n]]},PlotRange>All,PlotStyle^{Red,Thick},BaseStyle^{Bold,FontSize^16},AxesLabelT{x,error},AxesStyle^Thick]

Untitled-1*ln[101];=error=.=Plot[f[x]-pl[x],(x,xd[[l]],xd[[n]]},PlotRmigeM.1,Plotstylet{Red,Thick},BaseStylet{Bold,FontSize-*16}z13问冈AKesLabel->{x；error}；AxesStyle-jThick]13问冈软件工程实验报告I■!■■■■【实验编号】 2【实验名称】 数值积分【实验目的】掌握复合梯形公式和复合辛普森公式；掌握用mathematica编写数值积分的应用程序；【实验内容】给定函数f(x)二12+4(sin2x)ex/io,xe[0,11]，动画演示复合梯形公式和符合辛普森公式。【实验结果】1・ntos[n_Integer]:=Module[{s},s=IntegerDigits[n];If[n<10,"0"<>FromCharacterCode[s+48],FromCharacterCode[s+48]]]ntos[45]2・f[x_]:=12+4Sin[2x]Exp[x/10]t1=Plot[f[x],{x,0,11},AxesLabel->{x,y},PlotRange->{{0,12},{0,25}},PlotStyle->{Thickness[0・01],RGBColor[1,0,1]},AxesStyle->Thick]3.a=0;b=11;

毎软件工程 实验报告IIIIIIManipulate[dat=Table[{xp=a+((b-a)i)/n,f[xp]},{i,O,n}];l1=Interpolation[dat,InterpolationOrder->1];g1=Plot[l1[x],{x,a,b},PlotRange->{{0,12},{0,25}},Filling->Axis,FillingStyle->Green];Show[t1,g1,PlotLabelntos[n+1]<>"PointsFigure"],{n,2,40,2}]ntos[n_Integer]:=Module{s},s=IntegerDigits[n];If[n<10,"0"<>FromCharacterCode[s+48],FromCharacterCode[s+48]]]g[x_]:=12+4Sin[2x]Exp[x/10]t2=Plot[g[x],{x,0,11},AxesLabelD{x,g[x]},PlotRangeD{{0,12},{0,25}},PlotStyleD{Thickness[0・01],RGBColor[1,0,1]},AxesStyleDThick,软件工程 实验报告IBaseStyleDFontSizeD16]a=0;b=11;Manipulate[dat=Table[{xp=a+((b-a)i)/n,f[xp]},{i,O,n}];l2=Interpolation[dat,InterpolationOrdezEI2];g2=Plot[l2[x],{x,a,b},PlotRangeC{{0,12},{0,25}},FillingDAxis,FillingStyleDGreen];Show[t2,g2,PlotLabelEntos[n+1]V>"PointsFigure"],{n,2,30,2}]【实验编号】 3【实验名称】 线性方程组求解软件工程 实验报告IIIIIII【实验目的】掌握高斯类消去法求解线性方程组；掌握用mathematica编写直接求解线性方程组的应用程序;【实验内容】列主元高斯消去法求解线性代数方程组Ax=b，其中/21-3-1、(1「31072A=-124-2，b=-1<10-15要求：对矩阵A,找出第2列中从第2个元素开始的绝对值最大元素所在位置，并与第2行进行交换。使用列主元高斯消去法将矩阵A化为上三角形，矩阵A下三角部分存放消元过程中的变换元素。回代求出线性方程组的解，存放在向量b中。将求解Ax=b的列主元高斯消去法写成函数形式。【实验结果】A={{2,1,-3,-1},{3,1,0,7},{-1,2,4,-2},{1,0,-1,5}}i=2;m=First[Dimensions[A]];temp=Flatten[Abs[Take[A,{i,m},{i}]]];imax=Position[temp,Max[temp]][[1,1]]+i-1;{A[[i]],A[[imax]]}={A[[imax]],A[[i]]}A={{2,1,-3,-1},{3,1,0,7},{-1,2,4,-2},{1,0,-1,5}};b={1,2,-1,5};{m,n}=Dimensions[A];Do[temp=Flatten[Abs[Take[A,{k,m},{k}]]];kmax=Position[temp,Max[temp]][[1,1]]+k-1;If[k!=kmax,{A[[k]],A[[kmax]]}={A[[kmax]],A[[k]]};{b[[k]],b[[kmax]]}={b[[kmax]],b[[k]]}];If[A[[k,k]]!=0,Do[A[[i,k]]=A[[i,k]]/A[[k,k]];Do[A[[i,j]]=A[[i,j]]-A[[i,k]]*A[[k,j]],{j,k+1,n}];b[[i]]=b[[i]]-A[[i,k]]*b[[k]],{i,k+1,m}],Print["Error"];Break[]],{k,m-1}];b[[m]]=b[[m]]/A[[m,m]];

谿软件工程 实验报告IIIIIIIDo[b[[i]]=(b[[i]]-Sum[A[[i,j]]*b[[j]],{j,i+1,m}])/A[[i,i]],{i,m-1,1,-1}]；bgausspiv[A_,b_]:=MOdule[{m,n,temp,kmax,A0=A,b0=b},{m,n}=Dimensions[A0];Do[temp=Flatten[Abs[Take[A0,{k,m},{k}]]];kmax=Position[temp,Max[temp]][[1,1]]+k-1;If[旧kmax,{AO[[k]],AO[[kmax]]}={AO[[kmax]],AO[[k]]};{bO[[k]],bO[[kmax]]}={bO[[kmax]],bO[[k]]}];If[A0[[k,k]口0,Do[A0[[i,k]]=A0[[i,k]]/A0[[k,k]];Do[A0[[i,j]]=A0[[i,j]]-A0[[i,k]]*A0[[k,j]],{j,k+1,n}];bO[[i]]=bO[[i]]-AO[[i,k]]*bO[[k]],{i,k+1,m}],Print["Error"];Break[]],{k,m-1}];bO[[m]]=bO[[m]]/AO[[m,m]];Do[bO[[i]]=(bO[[i]]-Sum[AO[[i,j]]*bO[[j]],{j,i+1,m}])/AO[[i,i]],{i,m-1,1,-1}];b0]x=gausspiv[A,b]软件工程实验报告I软件工程【实验编号】 4【实验名称】 非线性方程求根【实验目的】掌握二分法和牛顿法求解非线性方程的根；掌握用mathematica编写非线性方程求根的应用程序;【实验内容】给定非线性方程sinlOx-xe-x二x,完成以下工作：⑴绘出函数f(x)二sinlOx—xe-x-x的图形，确定有根区间。使用二分法求出方程的根，要求精度＜10-8使用牛顿法求出方程的根，要求精度＜10-8【实验结果】(1)bisection[f,a,b,eps,kmax]:=Module[{k=0,v=1,a0=a,b0=b,va,vb,c},While[Sia:"];b=Input["Inputb:"]];Print["a=",a,"b=",b];va=f[a0];vb=f[b0];While[Abs[v]>eps&&k<kmax,c=(a0+b0)/2・0;v=f[c];If[Sign[va]EISign[v],a0=c;va=v,b0=c;vb=v];k++;Print["k=",k,"c=",c,"f(",c,")=",v]];{c,v}]f[x_]:=Sin[10x]-xExp[-x]-x;

res=bisection[f,0・1,0・5,10^-8,100]

SUntitled-1*In[18]：=bisectionb_feps_fkrtax_]:=Module[{k=0irv=l/a0=airl)0=2*irva^irt)』c}』THiile[Sign[f[ii]]==Sign[f[i]],a=Iiqiut["Iiqiuta:b=Iiqiut["In^utJ):11]];Print["a="』nJl)="』*];va=f[aO];irt)=f[l)OJ;THiile[v]kepsr&&k-=： c=(aO+l)0)/2・0;V=f[c];If[Sign[va]==Si(jn[v],raO=c;va=v,rl)0=c;=v];k+;Print ,k,11c=",c,11f(",c,11)=",v]];{C,V}]f[ji^J:=Sin[10x]-JirEiq)[-x]-x;res=bisection[f,r0.1』10f,-8f100]3=1」■丄}j=Li・5k=lc=U・E；f(U・工〕=一[1・3E；丄丄25k=2c=U・2f(U.2)=0-.5455.51k=3c=0.25f(0.25)=0.153772k=4c=U.275f(0・275；1=—匚1■丄02221k=5c=U.2625fi；U.262-5)=0.0295246k=6c=U・2£275f(U・2££75；1=—匚1.0355丄77k=7c=U・2£5£25fi；U・2石5£25；1=—0・匚1匚127了522k=8c=U fi；U・2石4匚|£3；1=匚1・0丄343丄？k=9c=U fi：：i：i・2石48E4；i=0 2也k=10c=0.265234f(0.265234)=0.00128696 2133326^fiHl

软件工程实验报告■Untitled-1*k=2c=U・2fi；n.2；i=Li-545551软件工程实验报告■Untitled-1*k=2c=U・2fi；n.2；i=Li-545551k=3c=U・25f(0・25)=0.153772k=4c=U.275f(0.275)=-0.102221k=.5c=0.2625f(0.2625)=0.0295246k=6c=U.26875f(0.26875)=-0.0355177k=7c=U・265625f(0.265625)=-0.00277528k=8c=Uf(0.264063)=0.0134317k=9c=U・2£EE；电4f(0.264844)=0.00534224k=10c=0・2£5234f(0.265234)=0.00128696 21333265(1^k=llc=0・26543f(0.26543)=-0.000743294k=12c=0・26.5332f(0.265332)=0.00027205k=13c=0.26S381f(0.265381)=-0.000235567k=14c=U・2653-56f(0.265356)=0.0000182552k=15c=0.265369f(0.265369)=-0.000108653k=16c=0・265363f(0.265363)=-0.0000451979k=17c=0.2E5359f(0.265359)=-0.0000134712k=18c=0・265358f(0.265358)=2.39204xlO-6k=19c=0・2£5359f(0.265359)=-5.53955x10-6k=20c=0・2653-58f(0.265358)=-1.57375x10-6k=21c=U・2653-58f(0.265358)=4.09146xl0_Tk=22c=U・2653-58f(0.265358)=-5.82303x10-Tk=23c=0・265358f(0.265358)=-8.65783x10_Sk=24c=0・265358f(0.265358)=1.61284xl0_Tk=25c=0・2653-58f(0.265358)=3.73528xl0_Sk=26c=0・2653-58f(0.265358)=-2.46127x10-0k=27c=0・2653-58f(0.265358)=6.37003xlO-3=A-]J100%]J(2)newton[f_,x0_,eps_,kmax_]:=Module[{k=0,v=1,x00=x0,x1,fd,t},fd[t_]=D[f[t]Print["k=",k,"x0=",x00,"f(",x00,")=",v];While[Abs[v]eps&&k<kmax,x1=xOO-v/fd[xOO];v=f[x1];x00=x1;k++;Print["k=",k,"x1=",x1,"f(",x1,")=",v]];If[Abs[v]>eps,Print["Newton'smethoddidnotproduceanaccurateappro{x1,v}]f[x_]:=Sin[10x]-xExp[-x]-x;res=newton[f,0.4,10A-8,100]OUntitled-1*区k=24c=0・265358k=25c=0・265358f(0.265358)=1.61284xlO_Tf(0.265358)=3.73528xlO_Sk=26c=0・2653-58k=27c=0・2653-58f(0.265358)=-2.46127x10"f(0.265358)=6.37003xlO-3]J]JOut[20]=(0.265358,6.37003x10-9}newtan[f_,xO_feps_f :=Module[{k=0,v=1^xOO=xOfxl,fd,t}「fd[t_]=D[f[t],t];v=f[x00]//U;Print["k=",k,"x0=",xOO,"f("rxOO,")=",v];Wile[Hbs[v]aepsr&&k< xl=xOO-v/fd[x00];V=f[xl];xOO=xl;k++;Print["k=",k,11xl=",xl,"f("rxl,")=",v]];Tf[11js[v]aepsfPrint["Newton'smethoddidnotproduceanaccuratea^iproxunation"]]:{xl,v}]f[x_]:=Sin[10x]-xEsp[-x]-x;res=neirton[f,0.4，10-8,100] 2133326孙ffilk=0XU=LI.4f(0.4)=-1.42493k=lf(0.220507)=0.408126k=2xl=0电555f(U 电555〕=-0・09742丄5xl=0・2£55丄5k=4k=5xl=U.265358f(U・2£55丄5；i=-0・00丄丘3478f(U・2£5358；i=-.5・£丄33丄:：■：:1U-7fi工i・265358：：i=-6・£39丄dxIO"14T33331Out[23]=[U.265358,-6.63913:=<lU"14]J|v|100%▲(3)res=newton[f,-1,10A-8,100]Untitled-1*口问区］f[x_]:=Sin[10x]-jcEsp[-x]-x;res=newton[f,0.4^10A-8,100]k=0x0=0.4f(0.4)=-l.42493k=lxl=0.220507f(0.220507)=0.408126k=2xl=0.274555f(0.274555)=-0.0974215k=3xl=0.265515k=4xl=0.265358k=5xl=0.265358f(0.265515)=-0.00163478f(0.265358)=-5.61331x10_Tf(0.265358)=-6.63913x10-14Out[23]={0.265358,-6.63913x10-14res=newton[f,-1,10A-8,100] 2133326Iftfflk=0XU=-1f(-1)=^1.2623k=lxl=-U・7丄25E；£f(-0.712536)=丄.41941k=2xl=一丄.36777fJ一丄.36777)=5.84217k=3xl=-U・37丄322f(一0・371322)=丄.45061k=4xl=-U・2E4D55f(-0・2也4055〕=-U・k=.5xl=-U.252799f(-0・252799〕二0・00248912k=6xl=-0.252568f(-0.252568)=1.61L34y10k=7xl=-U・2525£8f(-0.252568)=6.78235x10Out[24]={-0.252568,6.78235x10-12}A-|100%▲(4)res=newton[f,-0.1,10A-8,100]ED回冈E=i：iED回冈E=i：ixu=-l£(-lj=4・£五£3K=1xl=-0■了丄忌号36El；-□■了丄££3石〕二丄.41^41k=2xl=一丄.36777fJ-丄.36777)=5.84217k=3xl=-U・37丄322f(-LL371322)=1.45061k=4xl=-0・244055f(-0・244匚I55：i=-U・DE；944E；丄k=5xl=-U.252799fi；-U・252799；i=0・00248912k=6xl=-U.252568f(-U・2525£8；i二丄.61134:：■：:10k=7xl=-0・2525丘8f(-U・2525£8：i二丘・78235:：：:10res=ne^rton[f,r-丄』L0-8,r100]|100%▲HUntitled-1*k=5xl=O.Z6535S£(0.Z6535S)=-□.63913x10-14Out[23]={0.265358,-6.63913x10-11}Out[24]=[-0.252568,6.78235x10-13]res=newton[f,-0.1,10A-8,100] 2133326Shffl]k=0x0=-0.1f(-0.1)=-0.630954k=lxl=0.0979566f(0.0979566)=0.643483k=2xl=-0.0733665f(-0.0733665)=-0.517278k=3xl=0.0247483f(0.0247483)=0.196073k=4xl=-0.000571313ff-0.0005713131=-0.00457015k=s*i=4-F；Knnqyin"f(4.Rfinnq-x-in-s)=3.RARn7vi