数学建模数模第一次作业(章绍辉版)

1．（1）n=101;x1=linspace(-1,1,n);x2=linspace(-2,2,n);y1=[sqrt(1-x1.^2);-sqrt(1-x1.^2)];y2=[sqrt(4-x2.^2);-sqrt(4-x2.^2);sqrt(1-(x2.^2)/4);-sqrt(1-(x2.^2)/4)];plot(x1,y1)holdon;plot(x2,y2)title('椭圆x^2/4+y^2=1的内切圆和外切圆')axisequal-2.5-2-1.5-1-0.500.511.522.5-2-1.5-1-0.500.511.52椭圆x2/4+y2=1的内切圆和外切圆（2）x1=linspace(-2,2,101);x2=linspace(-2,8);axisequalplot(exp(x1),x1,x1,exp(x1),x2,x2)title('指数函数y=exp(x)和对数函数y=ln(x)关于y=x对称')-2-1012345678-2-1012345678指数函数y=exp(x)和对数函数y=ln(x)关于y=x对称（3）holdonq=input('请输入一个正整数q;')fori=1:qforj=1:iifrem(j,i)plot(j/i,1/i)endendend00.10.20.30.40.50.60.70.80.9100.050.10.150.20.250.30.350.40.450.53．代码如下：n=input('请输入实验次数n=')k=0;fori=1:nx=ceil(rand*6)+ceil(rand*6);ifx==3|x==11k=k+1;elseifx~=2&x~=7&x~=12y=ceil(rand*6)+ceil(rand*6);whiley~=x&y~=7y=ceil(rand*6)+ceil(rand*6);endify==7k=k+1;endendend更改试验次数n的值，打赌者赢得概率w随n的变化情况如下：试验次数n打赌者赢得概率w10000.528015000.495320000.499225000.515730000.501135000.511540000.511745000.511050000.5041从上表可看出打赌者赢的概率大约为0.5110。理论计算：掷一次骰子，得到点数及相应的概率点数23456789101112概率p1/362/363/364/365/366/365/364/363/362/361/36打赌者赢的情况有两种：（1）第一次就掷出3点或者11点；其概率P1=2/36+2/36=1/9；（2）当第1次掷出的点数之和是4，5，6，8，9或10，，继续掷骰子，直到掷出的点数之和是7或原来的值为止，先得到的点数之和是7；其概率P2=196/495则打赌者赢的概率P=P1+P2=0.50707.4.（1）（i）输入代码：f=@(r,t)3.9*exp(r.*(t-1790));t=1790:10:2000;c=[3.9,5.3,7.2,9.6,12.9,17.1,23.2,31.4,38.6,50.2,62.9,76.0,92.0,106.5,123.2,131.7,150.7,179.3,204.0,226.5,251.4,281.4];r0=0.036;r=nlinfit(t,c,f,r0)sse=sum((c-f(r,t)).^2)得到：r=0.02119sse=17418.48（ii）输入代码：f=@(k,t)k(1)*exp(k(2).*(t-1790));k0=[3.9,0.036];k=nlinfit(t,c,f,k0)sse=sum((c-f(k,t)).^2)得到：k=14.993959233344270.01422307528878sse=2.263917490357360e+003即：x0=14.994r=0.014223sse=2263.92(iii)输入代码：f=@(k,t)k(1)*exp(k(2).*(t-k(3)));k0=[3.9,0.036,1790];k=nlinfit(t,c,f,k0)sse=sum((c-f(k,t)).^2)得到：k=1.0e+003*0.007529449634170.000014223084281.74157000966658sse=2.263917490325156e+003即：x0=7.52945r=0.014223t0=1741.57sse=2263.92从误差平方和sse来看（ii）和（iii）的拟合效果较好。（2）对两边取对数得令y=,x=t-1790,A=r,B=，则原方程变为：y=Ax+B。用polyfit拟合参数A、B，代码如下：t=1790:10:2000;c=[3.9,5.3,7.2,9.6,12.9,17.1,23.2,31.4,38.6,50.2,62.9,76.0,92.0,106.5,123.2,131.7,150.7,179.3,204.0,226.5,251.4,281.4];x=t-1790;y=log(c);k=polyfit(x,y,1)r=k(1),x0=exp(k(2))sse=sum((c-exp(polyval(k,x))).^2)plot(t,c,'k+',t,exp(polyval(k,t-1790)),'k-')xlabel('t'),ylabel('x(t)'),legend('22个已知数据点’,2)title('polyfit函数拟合效果图')得到：k=0.0202193337618981.799226697714663r=0.020219333761898x0=6.044971066570838sse=3.489176713658117e+004175018001850190019502000050100150200250300350400450tx(t)polyfit函数拟合效果图22个已知数据点即：A=0.020219B=1.7992r=A=0.020219x0=6.04497误差平方和为34891.77。（3）比较（1），（2）两种拟合方式的误差平方和，显然线性拟合带来的误差比非线性拟合大得多。用函数plot绘制两种拟合方式的误差比较图，输入代码：f=@(k,t)k(1)*exp(k(2).*(t-1790));k0=[3.9,0.036];k1=nlinfit(t,c,f,k0);x=t-1790;y=log(c);k2=polyfit(x,y,1);plot(t,c-f(k1,t),'k+',t,c-exp(polyval(k2,x)),'kp')title('非线性拟合与线性拟合误差比较效果图')xlabel('t'),ylabel('x(t)'),legend('非线性拟合',2,'线性拟合',2)得到：175018001850190019502000-160-140-120-100-80-60-40-2002040非线性拟合与线性拟合误差比较效果图时间t误差非线性拟合线性拟合由图可知非线性拟合产生的误差比线性拟合均匀得多，可能是由于作了变换Y=ln(x(t))。因此，这里采用非线性拟合比较合理。（4）（i）输入代码：f=@(k,t)3.9*k(2)./(3.9+(k(2)-3.9)*exp(-k(1)*(t-1790)));t=1790:10:2000;c=[3.9,5.3,7.2,9.6,12.9,17.1,23.2,31.4,38.6,50.2,62.9,76.0,92.0,106.5,123.2,131.7,150.7,179.3,204.0,226.5,251.4,281.4];k0=[0.1,5.3000];k=nlinfit(t,c,f,k0)sse=sum((c-f(k,t)).^2)得到：k=0.0274342.4423sse=1.2249e+003即：r=0.0274N=342.44误差平方和为：1224.9(ii)输入代码：f=@(k,t)k(3)*k(2)./(k(3)+(k(2)-k(3))*exp(-k(1)*(t-1790)));k0=[0.1,5.3000,5.3000];k=nlinfit(t,c,f,k0)sse=sum((c-f(k,t)).^2)得到：Warning:Rankdeficient,rank=2,tol=2.4945e-014.Innlinfitat161k=0.0215446.57327.6981sse=457.7405即：r=0.0215N=446.5732x0=7.6981误差平方和为：457.7405（iii）输入代码：f=@(k,t)k(3)*k(2)./(k(3)+(k(2)-k(3))*exp(-k(1)*(t-k(4))));k0=[0.1,5.3000,5.3000,1790];k=nlinfit(t,c,f,k0)sse=sum((c-f(k,t)).^2)得到：Warning:Rankdeficient,rank=2,tol=2.5942e-014.Innlinfitat161k=1.0e+003*0.000021547329460.446572485004110.002670029180211.74033082470710sse=4.577406556557742e+002即：r=0.02155N=446.57x0=2.670t0=1740.33误差平方和为：457.7406此处（i）的误差平方和很大，（ii）和（iii）的误差平方和几乎一样且matlab都出现了警告（这点跟书上的参考答案不一样，书上是只有（iii）出现了警告，但拟合出来的参数值又跟书上一样。不知道为什么）所以这里我把（ii）和(iii)的拟合效果图都给出：175018001850190019502000050100150200250300(ii)拟合效果图175018001850190019502000050100150200250300(iii)拟合效果图