王樹然李天磊姜磊

1、2011高教社杯全國大學(xué)生數(shù)學(xué)建模競賽承 諾 書我們仔細(xì)閱讀了中國大學(xué)生數(shù)學(xué)建模競賽的競賽規(guī)則.我們完全明白，在競賽開始后參賽隊員不能以任何方式（包括電話、電子郵件、網(wǎng)上咨詢等）與隊外的任何人（包括指導(dǎo)教師）研究、討論與賽題有關(guān)的問題。我們知道，抄襲別人的成果是違反競賽規(guī)則的, 如果引用別人的成果或其他公開的資料（包括網(wǎng)上查到的資料），必須按照規(guī)定的參考文獻(xiàn)的表述方式在正文引用處和參考文獻(xiàn)中明確列出。我們鄭重承諾，嚴(yán)格遵守競賽規(guī)則，以保證競賽的公正、公平性。如有違反競賽規(guī)則的行為，我們將受到嚴(yán)肅處理。我們參賽選擇的題號是（從A/B/C/D中選擇一項填寫）： 我們的參賽報名號為（如果賽區(qū)設(shè)置報名

2、號的話）： 所屬學(xué)校（請?zhí)顚懲暾娜?參賽隊員 (打印并簽名) ：1. 2. 3. 指導(dǎo)教師或指導(dǎo)教師組負(fù)責(zé)人 (打印并簽名)： 日期： 年 月 日賽區(qū)評閱編號（由賽區(qū)組委會評閱前進(jìn)行編號）：2011高教社杯全國大學(xué)生數(shù)學(xué)建模競賽編 號 專 用 頁賽區(qū)評閱編號（由賽區(qū)組委會評閱前進(jìn)行編號）：賽區(qū)評閱記錄（可供賽區(qū)評閱時使用）：評閱人評分備注全國統(tǒng)一編號（由賽區(qū)組委會送交全國前編號）：全國評閱編號（由全國組委會評閱前進(jìn)行編號）：降落傘的選購模型摘要 本模型研究的是空投物資所用降落傘的選購問題，在滿足降落傘落地時速度不超過的前提條件下將物資全部安全空投到救災(zāi)區(qū)，要求所購買的降落傘的費用最少。

3、為方便研究，我們將傘面、繩索、物資看成整體，忽略傘和繩索的質(zhì)量。首先確定不同半徑的降落傘的價格，運用勾股定理計算單根繩長，根據(jù)已知傘面費用、半徑、繩索單價以及固定費用得出不同半徑傘的價格。分析物體的運動狀態(tài)，根據(jù)假設(shè)降落傘豎直向上只受空氣阻力和自身的重力作用，水平方向不受力。通過對降落傘受力分析建立起了高度與時間的方程，根據(jù)（高度）和（時間）的一階導(dǎo)得出的圖像可知，物體做向下加速度減小的加速運動，到開始向下做勻速運動，根據(jù)牛頓第二定律計算空氣阻力系數(shù)，首先得出（高度）和（時間）的二階導(dǎo)函數(shù)，利用加速度以及結(jié)合表二運用Matlab軟件擬合出空氣阻力系數(shù)值，根據(jù)已知物資在落地時速度必須小于，研究圖

4、像可知物資在空氣阻力和自身重力的作用下在落地之前達(dá)到平衡狀態(tài)，根據(jù)牛頓第二定律得知即就是，不同不同，這樣就可以根據(jù)得出對應(yīng)不同半徑傘的最大載重量。下面就是解決最優(yōu)解的問題，通過已知，空投物資必須大于等于以及降落傘的個數(shù)必須是正整數(shù)，得出目標(biāo)函數(shù)的約束條件，根據(jù)每種傘的價格和個數(shù)列出目標(biāo)函數(shù)。 最后運用Lingo軟件求解得: .即購買半徑為的降落傘所需費用最少，最大載重量為，最低費用為。關(guān)鍵字：空氣阻力 、勾股定理、線性規(guī)劃、數(shù)據(jù)擬合一、問題重述為向災(zāi)區(qū)空投救災(zāi)物資共2000kg，需選購一些降落傘。已知空投高度為500m，要求降落傘落地時的速度不能超過20m/s。降落傘面為半徑r的半球面，用每根

