清華數(shù)學(xué)實(shí)驗(yàn)復(fù)習(xí)試題八(蒙特卡羅方法-龍格-庫(kù)塔方法)

1、考試課程 數(shù)學(xué)實(shí)驗(yàn) 2005.6.15下午班級(jí) 姓名 學(xué)號(hào) 得分 說(shuō)明（1）第一、二、三題的答案直接填在試題紙上；（2）第四題將數(shù)學(xué)模型、簡(jiǎn)要解題過程和結(jié)果寫在試題紙上；卷面空間不夠時(shí)，請(qǐng)寫在背面；（3）除非特別說(shuō)明，所有計(jì)算結(jié)果小數(shù)點(diǎn)后保留4位數(shù)字。（4）考試時(shí)間為120分鐘。一、（10分）某廠生產(chǎn)A、B兩種產(chǎn)品，1千克原料在甲類設(shè)備上用12小時(shí)可生產(chǎn)3件A，可獲凈利潤(rùn)64元；在乙類設(shè)備上用8小時(shí)可生產(chǎn)4件B，可獲凈利潤(rùn)54元。該廠每天可獲得55千克原料，每天總的勞動(dòng)時(shí)間為480小時(shí)，且甲類設(shè)備每天至多能生產(chǎn)80件A。試為該廠制訂生產(chǎn)計(jì)劃使每天的凈利潤(rùn)最大。1） 以生產(chǎn)A、B產(chǎn)品所用原料的數(shù)

2、量x1、x2（千克）作為決策變量，建立的數(shù)學(xué)規(guī)劃模型是：決策變量：生產(chǎn)A原料x1；生產(chǎn)B原料x2目標(biāo)函數(shù)：y=64*x1+54*x2約束條件：x1+x2 5512*x1+8*x24803*x180x1,x20基本模型：max(y)= 64*x1+54*x2s.t. x1+x2 5512*x1+8*x24803*x180x1,x20c=64 54;A1= 1 1 ; 12 8 ; 3 0;b1=55;480;80; v1=0 0;x,z,ef,out,lag=linprog(-c,A1,b1,v1)lag.ineqlin輸出結(jié)果:x =10.000000004005848 44.99999999

3、3870908z =-3.069999999925403e+003ans = 33.999999998919357 2.500000000140441 0.0000000002784052） 每天的最大凈利潤(rùn)是_3070_元。若要求工人加班以增加勞動(dòng)時(shí)間，則加班費(fèi)最多為每小時(shí)_2.5_元。若A獲利增加到26元/件，應(yīng)否改變生產(chǎn)計(jì)劃？_不變_ c=78 54;A1= 1 1 ; 12 8 ; 3 0;b1=55;480;80; v1=0 0;x,z,ef,out,lag=linprog(-c,A1,b1,v1)x = 9.999999999999400 45.000000000000625z =

4、-3.209999999999987e+003二、(10分) 已知常微分方程組初值問題 試用數(shù)值方法求_ 1.73205_(保留小數(shù)點(diǎn)后5位數(shù)字)。你用的MATLAB命令是_ ode45(ff, ts,y0)_，其精度為_四階_。%待解常微分方程組函數(shù)M文件源程序：function dy=ff (x,y) dy=y(2);-y(2)./x-y(1)*(x.2-0.25)/(x.2);%應(yīng)用歐拉方法和龍格-庫(kù)塔方法求解該常微分方程：ts=pi/2:- pi/12:pi/6; ！步長(zhǎng)必須是可以整除步長(zhǎng)區(qū)間長(zhǎng)度的數(shù)y0=2,-2/pi;x,y=ode45(ff, ts,y0); %龍格-庫(kù)塔方法求數(shù)

