《tomxwd数学数学实验》由会员分享，可在线阅读，更多相关《tomxwd数学数学实验（18页珍藏版）》请在金锄头文库上搜索。
1、一、 实验目的：算；二、实验报告题：蒲丰（Buffon）的随机掷针法是计算的另一种随机模拟，其方法为：在白纸上画许多条等距为d的平行线，将一根长d/2的直针随机掷向白纸，若n次掷针中有m次与平行线相交，当n很大时的近似值为n/m。试证明这一结论，并用MATLAB软件做实验。二、 源代码：M文件(Buffon.m)：function kk=Buffon(n)d=10;x=unifrnd(0,d/2,n,1);c=unifrnd(0,pi,n,1);P=x K=Buffon(1000000)K =3.1434184675835 KK=Buffon(10000000)KK =3.1425146087
2、6479三、结果：投掷一百万次结果为3.1434184675835，投掷一千万次结果为3.14251460876479。小结：随着实验次数的增多，确实会越来越接近。三、 实验目的：学会用matlab实现插值与数值积分，并解决实际问题。二、实验报告题：1.在桥梁的一端每隔一段时间记录1min有几辆车过桥，得到数据，估计一天通过桥梁的车流量。（用线性插值和三次样条插值实现）2.给出机翼剖面的轮廓线数据，求机翼剖面的面积。（用梯形和辛普森实现）三、源代码：第一题：x1=0 2 4 5 6 7 8 9 10.5 11.5 12.5 14 16 17 18 19 20 21 22 23 24;y1=2
3、2 0 2 5 8 25 12 5 10 12 7 9 28 22 10 9 11 8 9 3;x1=x1.*60;x=0:0.1:24*60;%线性插值y=interp1(x1,y1,x);subplot(121),plot(x1,y1,*,x,y)title(线性插值),xlabel(时间),ylabel(车流量),legend(原始数据,线性插值)%三次样条插值y2=interp1(x1,y1,x,spline);subplot(122),plot(x1,y1,*,x,y2)title(三次样条插值),xlabel(时间),ylabel(车流量),legend(原始数据,三次样条插值)%
4、利用梯形公式求积分y11=trapz(x,y)y12=trapz(x,y2)第二题：x1=0 3 5 7 9 11 12 13 14 15;y1=0 1.8 2.2 2.7 3.0 3.1 2.9 2.5 2.0 1.6;y2=0 1.2 1.7 2.0 2.1 2.0 1.8 1.2 1.0 1.6;%三次样条插值x=0:0.1:15;y11=interp1(x1,y1,x,spline);y12=interp1(x1,y2,x,spline);plot(x1,y1,*,x1,y2,*);title(机翼剖面面积);hold on;plot(x,y11,x,y12);legend(原始数据y
5、1,原始数据y2,y1插值,y2插值);hold off;%复合梯形公式xx=0:0.02:15;n=length(x)-1;T1=0;T2=0;for k=1:n T1=T1+(y11(k)+y11(k+1)*(x(k+1)-x(k)/2; T2=T2+(y12(k)+y12(k+1)*(x(k+1)-x(k)/2;endT=T1-T2%复合辛普森公式m=(length(x)-1)/2;S1=0;S11=0;S2=0;S22=0;h=(15-0)/(2*m);for k=1:m S1=S1+(y11(2*k+1); S2=S2+(y12(2*k+1);endfor k=2:m S11=S11
6、+y11(2*k); S22=S22+y12(2*k);ends1=h/3*(y11(1)+y11(2*m+1)+4*S1+2*S11);s2=h/3*(y12(1)+y12(2*m+1)+4*S2+2*S22);S=s1-s2四、结果：第一题：答:用线性插值得到车流量为12990辆/天，用三次样条插值得到车流量约为12669辆/天。第二题：答：梯形公式求得面积约为11. 3444，复合辛普森公式求得面积约为11.3388。小结：插值与数值积分很实用，掌握了之后可以解决很多问题。四、 实验目的：了解并掌握改进欧拉法与龙格-库塔法解决数值微分问题，并用MATLAB实现。解决实际例题。 二、实验报
7、告题：三、源代码：创建一个M文件名为Untitled；format long;x=pi/2:0.01*pi:10*pi;h=0.01*pi;y=2,-2/pi;f=(x,y)(y(2); -(x*y(2)+(x2-0.52)*y(1)/x2);%改进欧拉法x=x(:);n=length(y);m=length(x);Y=zeros(m,n);Y(1,:)=y(:);for k=2:m Y(k,:)=Y(k-1,:)+h*feval(f,x(k-1),Y(k-1,:); YY=Y(k-1,:)+h*feval(f,x(k),Y(k,:); Y(k,:)=0.5*(Y(k,:)+YY);endfi
8、gure,subplot(311),plot(x, Y(:,1),title(改进欧拉法)%龙格-库塔法opt=odeset(reltol,1e-6,abstol,1e-9);x,L=ode45(f,x,y);subplot(312),plot(x,L(:,1),title(龙格-库塔法)%精确解J=sin(x).*sqrt(2*pi./x);subplot(313),plot(x,J),title(精确解)%误差分析disp(改进欧拉法最大误差为：,num2str(max(abs(Y(:,1)-J)disp(龙格-库塔法最大误差为：,num2str(max(abs(L(:,1)-J)四、结果
9、： Untitled改进欧拉法最大误差为：0.0021394龙格-库塔法最大误差为：0.0013725小结：用各种方法解决数值微分问题时，会得到不一样精度的解，相比之下，龙格-库塔法会比改进欧拉法精确一点。五、 实验目的：1.掌握用MATLAB软件求解非线性方程和方程组的基本用法，并对结果作初步分析。2.练习用非线性方程和方程组建立实际问题的模型并进行求解。二、实验报告题：P137 3.(1)小张夫妇以按揭方式贷款买了1套价值20万元的房子，首付了5万元，每月还款1000元，15年还清。问贷款利率是多少？(2)某人欲贷款50万元购房，他咨询了两家银行，第一家银行开出的条件是每月还4500元，1
10、5年还清；第二家银行开出的条件是每年还45000元，20年还清。从利率方面看，哪家银行比较优惠(简单地假设年利率=月利率*12，可参考(1)的方法）？设贷款月利率是r，总贷款为x0，第k次还款后剩余贷款为xk，每月固定还款为b元。三、源代码：1)解：若设贷款月利率是r，总贷款为x0，第k次还贷后剩余贷款为xk，每月固定还b元，可得到方程：x1=x0*(1+r)-bx2=x1*(1+r)-bx3=x2*(1+r)-b xk=xk-1*(1+r)-b则可得 x1=x0*(1+r)-bx2=x0*(1+r)2-b*(1+r)-bx3=x0*(1+r)3-b*(1+r)2+(1+r)+1)xk=x0*
11、(1+r)k-b*(1+r)(k-1)+(1+r)(k-2)+(1+r)+1)整理xk可得：xk=x0*(1+r)k+b*(1-(1+r)k)/r;又因为总贷款为15万，总月份为15*12个月，每月还款1000元，则x0=15，xk=0，k=180，b=0.1；带进方程xk，得到：0=15*(1+r)180+0.1*(1-(1+r)180)/r;用MATLAB实现非线性方程的求解：代码如下： r=fzero(inline(15*(1+r)180-0.1*(1-(1+r)180)/(-r),0.1)r = 0.0020811638894602)解：由第一题得到的方程xk，带入第一家和第二家的数据
12、，总贷款为50万，总月份分别为15*12和20*12个月，每月还款分别为4500和45000/12，则x0=50，xk1=xk2=0，k1=180，k2=240，b1=0.45，b2=0.375；带入xk可得：0=50*(1+r1)180+0.45*(1-(1+r1)180)/r1; 0=50*(1+r2)240+0.375*(1-(1+r2)240)/r2; 用MATLAB实现非线性方程的求解：代码如下： r1=fzero(inline(50*(1+r1)180+0.45*(1-(1+r1)180)/r1),0.1)r1 = 0.005850792582845 r2=fzero(inline
13、(50*(1+r2)240+0.375*(1-(1+r2)240)/r2),0.1)r2 = 0.005479234984688所以，从月利率方面看，第二家银行更优惠。小结：通过本次实验，掌握了用MATLAB软件求解非线性方程和方程组的基本用法，并对结果作初步分析，用非线性方程和方程组建立实际问题的模型并进行求解。获益良多。六、 实验目的：用迭代法求方程的根，解决数据的统计与分析问题。二、实验报告题：P273, 第7题。(用迭代法求方程的根时，取初值为 ，迭代误差为0.5)(用迭代法求上述方程的根，取初值n=2000)；三、源代码： format long;n=2000;miu=2000;sgm=50;A=0.5;K=50000;b=0.5;c=0.35;N=1-(A-2*A/K*n-c)/(b-c);N1=norminv(N,miu,sgm);NNN=;for i=1:1000 NNN=NNN n; if abs(n-N1)=0.5 n=N1; N=1-(A-2*A/K*n-c)/(b-c); N1=norminv(N,miu,sgm); else break endendNNNni 四、结果： baotongNNN = NNN = 1.0e+03 * 2.