5、長 L, 共16根繩索連接的載重m的物體位于球心正下方球面處，每個降落傘的價格由三部分組成。傘面費用C1由傘的半徑r決定，見表1；繩索費用C2由繩索總長度及單價4元/米決定；固定費用C3為200元。 2 2.5 3 3.5 4 費用（元） 65 170 350 660 1000 表1降落傘在降落過程中受到重力作用外還受到的空氣阻力，可以認(rèn)為與降落速度和傘的受力面積的乘積成正比。為了確定阻力系數(shù)，用半徑r=3m、載重m=300kg的降落傘從500m高度作降落試驗，測得各時刻的高度 ，見表2。時刻 0 3 6 9 12 15 18 21 24 27 30 高度 500 470 425 372 31

6、7 264 215 160 108 55 1 表2試根據(jù)以上條件確定降落傘的選購方案，即共需多少個，每個傘的半徑多大（在表1中選擇），在滿足空投要求的條件下，使費用最低。二、模型假設(shè)1假設(shè)物體下落時初速度為0；2假設(shè)降落傘在空投的瞬間打開；3假設(shè)傘和繩子的重量忽略不計；42000kg的物資可以任意分割；5假設(shè)降落傘水平方向不受力；6降落傘在下落的過程中只考慮空氣阻力對其的影響；三、符號說明: 購買半徑為r的降落傘數(shù)目；：空氣的阻力；：空氣的阻力系數(shù)；：降落傘的傘面費用；：降落傘懸掛物資的繩索費用；：降落傘的固定費用；：降落傘所用的總費用；：降落傘的半徑；：降落傘的最大載重量；：降落傘離地的高度

7、；：降落傘的下落的位移；：降落傘每根繩的長度；：降落傘下落的加速度；：重力加速度；：降落傘的傘面面積；四、問題分析1分析傘面費用：由題意可知每個傘的價格，有三部分構(gòu)成：傘面費用、繩索費用、固定費用，傘面費用由傘的半徑?jīng)Q定，題目中已經(jīng)給出傘面的具體費用，繩索費用由繩索的長度及單價所決定，根據(jù)題意“載重m的物體位于球心正下方球面處”（如圖一所示），根據(jù)勾股定理得繩索的長度，單價已經(jīng)得知，每個傘由16根繩索構(gòu)成，所以繩索費用，費用是題目中給出的固定費用200。已知題目中已經(jīng)給出每種傘形的半徑，所以每種傘的單價為定值。2分析物體的運動狀態(tài)：我們將降落傘、繩索和物資看成一整體，對它進(jìn)行受力分析（如圖二所

8、示），根據(jù)牛頓第二定律：，是物體的重量，是物體下落時遇到向上的阻力，表示降落傘的面積，所以,所以可以得出：由于需要降落傘對它的運動速度進(jìn)行控制，并且降落傘在運動的過程中受到空氣阻力的影響，所以可以得出物資在空中不是自由落體運動，物資在空投的時候它的運動狀態(tài)可以分為三種：2.1：物資在投出的時候阻力大于重力，它的運動狀態(tài)為向上變速運動；2.2：物資在投出的時候阻力等于重力，它的運動狀態(tài)為向上勻速運動或向下勻速運動或靜止不動；2.3：物資在投出的時候阻力小于重力，它的運動狀態(tài)為向下變速運動；根據(jù)題意可判斷出這三種都有可能，要具體得知是哪種,就需要研究物資從投出到落地一瞬間的運動狀態(tài)，題意中已給出物