5、值解x, y(:,1)輸出結(jié)果：0.523598775598299 1.732050795523993三、 (10分) 已知線性代數(shù)方程組Ax=b, 其中, , ,若方程組右端項(xiàng)有小擾動(dòng)，試根據(jù)誤差估計(jì)式估計(jì)_0.0743_（分別表示原問題的解和右端項(xiàng)小擾動(dòng)后對(duì)應(yīng)的解的變化量）；若取初值，則用高斯-賽德爾迭代法求解Ax=b時(shí)，_(1.7160, 0.3926, -0.1306, 0.1381)_；對(duì)本題而言，此迭代方法是否收斂_是_，原因是_譜半徑(B)=0.397<1_。線性代數(shù)方程組解的誤差分析：故其誤差上限為：A=5 -7 0 1 ;-3 22 6 2 ;5 -1 31 -1 ;2

6、 1 0 23; b=6 3 4 7;db=0 0 0 0.1;d=cond(A,1)*norm(db,1)/norm(b,1)輸出結(jié)果:d =0.074339065208930A=5 -7 0 1 ;-3 22 6 2 ;5 -1 31 -1 ;2 1 0 23;D=diag(diag(A); %從稀疏矩陣A中提取DL=-tril(A,-1); %從稀疏矩陣A中提取LU=-triu(A,1); %從稀疏矩陣A中提取U b=6 3 4 7' %設(shè)定方程組右端項(xiàng)向量x= zeros(4,1); %設(shè)定方程組初始向量 m= inv(D-L)*U;n= inv(D-L)*b; %高斯-賽德爾

7、迭代法for j2=1:5 y=m*(x(:,j2); for i=1:4 x(i,j2+1)=y(i,:)+n(i,:); end endt2=x(:,end) %輸出迭代法最終結(jié)果j2 輸出結(jié)果:t2 = 1.715972347226445 0.392646824062879 -0.130571100623047 0.138061238325401判斂:lamda=eig(inv(D-L)*U)pubanjing=max(abs(lamda)輸出結(jié)果:pubanjing =0.396832340862002四、（20分）炮彈射擊的目標(biāo)為一橢圓形區(qū)域，在X方向半軸長(zhǎng)110m，Y方向半軸長(zhǎng)90

8、m.當(dāng)瞄準(zhǔn)目標(biāo)的中心發(fā)射炮彈時(shí)，在眾多隨機(jī)因素的影響下，彈著點(diǎn)服從以目標(biāo)中心為均值的二維正態(tài)分布，設(shè)彈著點(diǎn)偏差的均方差在X方向和Y方向分別為70m和50m。今測(cè)得一組彈著點(diǎn)的橫縱坐標(biāo)如下：X-6.3-71.665.6-79.2-49.7-81.974.6-47.6-120.856.9Y 28.91.661.7-68-41.3-30.58717.3-17.81.2X100.9479.7-60.1-52.78680.6-42.656.415.2Y-12.639.18532.728.1-9.3-4.55.1-32-9.51) 根據(jù)這組數(shù)據(jù)對(duì)X方向和Y方向的均值和均方差進(jìn)行假設(shè)檢驗(yàn)（設(shè)顯著性水平為0.

9、05）。均值假設(shè)檢驗(yàn):H0:=0;H1:0；x=-6.3 -71.6 65.6 -79.2 -49.7 -81.9 74.9 -47.6 -120.8 56.9 100.9 47 9.7 -60.1 -52.7 86 80.6 -42.6 56.4 15.2;y=28.9 1.6 61.7 -68 -41.3 -30.5 87 17.3 -17.8 1.2 -12.6 39.1 85 32.7 28.1 -9.3 -4.5 5.1 -32 -9.5;h1=ztest(x,0,70)h2=ztest(y,0,50)輸出結(jié)果：h1 =0h2 =0方差假設(shè)檢驗(yàn):H0：2=02；H1：202；x=-6

10、.3 -71.6 65.6 -79.2 -49.7 -81.9 74.9 -47.6 -120.8 56.9 100.9 47 9.7 -60.1 -52.7 86 80.6 -42.6 56.4 15.2;y=28.9 1.6 61.7 -68 -41.3 -30.5 87 17.3 -17.8 1.2 -12.6 39.1 85 32.7 28.1 -9.3 -4.5 5.1 -32 -9.5;function h=ktest(x,s0,alpha,tail)n=length(x);k=(n-1)*var(x)/(s02) %2分布檢驗(yàn)方差if tail=0k1=chi2inv(alpha

11、/2,n-1)k2=chi2inv(1-alpha/2,n-1)if k>=k1&k<=k2 h=0;else h=1;endendif tail=1k0=chi2inv(1-alpha,n-1)if k<=k0 h=0;else h=1;endendif tail=-1k0=chi2inv(alpha,n-1)if k>=k0 h=0;else h=1;endh1=ktest(x,70,0.05,0)h2=ktest(y,50,0.05,0)輸出結(jié)果：h1 =0h2 =02) 根據(jù)這組數(shù)據(jù)給出隨機(jī)變量X和Y相關(guān)系數(shù)的一個(gè)點(diǎn)估計(jì)。相關(guān)系數(shù)點(diǎn)估計(jì)：x=-6.3 -

12、71.6 65.6 -79.2 -49.7 -81.9 74.9 -47.6 -120.8 56.9 100.9 47 9.7 -60.1 -52.7 86 80.6 -42.6 56.4 15.2;y=28.9 1.6 61.7 -68 -41.3 -30.5 87 17.3 -17.8 1.2 -12.6 39.1 85 32.7 28.1 -9.3 -4.5 5.1 -32 -9.5;r=corrcoef(x,y)輸出結(jié)果：r= 0.3134123051021973) 用蒙特卡羅方法求炮彈落在橢圓形區(qū)域內(nèi)的概率（取10000個(gè)數(shù)據(jù)點(diǎn)；請(qǐng)附程序）。%炮彈命中橢圓形區(qū)域概率源程序：a=110