9、資的時間與位移的數(shù)據(jù)，可由得出它的一階導(dǎo)速度的圖像，根據(jù)物資的速度運動狀態(tài)圖像（如圖三所示）可得出第一種運動狀態(tài)不可能發(fā)生，空投物資是需要將物資從高空向下投到地面，所以，具體就是分析物資向下運動的狀態(tài)是向下勻速運動還是向下變加速運動，根據(jù)分析可以得出物資先是做的是向下的加速度減小的加速運動，到差不多快10秒的時間的時候，物體開始做勻速運動。物體下落速度變化圖現(xiàn)在知道了物體的運動狀態(tài)，可以果斷排除前面兩種情況，即就是物體在空投后是做向下加速度減小的變加速運動。3分析空氣阻力系數(shù)：空氣阻力系數(shù)決定著降落傘的最大載重量，根據(jù)牛頓第二定律可以得出，所以要想知道空氣阻力系數(shù)，就必須得出物資下落的加速度，

10、現(xiàn)在從題目中可以已知（時間）和 (高度)的對應(yīng)關(guān)系，我們在數(shù)學(xué)中學(xué)習(xí)過和的一階導(dǎo)是（速度），而速度的一階導(dǎo)就是(加速度)，首先根據(jù)和的對應(yīng)數(shù)據(jù)再結(jié)合題目中的實驗說明（半徑、載重的降落傘從高度作降落試驗）利用Matlab軟件計算出H=1/2*g/k/pi/r2*m*t+1/4*g/k2/pi2/r4*m2*exp(-2*k*pi*r2/m*t)的函數(shù)關(guān)系，再得出（高度）和（時間）的二階導(dǎo)，利用Matlab軟件擬合出空氣阻力系數(shù)的值。4：分析降落傘的最大承載重量：根據(jù)題目給出的條件（已知空投高度為，要求降落傘落地時的速度不能超過），要想在物資落地的一瞬間速度不超過，就需要降落傘的承載重量不能太大，

11、那么這個就是一個零界狀態(tài)，根據(jù)物體下落速度變化圖可以得出物體下落時間在10s的時候（加速度）=0，所以說從10s開始物體做的是向下的勻速運動，根據(jù)牛頓第二定律，即就是，得出解得，表1中已經(jīng)給出每種降落傘的半徑，利用Matlab軟件解得每種半徑的降落傘的最大承重量。5利用線性規(guī)劃，使用lingo軟件得出最優(yōu)解：前面已經(jīng)解出每種降落傘的價格，目標(biāo)函數(shù)最低費用就是每種傘的 (價格)乘以它的(個數(shù))相加，再根據(jù)題目給出的物資共得出，降落傘的個數(shù)是正整數(shù)且空投總物資重量必須大于等于，所以得出約束條件，再使用Lingo軟件求得最優(yōu)解。圖一 圖二五、建模與求解1首先確定阻力系數(shù)k為了方便對物資進(jìn)行受力分析，

12、我們把降落傘和物資看作一個整體如圖二。由假設(shè)5可知物體A只受到豎直向上的空氣阻力和豎直向下的重力作用。又由題可知空氣阻力與降落速度v和傘的受力面積S的乘積成正比。則物體A在豎直方向上受到的合外力為：由運動學(xué)方程：得由物體位移H和時間(t)的二次微分等于加速度建立方程得: 用MATLAB解微分方程得：(程序見附錄1)則題目已經(jīng)給t-h數(shù)據(jù)為：時刻t（s）036912151821242730高度h（m）500470425372317264215160108551時間高度對應(yīng)關(guān)系表對給定的數(shù)據(jù)以為擬合函數(shù)進(jìn)行擬合，r=3m,m=300kg,g=9.8,得出k=2.9377 。（程序見附錄2）2求解最