13、;b=90;sx=70;sy=50;r=0.3134123;z=0;n=10000;x=unifrnd(-a,a,1,n);y=unifrnd(-b,b,1,n);for i=1:n if (x(i)2)/(a2)+y(i)2/(b2)<=1 u=exp(-0.5/(1-r2)*(x(i)2/sx2-2*r*x(i)*y(i)/(sx*sy)+y(i)2/sy2); z=z+u; endendP=4*a*b*z/(2*pi*sx*sy*sqrt(1-r2)*n)輸出結(jié)果：P = 0.761272218724379考試課程 數(shù)學(xué)實(shí)驗(yàn) 2005.6.15下午班級(jí) 學(xué)號(hào) 姓名 得分 說(shuō)明（1）

14、第一、二、三題的答案直接填在試題紙上；（2）第四題將數(shù)學(xué)模型、簡(jiǎn)要解題過程和結(jié)果寫在試題紙上；卷面空間不夠時(shí)，請(qǐng)寫在背面；（3）除非特別說(shuō)明，所有計(jì)算結(jié)果小數(shù)點(diǎn)后保留4位數(shù)字。（4）考試時(shí)間為120分鐘。一、（10分）某廠生產(chǎn)A、B兩種產(chǎn)品，1千克原料在甲類設(shè)備上用12小時(shí)可生產(chǎn)3件A，可獲凈利潤(rùn)64元；在乙類設(shè)備上用8小時(shí)可生產(chǎn)4件B，可獲凈利潤(rùn)56元。該廠每天可獲得55千克原料，每天總的勞動(dòng)時(shí)間為480小時(shí)，且甲類設(shè)備每天至多能生產(chǎn)80件A。試為該廠制訂生產(chǎn)計(jì)劃使每天的凈利潤(rùn)最大。1） 以生產(chǎn)A、B產(chǎn)品所用原料的數(shù)量x1、x2（千克）作為決策變量，建立的數(shù)學(xué)規(guī)劃模型是：2） 每月的最大凈利

15、潤(rùn)是_元。若要求工人加班以增加勞動(dòng)時(shí)間，則加班費(fèi)最多為每小時(shí)_元。若A獲利增加到27元/件，應(yīng)否改變生產(chǎn)計(jì)劃？_ 二、(10分) 已知常微分方程組初值問題試用數(shù)值方法求_(保留小數(shù)點(diǎn)后5位數(shù)字)。你用的MATLAB命令是_，其精度為_。 三、(10分) 已知線性代數(shù)方程組Ax=b, 其中, , , 若方程組右端項(xiàng)有小擾動(dòng)，試根據(jù)誤差估計(jì)式估計(jì)_ （分別表示原問題的解和右端項(xiàng)小擾動(dòng)后對(duì)應(yīng)的解的變化量）；若取初值，則用高斯-賽德爾迭代法求解Ax=b時(shí)，_；對(duì)本題而言，此迭代方法是否收斂_，原因是_。四、（20分）炮彈射擊的目標(biāo)為一橢圓形區(qū)域，在X方向半軸長(zhǎng)100m，Y方向半軸長(zhǎng)80m.當(dāng)瞄準(zhǔn)目標(biāo)的

16、中心發(fā)射炮彈時(shí)，在眾多隨機(jī)因素的影響下，彈著點(diǎn)服從以目標(biāo)中心為均值的二維正態(tài)分布，設(shè)彈著點(diǎn)偏差的均方差在X方向和Y方向分別為70m和50m。今測(cè)得一組彈著點(diǎn)的橫縱坐標(biāo)如下：X-6.3-71.665.6-79.2-49.7-81.974.6-47.6-100.856.9Y 28.91.661.7-68-41.3-30.58717.3-17.81.2X100.9479.7-60.1-52.78680.6-42.656.415.2Y-12.639.18532.728.1-9.3-4.55.1-32-9.51） 根據(jù)這組數(shù)據(jù)對(duì)X方向和Y方向的均值和均方差進(jìn)行假設(shè)檢驗(yàn)（設(shè)顯著性水平為0.05）。2） 根

17、據(jù)這組數(shù)據(jù)給出隨機(jī)變量X和Y相關(guān)系數(shù)的一個(gè)點(diǎn)估計(jì)。3） 用蒙特卡羅方法求炮彈落在橢圓形區(qū)域內(nèi)的概率（取10000個(gè)數(shù)據(jù)點(diǎn)；請(qǐng)附程序）?？荚囌n程 數(shù)學(xué)實(shí)驗(yàn) 參考答案與評(píng)分標(biāo)準(zhǔn) 2005.6.15A卷（班級(jí)-姓名-學(xué)號(hào)-得分）一、1） （如果進(jìn)一步要求3x1和4x2為非負(fù)整數(shù)，不扣分） 2）3070元，2.5元；不變二、1.73203(或1.73205), ode23(或ode45)，3級(jí)2階(或5級(jí)4階) （不寫3級(jí)(或5級(jí)）不扣分；個(gè)別同學(xué)可能用其他命令，則結(jié)果相應(yīng)略有變化)三、0.0743, 1.7160, 0.3926, -0.1306, 0.1381, 收斂,譜半徑為0.3968<

18、1（不寫出譜半徑的具體數(shù)值不扣分，但寫錯(cuò)要扣分）四、1）對(duì)均值做的假設(shè)為(X,Y方向相同) X,Y方向均接受H0對(duì)X方向的方差做的假設(shè)為 H0: =4900, H0: 4900（如果做單側(cè)檢驗(yàn)，可不扣分） 接受H0對(duì)Y方向的方差做的假設(shè)為H0: =2500, H0: 2500（如果做單側(cè)檢驗(yàn)，可不扣分）接受H02) 相關(guān)系數(shù)的點(diǎn)估計(jì)為0.313 （用r=corrcoef(x,y)命令）3)大約0.76，結(jié)果具有隨機(jī)性附主要程序示例：%1)2)x=-6.3 -71.6 65.6 -79.2 -49.7 -81.9 74.6 -47.6 -120.8(%B-100.8) 56.9 100.9 47

19、 9.7 -60.1 -52.7 86 80.6 -42.6 56.4 15.2;h1=ztest(x,0,70), %x方向均值檢驗(yàn)y=28.9 1.6 61.7 -68 -41.3 -30.5 87 17.3 -17.8 1.2 -12.6 39.1 85 32.7 28.1 -9.3 -4.5 5.1 -32 -9.5;h2=ztest(y,0,50), %Y方向均值檢驗(yàn)r=corrcoef(x,y) %相關(guān)系數(shù)的點(diǎn)估計(jì)pausen=20;alpha=0.05;sx2=var(x),sx0=70;chi2=(n-1)*sx2/(sx02)chi2alpha=chi2inv(1-alpha

20、,n-1)if chi2<=chi2alpha H0=0else H0=1endsy2=var(y),sy0=50;chi2=(n-1)*sy2/(sy02)chi2alpha=chi2inv(1-alpha,n-1)if chi2<=chi2alpha H0=0else H0=1end%3) a=0.7;b=0.5;m=0;z=0;p=0.313;c=1.1;d=0.9;%A % p=0.311;c=1; d=0.8; %Bn=10000;for i=1:n x=2*rand(1,2)-1; y=0; if x(1)2+x(2)2<=1 y=exp(-0.5/(1-p*p)

21、*(c2*x(1)2/a2+d2*x(2)2/b2-2*p*c*d*x(1)*x(2)/a/b); z=z+y; m=m+1; endendP=4*c*d*z/2/pi/a/b/sqrt(1-p*p)/n,mB卷（班級(jí)-學(xué)號(hào)-姓名）一、1） 2）3160元； 2元；不變二、1.53077 (或1.53073), 龍格-庫(kù)塔方法ode23(或ode45)，3級(jí)2階(或5級(jí)4階)三、0.0826, 2.3416,0.5359,-0.2961,-0.0095, 收斂, 譜半徑為0.3968<1四、1)同A卷。2)0.311; 3) 大約0.69，結(jié)果具有隨機(jī)性評(píng)分標(biāo)準(zhǔn)：一、1)5分(每個(gè)式子一

22、分)，2)前兩空每空2分，最后一空1分。二、第一空6分，后兩空每空2分。三、第一空3分，中間一空5分，最后兩空每空1分。四、1)8分（4個(gè)檢驗(yàn)每個(gè)2分）；2）2分；3）10分（看程序及結(jié)果分析：積分表達(dá)式2分，程序5分，結(jié)果3分）考試課程 數(shù)學(xué)實(shí)驗(yàn) 參考答案與評(píng)分標(biāo)準(zhǔn) 2005.6.15A卷（班級(jí)-姓名-學(xué)號(hào)-得分）一、1） （如果進(jìn)一步要求3x1和4x2為非負(fù)整數(shù)，不扣分） 2）3070元，2.5元；不變二、1.73203(或1.73205), ode23(或ode45)，3級(jí)2階(或5級(jí)4階) （不寫3級(jí)(或5級(jí)）不扣分；個(gè)別同學(xué)可能用其他命令，則結(jié)果相應(yīng)略有變化)三、0.0743, 1.