13、大承載量用速度對時間的微分等于加速度，且v0=0建立方程組得：用MATLAB解得(程序見附錄3)因為降落傘在下落過程中其質(zhì)量是不變的，所以我們把關(guān)系式中t看做一個定值，則關(guān)于m的方程為從上式我們可以知道是關(guān)于m的單調(diào)遞增函數(shù):證明過程如下：由數(shù)學(xué)知識可知：函數(shù)的一階導(dǎo)數(shù)大于零，則原函數(shù)是單調(diào)遞增的。一階導(dǎo)數(shù)小于零，則原函數(shù)是單調(diào)遞減的。對求一階導(dǎo)數(shù)得：由上式分析可知無法確定其是否大于零，在對其求二階導(dǎo)數(shù)為：則一階導(dǎo)數(shù)為單調(diào)遞減函數(shù)，當(dāng)m趨近于無窮大時對一階導(dǎo)數(shù)求極限可知由此可得：則原函數(shù)是單調(diào)遞增函數(shù)，即速度v和m是成正比關(guān)系的。用Matlab畫出V-t圖，（程序見附表4）識圖可知物資在下過程

14、中：開始的加速度減小的加速運動，隨著速度的增大, 由 得出阻力不斷增大。當(dāng)，物質(zhì)做勻速運動。物體下落速度隨時間的變化圖因此，求最大載重量取傘在下降到地面的瞬間達(dá)到最大速度,此時m，由方程調(diào)用MATLAB分別解得半徑為r的降落傘在滿足空投條件下的最大載重量如下表：（程序見附錄5）降落傘在滿足空投條件下的最大承重量表3線性規(guī)劃求解數(shù)量和費用由分析可知每種傘的單價：由題可知為:r（m）22.533.54C1（元）651703506601000不同半徑的降落傘的傘面費用表為： 為固定值即： 由以上數(shù)據(jù)求得每種傘的單價見下表：則購買每把不同半徑的降落傘的各需總費用C如下:r(m)22.533.54總費用

15、（元）446596.3821.51176.81562.0表3我們設(shè)每種傘分別取n2,n2.5,n3,n3.5,n4個,則其目標(biāo)函數(shù)為：min=446*n1+596.3*n2+821.5*n3+1176.8*n4+1562*n5;150*n1+235*n2+339*n3+461*n4+602*n5>=2000;n1>=0;n2>=0;n3>=0;n4>=0;n5>=0;gin(n1);gin(n2);gin(n3);gin(n4);gin(n5);對其進(jìn)行優(yōu)化求解C的最小值，就是所需的最小費用。用LINGO求解得（程序見附件6） n2=0,n2.5=0,n3=

16、6,n3.5=0,n4=0購買半徑為3m的傘六個所需費用最少，最少總費用C=6*821.5=4929元。六、模型的評價優(yōu)點：1 本模型的求解過程中用了大量的電腦軟件，使得計算數(shù)據(jù)更加可靠，方法簡便，。2當(dāng)降落傘的半徑跟題目中所給的半徑相同時，這個模型同樣適用。缺點：1本模型太過理想化2. 本模型假設(shè)的是理想的無風(fēng)狀態(tài)下并沒有過多的考慮外界的其他因素，3. 本模型未考慮降落傘打開的時間，只是假設(shè)在空投的瞬間打開，與實際情況并不相符。七、參考文獻(xiàn)1薛定宇，<<高等應(yīng)用數(shù)學(xué)問題的Matlab求解>>，北京：清華大學(xué)出版社，2008、10、1234附錄一H=dsolve(

17、9;m*D2H+k*S*DH=m*g','H(0)=0,DH(0)=0','t')得：g/k2/S2*m2*exp(-k*S/m*t)+g/k/S*m*t-1/k2/S2*m2* 附錄二擬合k 利用matlab軟件，首先建立一個m文件： function H=fl(k,t);m=300,g=9.8,r=3H=1/2*g/k/pi/r2*m*t+1/4*g/k2/pi2/r4*m2*exp(-2*k*pi*r2/m*t) 然后在Command window中輸入：h=0 30 75 128 183 236 285 340 392 445 499;>&

18、gt; t=0 3 6 9 12 15 18 21 24 27 30;>> k=lsqcurvefit('f1',2,t,h)得出：m =300g =9.8000r =3H =Columns 1 through 6 68.9547 100.2389 163.1533 236.2754 312.6916 390.1710 Columns 7 through 11 467.9934 545.9266 623.8955 701.8759 779.8601m = 300g = 9.8000r = 3H = Columns 1 through 6 68.9547 100.23

19、89 163.1533 236.2754 312.6916 390.1710 Columns 7 through 11 467.9934 545.9266 623.8955 701.8759 779.8601m =300g =9.8000r =3H = Columns 1 through 6 36.3091 64.2311 114.7885 170.1093 226.4325 282.9666 Columns 7 through 11 339.5450 396.1329 452.7227 509.3129 565.9031m = 300g = 9.8000r = 3H = Columns 1

20、through 6 36.3091 64.2311 114.7885 170.1093 226.4325 282.9665 Columns 7 through 11 339.5450 396.1329 452.7226 509.3128 565.9031m = 300g = 9.8000r = 3H = Columns 1 through 6 27.7320 54.1175 99.6966 148.5016 197.8487 247.2869 Columns 7 through 11 296.7404 346.1966 395.6531 445.1097 494.5664m = 300g =

21、9.8000r = 3H = Columns 1 through 6 27.7320 54.1175 99.6966 148.5016 197.8487 247.2869 Columns 7 through 11 296.7404 346.1966 395.6531 445.1097 494.5664m = 300g = 9.8000r = 3H = Columns 1 through 6 26.6293 52.7783 97.6260 145.5035 193.8719 242.3200 Columns 7 through 11 290.7809 339.2439 387.7072 436.

22、1706 484.6340m = 300g = 9.8000r = 3H = Columns 1 through 6 26.6293 52.7783 97.6260 145.5035 193.8719 242.3200 Columns 7 through 11 290.7809 339.2438 387.7072 436.1706 484.6340m = 300g = 9.8000r = 3H = Columns 1 through 6 26.6101 52.7550 97.5897 145.4509 193.8022 242.2328 Columns 7 through 11 290.676

23、3 339.1219 387.5678 436.0138 484.4597m = 300g = 9.8000r = 3H = Columns 1 through 6 26.6101 52.7550 97.5897 145.4509 193.8022 242.2328 Columns 7 through 11 290.6763 339.1219 387.5678 436.0138 484.4597m = 300g = 9.8000r = 3H = Columns 1 through 6 26.6102 52.7550 97.5898 145.4511 193.8024 242.2331 Colu

24、mns 7 through 11 290.6766 339.1222 387.5682 436.0142 484.4602m = 300g = 9.8000r = 3H = Columns 1 through 6 26.6102 52.7550 97.5898 145.4510 193.8024 242.2331 Columns 7 through 11 290.6766 339.1222 387.5682 436.0142 484.4602Optimization terminated: relative function value changing by less than OPTION

25、S.TolFun.k = 2.9377附錄三 v=dsolve('m*Dv+k*S*v-m*g=0','v(0)=0','t')得：附錄四畫出V-t圖t=0:3:30;g=9.8;m=300;s=2*pi*3*3;k=3.2145;v=g/k/s*m-exp(-k*s/m*t)*g/k/s*m;plot(t,v) 附錄五 求解每個降落傘的最大載重量，利用matlab軟件先建立一個m文件：H=1/2*g/k/pi/r2*m*t+1/4*g/k2/pi2/r4*m2*exp(-2*k*pi*r2/m*t);然后再 Command window 中輸入：>> k=2.9377;>> g=9.8;>> for r=2:0.5:4m=k*40*pi*r2/gend得出：m = 150.6785m = 235.4351m = 339.0266m = 461.4528m = 602.