23、7160, 0.3926, -0.1306, 0.1381, 收斂,譜半徑為0.3968<1（不寫出譜半徑的具體數(shù)值不扣分，但寫錯(cuò)要扣分）四、1）對(duì)均值做的假設(shè)為(X,Y方向相同) X,Y方向均接受H0對(duì)X方向的方差做的假設(shè)為 H0: =4900, H0: 4900（如果做單側(cè)檢驗(yàn)，可不扣分） 接受H0對(duì)Y方向的方差做的假設(shè)為H0: =2500, H0: 2500（如果做單側(cè)檢驗(yàn)，可不扣分）接受H02) 相關(guān)系數(shù)的點(diǎn)估計(jì)為0.313 （用r=corrcoef(x,y)命令）3)大約0.76，結(jié)果具有隨機(jī)性附主要程序示例：%1)2)x=-6.3 -71.6 65.6 -79.2 -49.7 -81.9 74.6 -47.6 -120.8(%B-100.8) 56.9 100.9 47 9.7 -60.1 -52.7 86 80.6 -42.6 56.4 15.2;h1=ztest(x,0,70), %x方向均值檢驗(yàn)y=28.9 1.6 61.7 -68 -41.3 -30.5 87 17.3 -17.8 1.2 -12.6 39.1 85 32.7 28.1 -9.3 -4.5 5.1 -32 -9.5;h2=ztest(y,0,50), %Y方向均值檢驗(yàn)r=corrcoef(x,y) %相關(guān)系數(shù)的點(diǎn)估計(jì)paus