數(shù)學(xué)實(shí)驗(yàn)作業(yè)

1、精選優(yōu)質(zhì)文檔-傾情為你奉上精選優(yōu)質(zhì)文檔-傾情為你奉上專心-專注-專業(yè)專心-專注-專業(yè)精選優(yōu)質(zhì)文檔-傾情為你奉上專心-專注-專業(yè) 數(shù)學(xué)實(shí)驗(yàn)作業(yè) 學(xué)院：理學(xué)院 班級(jí)：統(tǒng)計(jì)11-1 姓名：吳 學(xué)號(hào)：6 實(shí)驗(yàn)1 一、問(wèn)題的提出已知方程組，其中,定義為通過(guò)迭代法求解方程組。（1）選取不同的初始向量，和不同的方程組右端向量，給定迭代誤差要求，用雅克比和高斯賽德爾迭代法計(jì)算，觀察得到的迭代向量序列是否收斂？（2）取定右端向量和初始向量，將的主對(duì)角元素成倍增長(zhǎng)若干次，非主對(duì)角元素不變，每次用雅克比迭代法計(jì)算，要求迭代誤差滿足，比較收斂速度。二、問(wèn)題分析問(wèn)題（1）要求針對(duì)不同的初始向量，和不同的方程組右端向量，

2、在迭代誤差一定的情況下，分別用用雅克比和高斯賽德爾迭代法計(jì)算方程組，并分析迭代向量序列的收斂性與迭代次數(shù)。問(wèn)題（2）要求在右端向量和初始向量一定的條件下，將的主對(duì)角元素成倍增長(zhǎng)若干次，非主對(duì)角元素不變。在迭代誤差滿足的條件下，用雅克比迭代法比較不同的收斂速度。三、模型建立對(duì)于一般的線性方程組，假設(shè)，雅克比迭代公式是如果將分解為,，迭代公式等價(jià)于如下的矩陣形式：或類似地，線性方程組的高斯賽德爾迭代公式是：等價(jià)于如下的矩陣形式：依據(jù)分析可知問(wèn)題（2）要求迭代誤差滿足，在此，不妨問(wèn)題（1）也采用相同的迭代誤差。針對(duì)問(wèn)題（1），出于簡(jiǎn)化模型的目的，不妨初始向量分別取和，方程組右端向量也分別取和。針對(duì)問(wèn)

3、題（2），出于簡(jiǎn)化模型的目的，分別取原的主對(duì)角元素的1到5倍，初始向量和方程組右端向量分別取和，在迭代誤差滿足的條件下，比較五次得到的迭代結(jié)果進(jìn)行分析。四、模型求解針對(duì)問(wèn)題（1），依據(jù)模型在初始向量為，方程組右端向量為時(shí)編寫如下程序：A=zeros(20);for i=1:20 for j=1:20 if i=j A(i,j)=3; end if i=j-1 A(i,j)=-1/2; end if j=i-1 A(i,j)=-1/2; end if i=j-2 A(i,j)=-1/4; end if j=i-2 A(i,j)=-1/4; end endendL=-tril(A,-1);U=-t

4、riu(A,1);D=diag(diag(A);b=ones(20,1); %右端向量以下的程序?qū)Υ顺鲎鱿鄳?yīng)的改動(dòng)%雅克比B1=D(L+U);f1=Db;x0=zeros(20,1); %初始向量以下的程序?qū)Υ顺鲎鱿鄳?yīng)的改動(dòng)i=1;x=B1*x0+f1; while norm(x-x0,inf)=1e-6 x0=x; %設(shè)置迭代精度為10e-5 x=B1*x0+f1; i=i+1;endx i%高斯賽德爾B2=(D-L)U;f2=(D-L)b;y0=zeros(20,1); %初始向量以下的程序?qū)Υ顺鲎鱿鄳?yīng)的改動(dòng)j=1;y=B2*y0+f2; while norm(y-y0,inf)=1e-6

5、 y0=y; %設(shè)置迭代精度為10e-5 y=B2*y0+f2; j=j+1;endy j運(yùn)行程序得到雅克比算法為了達(dá)到的迭代精度，迭代次數(shù)為20。高斯賽德爾算法為了達(dá)到的迭代精度，迭代次數(shù)為13。其結(jié)果為：0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.改變初始向量為，保持方程組右端向量不變?yōu)?，得到雅克比算法為了達(dá)到的迭代精度，迭代次數(shù)為19。高斯賽德爾算法為了達(dá)到的迭代精度，迭代次數(shù)為13。保持初始向量為不變，改變方程組右端向量為時(shí)，得到雅克比算法為了達(dá)到的迭代精度，迭代次數(shù)為1。高斯賽德爾算法為了達(dá)到的迭代精度，迭代次數(shù)為1。從結(jié)果看該取值很有可能使得方程

6、組出現(xiàn)了病態(tài)解。同時(shí)改變初始向量（）和方程組右端向量（）時(shí)，得到雅克比算法為了達(dá)到的迭代精度，迭代次數(shù)為20。高斯賽德爾算法為了達(dá)到的迭代精度，迭代次數(shù)為14。當(dāng)初始向量和方程組右端向量同時(shí)為時(shí)，方程求解的出現(xiàn)了特殊情況。因此在方程組右端向量為時(shí)，改變初始向量為，得到克比算法為了達(dá)到的迭代精度，迭代次數(shù)為21。高斯賽德爾算法為了達(dá)到的迭代精度，迭代次數(shù)為14。綜合上述程序運(yùn)行結(jié)果，可以看出，在相同的精度要求下，高斯賽德爾的收斂性更好一些。針對(duì)問(wèn)題（2），依據(jù)模型編寫如下程序：A=zeros(20);for i=1:20 for j=1:20 if i=j A(i,j)=3; end if i=

7、j-1 A(i,j)=-1/2; end if j=i-1 A(i,j)=-1/2; end if i=j-2 A(i,j)=-1/4; end if j=i-2 A(i,j)=-1/4; end endendT=cell(1,5);for k=1:5 for i=1:20 for j=1:20 if i=j A(i,j)=3*k; T1,k=A; end end endend %取原對(duì)角元素的1到5倍for i=1:5 L=-tril(T1,i,-1); U=-triu(T1,i,1); D=diag(diag(T1,i); b=ones(20,1);B=D(L+U);f=Db;x0=zer

8、os(20,1); i=1;x=B*x0+f; while norm(x-x0,inf)=1e-6 x0=x; %設(shè)置迭代精度為10e-5 x=B*x0+f; i=i+1;endx iend 運(yùn)行程序得到在迭代誤差滿足的要求下，原系數(shù)矩陣的主對(duì)角元素的1到5倍的迭代次數(shù)分別為20,10,8,7,6。由此可以看出當(dāng)?shù)炔蛔儠r(shí)，主對(duì)角元素成倍增大，達(dá)到精度所需的迭代次數(shù)逐漸下降。五、模型的評(píng)價(jià)與改進(jìn) 該模型充分運(yùn)用了分解，充分體現(xiàn)了雅克比迭代與高斯賽德爾迭代的特征。在相互對(duì)照的模式下，體現(xiàn)了雅克比迭代與高斯賽德爾迭代的收斂速度。也體現(xiàn)出了在系數(shù)矩陣的主對(duì)角元素增大的情況下雅克比迭代的收斂性。

9、實(shí)驗(yàn)2 一、問(wèn)題的提出設(shè)國(guó)民經(jīng)濟(jì)由農(nóng)業(yè)、制造業(yè)和服務(wù)業(yè)三個(gè)部門構(gòu)成，已知某年他們之間的投入產(chǎn)出關(guān)系、外部需求、初始投入等如表1所示。根據(jù)表回答下列問(wèn)題：農(nóng)業(yè)制造業(yè)服務(wù)業(yè)外部需求總產(chǎn)出農(nóng)造業(yè)301045115200服務(wù)業(yè)2060070150初始投入3511075總投入100200150表1：國(guó)民經(jīng)濟(jì)三個(gè)部門之間的投入產(chǎn)出表（1）如果今年對(duì)農(nóng)業(yè)、制造業(yè)和服務(wù)業(yè)的外部需求分別為50,150,100億元，問(wèn)這三個(gè)部門的總產(chǎn)出應(yīng)分別為多少？（2）如果三個(gè)部門的外部需求分別增加一個(gè)單位，問(wèn)他們的總產(chǎn)出應(yīng)分別增加多少？二、問(wèn)題分析題目中給出了農(nóng)業(yè)、制造業(yè)、服務(wù)業(yè)之間的投入產(chǎn)出關(guān)系、

10、外部需求、初始投入等。描述了國(guó)民經(jīng)濟(jì)三個(gè)部門之間的生產(chǎn)消耗和投入產(chǎn)出的數(shù)量關(guān)系。要求建立投入產(chǎn)出的數(shù)學(xué)模型，并運(yùn)用模型求出（1）農(nóng)業(yè)、制造業(yè)和服務(wù)業(yè)的外部需求分別為50,150,100億元時(shí)三個(gè)部門的總產(chǎn)出（2）三個(gè)部門的外部需求分別增加一個(gè)單位時(shí)他們的總產(chǎn)出應(yīng)分別增加多少。三、模型建立首先，依據(jù)題目所給的數(shù)據(jù)可以求出消耗系數(shù)（第個(gè)部門1單位的產(chǎn)出對(duì)第個(gè)部門的世界消耗量） （1）在技術(shù)水平?jīng)]有明顯提高的情況下，假設(shè)直接消耗系數(shù)不變，在這個(gè)假設(shè)下建立投入產(chǎn)出的數(shù)學(xué)模型。設(shè)有n個(gè)部門，記一定時(shí)期內(nèi)第個(gè)部門的總產(chǎn)出為，其中對(duì)第個(gè)部門的投入為，外部需求為，則 （2）由于每個(gè)部門的總產(chǎn)出等于總投入，所以

11、可以將（1）式中的視為第列的總投入，由（1）、（2）式得到 （3）記直接消耗系數(shù)矩陣，產(chǎn)出向量，需求向量，則式（3）可以寫作（其中是單位矩陣）： 或 從方程組看得出，表明總產(chǎn)出對(duì)外部需求是線性的，所以當(dāng)增加一個(gè)單位（記作）時(shí)的增量為。四、模型求解當(dāng)農(nóng)業(yè)、制造業(yè)和服務(wù)業(yè)的外部需求分別為50,150,100億元時(shí)，運(yùn)用建立的模型編寫以下程序：a=15 20 30;30 10 45;20 60 0;b=100 200 150;for i=1:3 for j=1:3 A(i,j)=a(i,j)/b(1,j); endend d=50;150;100;B=eye(3)-A;x=Bd運(yùn)行得到當(dāng)農(nóng)業(yè)、制造業(yè)

12、和服務(wù)業(yè)的外部需求分別為50,150,100億元時(shí)三個(gè)部門的總產(chǎn)出分別為139.2801，267.6056，208.1377。三個(gè)部門的外部需求分別增加一個(gè)單位時(shí)，承接上面程序用矩陣求逆命令計(jì)算：dx=inv(B)得到三個(gè)部門的外部需求分別增加一個(gè)單位時(shí)他們的總產(chǎn)出應(yīng)分別增加1.3459，0.5634，0.4382。五、模型的評(píng)價(jià)與改進(jìn) 該模型中充分運(yùn)用了題目中已知的數(shù)據(jù)，考慮到消耗系數(shù)在投入產(chǎn)出關(guān)系中的作用，并以消耗系數(shù)為起點(diǎn)建立了相應(yīng)的模型。合理分析了矩陣元素所代表的實(shí)際經(jīng)濟(jì)意義，恰當(dāng)?shù)牟捎昧四婢仃囀沟眠\(yùn)算大大簡(jiǎn)化。 實(shí)驗(yàn)3 一、問(wèn)題的提出有三個(gè)節(jié)點(diǎn)的鋼架結(jié)構(gòu)如圖1所示，點(diǎn)1收到100kg

13、外力的作用，點(diǎn)2是固定支點(diǎn)，點(diǎn)3是滑動(dòng)支點(diǎn)。利用力的平衡原理建立模型，討論外力變化1kg時(shí)對(duì)各個(gè)力的影響，求出力。圖1：鋼架結(jié)構(gòu)（左）和三個(gè)支點(diǎn)的受力分析圖（右）二、問(wèn)題分析 題中給出了三個(gè)節(jié)點(diǎn)的鋼架結(jié)構(gòu)，其中點(diǎn)2是固定支點(diǎn)，點(diǎn)3是滑動(dòng)支點(diǎn)。要求利用力的平衡原理建立模型求出點(diǎn)1收到100kg外力作用時(shí)，力的值，并討論外力變化1kg時(shí)對(duì)各個(gè)力的影響。三、模型建立由于此鋼架結(jié)構(gòu)構(gòu)成了一個(gè)三角形，因三角形是一種穩(wěn)定的結(jié)構(gòu)，故知其三邊不會(huì)發(fā)生垂直于桿方向的形變，即三邊只會(huì)在沿桿的方向上產(chǎn)生力的作用?？紤]到整個(gè)系統(tǒng)是處于平衡的，由力學(xué)知識(shí)知，系統(tǒng)在水平方向上與豎直方向上所受合力為零。對(duì)于每一個(gè)分點(diǎn)受力分

14、析也可得到各點(diǎn)所受合力為零。點(diǎn)1處受到一個(gè)豎直向下100kg的力，并設(shè)為力F，由點(diǎn)2指向點(diǎn)1由桿1產(chǎn)生的力F1，由點(diǎn)3指向點(diǎn)1由桿3產(chǎn)生的力F3，并且，F(xiàn)，F(xiàn)1，F(xiàn)3三個(gè)力的合力為0；點(diǎn)2處，由于點(diǎn)2是一個(gè)固定點(diǎn)，故可以產(chǎn)生一水平方向上的力H2，和一個(gè)豎直向上的力V2，以及一個(gè)由點(diǎn)1指向點(diǎn)2的力F1，水平向有的力F2，并且H2，V2，F(xiàn)2，F(xiàn)1四個(gè)力的合力為0；點(diǎn)3處，由于點(diǎn)3是一個(gè)在水平方向可移動(dòng)的點(diǎn)，故點(diǎn)3處只能產(chǎn)生一個(gè)豎直向上的力V3，并且還受到水平向左的力F2，由點(diǎn)1指向點(diǎn)3的力F3，其中，V3，F(xiàn)2，F(xiàn)3合力為0。其中，F(xiàn)1,F1，F(xiàn)2,F2，F(xiàn)3,F3分別是由桿1，桿2，桿3的內(nèi)力

15、，大小相等，方向相反。其受力分析如圖1所示。由點(diǎn)1的受力分析可得（以下分析式中力均為矢量）：由點(diǎn)2的受力分析可得：由點(diǎn)3的受力分析可得： 將6個(gè)矢量轉(zhuǎn)化為標(biāo)量，并將標(biāo)量式轉(zhuǎn)化成線性方程組Ax=b的形式如下：其中， ， 。解得：x=Ab四、模型求解依據(jù)上述建立的模型編寫程序如下：a=sqrt(3);A=1 0 a 0 0 0;0 -2 1 0 0 0;a 0 -1 0 0 0;-a 2 0 2 0 0;-1 0 0 0 2 0;0 0 a 0 0 2;b=200 0 0 0 0 0;x=Ab.運(yùn)行程序得到結(jié)果：的值分別為50.0000，43.3013，86.6025， 0，25.0000，-75

16、.0000。 改變初始值，當(dāng)力改為99kg時(shí)，的值分別為49.5000 42.8683 85.7365 0 24.7500 -74.2500；當(dāng)力改為101kg時(shí)的值分別為50.5000 43.7343 87.4686 0 25.2500 -75.7500?？梢缘贸霎?dāng)F增加1kg x的值的變化為：0.5000 0.43301 0.8660 0 0.2500 0.7500當(dāng)F減少1kg x的值的變化為：-0.5000 -0.4330 -0.8660 0 -0.2500 -0.7500。力F由90kg變化到110kg時(shí)各個(gè)力的大小的變化分別如表2和圖2所示：F9091929394959697989

17、9100F14545.54646.54747.54848.54949.550F238.97139.40439.83740.2740.70341.13641.56942.00242.43542.86843.301F377.94278.80879.67480.5481.40682.27283.13884.00484.8785.73786.603H200000000000V222.522.752323.2523.523.752424.2524.524.7525V367.568.256969.7570.571.257272.7573.574.2575dF100101102103104105106107

18、108109110F15050.55151.55252.55353.55454.555F243.30143.73444.16744.645.03345.46645.89946.33246.76547.19847.631F386.60387.46988.33589.20190.06790.93391.79992.66593.53194.39795.263H200000000000V22525.2525.525.752626.2526.526.752727.2527.5V37575.7576.577.257878.7579.580.258181.7582.5表2：F由90kg變化到110kg時(shí)各個(gè)

19、力的大小圖2：各個(gè)力的大小的變化圖像五、模型的評(píng)價(jià)與改進(jìn) 在本模型中合理的分析了三角形鋼架的受力情況，引用矩陣求解得到力F作用下的值以及F增加1kg和減少1kg時(shí)各個(gè)力的變化值，以此為基礎(chǔ)推出了F在90kg到110kg范圍內(nèi)對(duì)應(yīng)的值。 實(shí)驗(yàn)4 一、問(wèn)題的提出給定四種物質(zhì)對(duì)應(yīng)的參數(shù)和交互作用矩陣如下：在壓強(qiáng)下，為了形成共沸混合物，溫度和組分分別是多少？請(qǐng)盡量找出所有可能的解。二、問(wèn)題分析 共沸化合物是指由兩種或兩種以上物質(zhì)組成的液體混合物，當(dāng)在某種壓力下裝備蒸餾或者局部汽化時(shí)，在氣體狀態(tài)下和液體狀態(tài)下保持相同的組分。題目中給出了四種物質(zhì)對(duì)應(yīng)的參數(shù)和交互作用矩陣，要求在壓強(qiáng)下，找出能夠形成形成共沸

20、混合物的溫度和組分。三、模型建立設(shè)該混合物由n個(gè)可能的組份組成（），組份，所占的比例為，則 （1）依據(jù)均相共沸混合物應(yīng)該滿足的條件，即共沸混合物的每個(gè)組分在氣體狀態(tài)和液體狀態(tài)下具有相同的化學(xué)勢(shì)能。建立如下模型： （2）在題目給定了四種物質(zhì)對(duì)應(yīng)的參數(shù)和交互作用矩陣以及壓強(qiáng)的條件下，模型（1）、（2）式描述了確定共享混合物組份的條件，在不考慮（1）式非負(fù)限制時(shí)，該模型為含有5個(gè)未知數(shù)和五個(gè)方程的非線性方程組，因此用數(shù)值方法求解。注意到（1）式是簡(jiǎn)單線性等式，因此我們可以從中消去一個(gè)未知數(shù)達(dá)到更好的求解效果，即 （3）將它帶入（2）式得到含有n個(gè)未知數(shù)的非線性方程組。四、模型求解用建立的模型首先編寫

21、如下M文件：function f=azeofun(XT,n,p,a,b,c,Q)x(n)=1;for i=1:n-1 x(i)=XT(i); x(n)=x(n)-x(i);endT=XT(n);p=log(p);for i=1:n d(i)=x*Q(i,1:n); dd(i)=x(i)/d(i);endfor i=1:n f(i)=x(i)*(b(i)/(T+c(i) + log(x*Q(i,1:n) + dd*Q(1:n,i) - a(i) - 1 + p);end然后用題中所給的數(shù)據(jù)，取初始值四種物質(zhì)各占1/4，溫度為50作如下計(jì)算：n=4;p=760;a=18.607 15.841 20

22、.443 19.293;b=2643.31 2755.64 4628.96 4117.07;c=239.73 219.16 252.64 227.44;Q=1.0 0.192 2.169 1.611 0.316 1.0 0.477 0.524 0.377 0.360 1.0 0.296 0.524 0.282 2.065 1.0;XT0=0.3,0.25,0.3,75;XT,Y=fsolve(azeofun,XT0,n,p,a,b,c,Q)得到XT =0.0000 0.5858 0.4142 71.9657Y =1.0e-006 *-0.0445 -0.0262 0.0055 0.1266即四

23、種物質(zhì)組成均相混合物時(shí)的比例分別為：0,58.58%，41.42%，0，溫度為71.9657。五、模型的評(píng)價(jià)與改進(jìn) 在該模型中結(jié)合了化學(xué)中共沸化合物的特性，建立了非線性模型，運(yùn)用MATLAB優(yōu)化工具包fsolve求解，并結(jié)合（1）式是簡(jiǎn)單線性等式的特征消去了非線性模型中的一個(gè)未知數(shù)，使得計(jì)算得到簡(jiǎn)化。 實(shí)驗(yàn)5 一、問(wèn)題的提出假設(shè)商品在時(shí)期的市場(chǎng)價(jià)格為，當(dāng)需求函數(shù)為，而生產(chǎn)方的期望價(jià)格為，供應(yīng)函數(shù)為，當(dāng)供銷平衡時(shí)。若期望價(jià)格與市場(chǎng)價(jià)格不符，商品市場(chǎng)不均衡，生產(chǎn)方t+1時(shí)期的期望價(jià)格將會(huì)調(diào)整，方式為。以帶入，得到關(guān)于的遞推方程。設(shè)，以c為可變參數(shù)，討論期望價(jià)格的變化規(guī)律，是否有混沌現(xiàn)象出現(xiàn)，并找出

24、前幾個(gè)分岔點(diǎn)，觀察分岔點(diǎn)的極限趨勢(shì)是否符合Feigenbaum常數(shù)揭示的規(guī)律。二、問(wèn)題分析為了保證商品能滿足供銷平衡，在期望價(jià)格與市場(chǎng)價(jià)格不符時(shí)，生產(chǎn)方將會(huì)在第t+1時(shí)期調(diào)整期望價(jià)格，并得到關(guān)于市場(chǎng)價(jià)格的遞推方程。題中已知供應(yīng)函數(shù)與需求函數(shù)的參數(shù)d，以c為可變參數(shù)，討論期望價(jià)格的變化規(guī)律，是否有混沌現(xiàn)象出現(xiàn)，并找出前幾個(gè)分岔點(diǎn)，觀察分岔點(diǎn)的極限趨勢(shì)是否符合Feigenbaum常數(shù)揭示的規(guī)律。三、模型建立依據(jù)題目所給出的信息，記第t時(shí)期的期望價(jià)格為，得到t+1時(shí)期的期望價(jià)格的非線性差分方程：其中。不妨取初始值，以c為可變參數(shù)，討論期望價(jià)格的變化規(guī)律。觀察期望價(jià)格是否有出現(xiàn)混沌現(xiàn)象出現(xiàn)。找出其前幾

25、個(gè)分岔點(diǎn)，討論分岔點(diǎn)的極限趨勢(shì)是否符合Feigenbaum常數(shù)揭示的規(guī)律。四、模型求解依據(jù)建立的模型編寫如下MATLAB程序:c=1:8;q=100;n=50;for j=1:8; C=c(j); for t=1:n q(t+1)=0.7*q(t)- tand(4.8*q(t)+0.3*C; end Q(:,j)=q;endk=(0:50);k,Qsubplot(4,2,1),plot(k,Q(:,1),subplot(4,2,2),plot(k,Q(:,2),subplot(4,2,3),plot(k,Q(:,3),subplot(4,2,4),plot(k,Q(:,4),subplot(4

26、,2,5),plot(k,Q(:,5),subplot(4,2,6),plot(k,Q(:,6),subplot(4,2,7),plot(k,Q(:,7),subplot(4,2,8),plot(k,Q(:,8)得到結(jié)果如表3和圖3所示：k（c=1）（c=2）（c=3）（c=4）（c=5）（c=6）（c=7）（c=8）0 100.0000 100.0000 100.0000 100.0000 100.0000 100.0000 100.0000 100.0000 1 72.0321 72.3321 72.6321 72.9321 73.2321 73.5321 73.8321 74.1321

27、2 50.9763 51.4597 51.9435 52.4274 52.9116 53.3960 53.8806 54.3653 3 33.8692 34.2652 34.6100 34.8840 35.0563 35.0748 34.8449 34.1751 4 24.3224 24.8635 25.3740 25.8415 26.2470 26.5583 26.7176 26.6085 21 21.4611 19.5144 26.6675 19.9913 23.3474 22.5641 21.6621 21.6251 22 19.6497 29.8544 20.8470 24.7757

28、20.3099 20.6172 21.2808 21.6087 23 27.2971 22.2434 21.1265 20.3527 23.3255 22.5725 21.6422 21.6215 24 20.5573 19.4892 20.6448 22.8498 20.3076 20.6160 21.2956 21.6115 25 21.2443 30.3707 21.5979 19.9910 23.3352 22.5759 21.6247 21.6193 46 23.2879 32.4736 23.2521 24.7772 20.3083 20.6152 21.3922 21.6158

29、47 19.1040 23.7794 19.7009 20.3532 23.3322 22.5782 21.5182 21.6159 48 47.3819 19.4768 27.2173 22.8477 20.3083 20.6152 21.3972 21.6159 49 32.3786 30.6367 21.1170 19.9911 23.3322 22.5782 21.5130 21.6159 50 23.4225 22.6937 20.6585 24.7772 20.3083 20.6152 21.4018 21.6159 表3：不同參數(shù)c下期望價(jià)格的變化規(guī)律圖3： 不同參數(shù)c下期望價(jià)格

30、的變化規(guī)律從圖3和表3可以看出當(dāng)c=1,2，3時(shí)，沒(méi)有任何收斂的子序列，期望價(jià)格沒(méi)有任何變化規(guī)律；當(dāng)c=4時(shí)有三個(gè)收斂的子序列，分別趨向于20、22和24；當(dāng)c=5時(shí)有兩個(gè)收斂的子序列，分別趨向于20和22；當(dāng)c=6和7時(shí)有趨向于21，可以看做期望價(jià)格期期收斂。為了分析非線性差分方程出現(xiàn)的混沌現(xiàn)象，用MATLAB編程畫出非線性迭代序列隨著參數(shù)變化的收斂，分岔和混沌現(xiàn)象圖。首先編寫如下程序：function chaos(iter_fun,x0,r,n) % 該函數(shù)沒(méi)有返回值；iter_fun是迭代函數(shù)；kr=0; for rr=r(1):r(3):r(2) % 輸入中r(1),r(2) 是參數(shù)變

31、化的范圍，r（3）是步長(zhǎng) kr=kr+1; y(kr,1)=feval(iter_fun,x0,rr); for i=2:n(2) y(kr,i)=feval(iter_fun,y(kr,i-1),rr); endendplot(r(1):r(3):r(2),y(:,n(1)+1:n(2),k.); 針對(duì)本題的迭代函數(shù)寫如下M文件：function y=iter01(q,c)y=0.7*q-tand(4.8*q)+0.3*c;輸入命令chaos(iter01,100,2.5,7.1,0.01,50,300)可以得到圖4所示的分岔和混沌圖（橫坐標(biāo)寶石變化的參數(shù)，縱坐標(biāo)表示迭代的極限）。圖4：期望

32、價(jià)格解的收斂分岔混沌現(xiàn)象從圖中可以看出隨著c的減小序列收斂性出現(xiàn)分叉：一分二，二分四，直至出現(xiàn)混沌。因此，分岔點(diǎn)的極限趨勢(shì)是否符合Feigenbaum常數(shù)揭示的規(guī)律。五、模型的評(píng)價(jià)與改進(jìn) 該模型中依據(jù)期望價(jià)格的變化規(guī)律合理的應(yīng)用了非線性差分方程，對(duì)c先抽取具有代表性的8個(gè)值用MATLAB求解期望價(jià)格的遞增值觀察變化規(guī)律，驗(yàn)證其中是有混沌現(xiàn)象出現(xiàn)，而后對(duì)出現(xiàn)的混沌現(xiàn)象作圖觀察其分岔點(diǎn)出現(xiàn)位置，以及分叉點(diǎn)的極限。 實(shí)驗(yàn)6 一、問(wèn)題的提出取不同的初值計(jì)算非線性規(guī)劃，盡可能求出所有局部極小點(diǎn)，進(jìn)而找出全局極小點(diǎn)，并對(duì)不同算法的結(jié)果進(jìn)行分析比較。二、問(wèn)題分析題目要求對(duì)非線性規(guī)劃在求出其所有局部極小點(diǎn)的基

33、礎(chǔ)上，找出全局極小點(diǎn)。并比較不同初值以及不同算法的計(jì)算結(jié)果。三、模型建立對(duì)于類似題目中所給出的無(wú)約束優(yōu)化問(wèn)題，通常先尋找問(wèn)題的局部最優(yōu)解，然后將局部最優(yōu)解進(jìn)行對(duì)比得到全局最優(yōu)解。在此采用擬牛頓法和單純型搜索法（直接法），選取不同的初值求解。其搜索最優(yōu)解的基本思想是用迭代法搜索最優(yōu)解，在迭代的第k步，即對(duì)n維空間中的一點(diǎn)，確定一個(gè)搜索方向和步長(zhǎng)，沿此方向、按此步長(zhǎng)走一步到達(dá)下一點(diǎn)時(shí)，函數(shù)值下降。基本步驟為：1、選初始解。2、對(duì)于第次迭代解，確定搜索方向，并在此方向確定搜索步長(zhǎng)，令，使。3、若符合給定的迭代終止原則，迭代停止，最優(yōu)解；否則，轉(zhuǎn)步驟2。四、模型求解依據(jù)建立的模型，首先建立M文件計(jì)算函

34、數(shù)值：function x=shiyan6(x,c1,c2,a1,a2);f1=-(x-a1)*(x-a1)+c1;f2=-(x-a2)*(x-a2)+c2;z=1/f1+1/f2;取初始值1,1，輸入程序：x0=1,1; % 初值以下更改初值的過(guò)程中更改此處for i=1:2 a1=4,4; a2=2.5,3.8; z1,f1,e1,out1=fminunc(shiyan6,x0(i),0.7,0.73,a1(i),a2(i) z2,f2,e2,out2=fminsearch(shiyan6,x0(i),0.7,0.73,a1(i),a2(i)end分別計(jì)算得到該初始值下擬牛頓法和單純型搜索

35、法的計(jì)算結(jié)果分別為-4.3585e+008, -6.3383e+028和-4.3585e+008, -6.3383e+028。兩種算法最終得到的計(jì)算結(jié)果是一致的。但是從輸出結(jié)果的out一項(xiàng)中分析可以看出單純型搜索法的迭代次數(shù)（100）遠(yuǎn)大于擬牛頓法（2），相比之下當(dāng)計(jì)算量大時(shí)擬牛頓法更合適。改變初值，多次運(yùn)行程序可以發(fā)現(xiàn)，在初值的改變對(duì)最終的計(jì)算結(jié)果沒(méi)有影響，初值的選取僅導(dǎo)致迭代次數(shù)不一樣。從以上結(jié)果來(lái)看，運(yùn)算結(jié)果不受選取的計(jì)算方法和初始值的影響，但是初值以及計(jì)算方法的選取將影響迭代次數(shù)。五、模型的評(píng)價(jià)與改進(jìn) 模型中應(yīng)題目要求，采用了擬牛頓法和單純型搜索法以及在不同算法下選取了不同的初始值，并

36、對(duì)得到的結(jié)果經(jīng)行了對(duì)比分析。 實(shí)驗(yàn)7 一、問(wèn)題的提出已知某分子由25個(gè)原子組成，并且已經(jīng)通過(guò)實(shí)驗(yàn)測(cè)量得到了其中某些原子對(duì)之間的距離如表4所示（假設(shè)在平面結(jié)構(gòu)上討論），請(qǐng)確定每個(gè)原子的位置關(guān)系。原子對(duì)距離原子對(duì)距離原子對(duì)距離原子對(duì)距離（4,1）0.9607(5,4)0.4758(18,8)0.8363(15,13)0.5725（12,1）0.4399(12,4)1.3402(13,9)0.3208(19,13)0.7660(13,1)0.8143(24,4)0.7006(15,9)0.1574(15,14)0.4394(17,1)1.3765(8,6)0.4945(22,9)1.2736(16,

37、14)1.0952(21,1)1.2722(13,6)1.0559(11,10)0.5781(20,16)1.0422(5,2)0.5294(19,6)0.6810(13,10)0.9254(23,16)1.8255(16,2)0.6144(25,6)0.3587(19,10)0.6401(18,17)1.4325(17,2)0.3766(8,7)0.3351(20,10)0.2467(19,17)1.0851(25,2)0.6893(14,7)0.2878(22,10)0.4727(20,19)0.4995(5,3)0.9488(16,7)1.1346(18,11)1.3840(23,19)

38、1.2277(20,3)0.8000(20,7)0.3870(25,11)0.4366(24,19)1.1271(21,3)1.1090(21,7)0.7511(15,12)1.0307(23,21)0.7060(24,3)1.1432(14,8)0.4439(17,12)1.3904(23,22)0.8025表4：分子中某些原子對(duì)之間的距離二、問(wèn)題分析 依據(jù)分子中已知的原子對(duì)之間的距離，在假設(shè)所有原子均在一個(gè)平面上的基礎(chǔ)上，建立相應(yīng)的模型確定每個(gè)原子的位置關(guān)系。三、模型建立記為第個(gè)原子的橫坐標(biāo)，為第個(gè)原子的縱坐標(biāo)，為第個(gè)原子與第個(gè)原子之間的距離，用歐氏距離紀(jì)錄任意兩個(gè)原子之間的距離有：現(xiàn)已知

39、某些原子對(duì)之間的距離，求解每個(gè)原子的位置關(guān)系，在此不妨假設(shè)第一個(gè)原子的坐標(biāo)為（0,0）。建立無(wú)約束最優(yōu)化模型求解每個(gè)原子的坐標(biāo)?？梢詺w結(jié)到無(wú)約束規(guī)劃模型，令，令，調(diào)用lsqnonlin命令來(lái)做。 四、模型求解依據(jù)上述建立的模型，首先編寫對(duì)應(yīng)條件第和第個(gè)原子間距離的M文件如下：function f=distance(x,d) f(1)=(x(1,3)-0)2+(x(2,3)-0)2-d(1)2 f(2)=(x(1,11)2+(x(2,11)2-d(2)2; f(3)=(x(1,12)2+(x(2,12)2-d(3)2; f(4)=(x(1,16)2+(x(2,16)2-d(4)2;f(5)=(x

40、(1,20)2+(x(2,20)2-d(5)2; f(6)=(x(1,4)-x(1,1)2+(x(2,4)-x(2,1)2-d(6)2; f(7)=(x(1,15)-x(1,1)2+(x(2,15)-x(2,1)2-d(7)2; f(8)=(x(1,16)-x(1,1)2+(x(2,16)-x(2,1)2-d(8)2; f(9)=(x(1,24)-x(1,1)2+(x(2,24)-x(2,1)2-d(9)2; f(10)=(x(1,4)-x(1,2)2+(x(2,4)-x(2,2)2-d(10)2; f(11)=(x(1,19)-x(1,2)2+(x(2,19)-x(2,2)2-d(11)2;

41、 f(12)=(x(1,20)-x(1,2)2+(x(2,20)-x(2,2)2-d(12)2; f(13)=(x(1,23)-x(1,2)2+(x(2,23)-x(2,2)2-d(13)2; f(14)=(x(1,4)-x(1,3)2+(x(2,4)-x(2,3)2-d(14)2; f(15)=(x(1,11)-x(1,3)2+(x(2,11)-x(2,3)2-d(15)2; f(16)=(x(1,23)-x(1,3)2+(x(2,23)-x(2,3)2-d(16)2; f(17)=(x(1,7)-x(1,5)2+(x(2,7)-x(2,5)2-d(17)2; f(18)=(x(1,12)-

42、x(1,5)2+(x(2,12)-x(2,5)2-d(18)2; f(19)=(x(1,18)-x(1,5)2+(x(2,18)-x(2,5)2-d(19)2; f(20)=(x(1,24)-x(1,5)2+(x(2,24)-x(2,5)2-d(20)2; f(21)=(x(1,7)-x(1,6)2+(x(2,7)-x(2,6)2-d(21)2; f(22)=(x(1,13)-x(1,6)2+(x(2,13)-x(2,6)2-d(22)2; f(23)=(x(1,15)-x(1,6)2+(x(2,15)-x(2,6)2-d(23)2; f(24)=(x(1,19)-x(1,6)2+(x(2,1

43、9)-x(2,6)2-d(24)2; f(25)=(x(1,20)-x(1,6)2+(x(2,20)-x(2,6)2-d(25)2; f(26)=(x(1,13)-x(1,7)2+(x(2,13)-x(2,7)2-d(26)2;f(27)=(x(1,17)-x(1,7)2+(x(2,17)-x(2,7)2-d(27)2;f(28)=(x(1,12)-x(1,8)2+(x(2,12)-x(2,8)2-d(28)2; f(29)=(x(1,14)-x(1,8)2+(x(2,14)-x(2,8)2-d(29)2; f(30)=(x(1,21)-x(1,8)2+(x(2,21)-x(2,8)2-d(3

44、0)2; f(31)=(x(1,10)-x(1,9)2+(x(2,10)-x(2,9)2-d(31)2;f(32)=(x(1,12)-x(1,9)2+(x(2,12)-x(2,9)2-d(32)2; f(33)=(x(1,18)-x(1,9)2+(x(2,18)-x(2,9)2-d(33)2; f(34)=(x(1,19)-x(1,9)2+(x(2,19)-x(2,9)2-d(34)2; f(35)=(x(1,21)-x(1,9)2+(x(2,21)-x(2,9)2-d(35)2; f(36)=(x(1,17)-x(1,10)2+(x(2,17)-x(2,10)2-d(36)2; f(37)=

45、(x(1,24)-x(1,10)2+(x(2,24)-x(2,10)2-d(37)2; f(38)=(x(1,14)-x(1,11)2+(x(2,14)-x(2,11)2-d(38)2; f(39)=(x(1,16)-x(1,11)2+(x(2,16)-x(2,11)2-d(39)2; f(40)=(x(1,14)-x(1,12)2+(x(2,14)-x(2,12)2-d(40)2; f(41)=(x(1,18)-x(1,12)2+(x(2,18)-x(2,12)2-d(41)2; f(42)=(x(1,14)-x(1,13)2+(x(2,14)-x(2,13)2-d(42)2; f(43)=

46、(x(1,15)-x(1,13)2+(x(2,15)-x(2,13)2-d(43)2; f(44)=(x(1,19)-x(1,15)2+(x(2,19)-x(2,15)2-d(44)2; f(45)=(x(1,22)-x(1,15)2+(x(2,22)-x(2,15)2-d(45)2; f(46)=(x(1,17)-x(1,16)2+(x(2,17)-x(2,16)2-d(46)2; f(47)=(x(1,18)-x(1,16)2+(x(2,18)-x(2,16)2-d(47)2; f(48)=(x(1,19)-x(1,18)2+(x(2,19)-x(2,18)2-d(48)2; f(49)=

47、(x(1,22)-x(1,18)2+(x(2,22)-x(2,18)2-d(49)2; f(50)=(x(1,23)-x(1,18)2+(x(2,23)-x(2,18)2-d(50)2; f(51)=(x(1,22)-x(1,20)2+(x(2,22)-x(2,20)2-d(51)2; f(52)=(x(1,22)-x(1,21)2+(x(2,22)-x(2,21)2-d(52)2; 其次，調(diào)用lsqnonlin命令，編寫如下程序求解任意兩個(gè)原子的位置關(guān)系：x0=zeros(1,24);ones(1,24);d=0.9607,0.4399,0.8143,1.3765,1.2722,0.5294

48、,0.6144,0.3766,0.6893,. 0.9488,0.8000,1.1090,1.1432,0.4758,1.3402,0.7006,0.4945,1.0559,. 0.6810,0.3587,0.3351,0.2878,1.1346,0.3870,0.7511,0.4439,0.8363,. 0.3208,0.1574,1.2736,0.5781,0.9254,0.6401,0.2467,0.4727,1.3840,. 0.4366,1.0307,1.3904,0.5725,0.7660,0.4394,1.0952,1.0422,1.8255,. 1.4325,1.0851,0.

49、4995,1.2277,1.1271,0.7060,0.8052;%設(shè)定初值x,norms=lsqnonlin(distance,x0,d);p=x; %p第一行為第二個(gè)原子的橫縱坐標(biāo) a=0,x(1,:); b=0,x(2,:); plot(a,b,*); %畫散點(diǎn)圖表示出原子的位置執(zhí)行程序得到任意兩個(gè)原子的位置關(guān)系如表5和圖5所示序號(hào)橫坐標(biāo)縱坐標(biāo)序號(hào)橫坐標(biāo)縱坐標(biāo)序號(hào)橫坐標(biāo)縱坐標(biāo)20.40951.1753100.15201.4901181.07910.745930.74761.639311-0.26341.0793190.75791.38124-0.56560.8388120.51860.1

50、088200.16571.14965-0.10521.1593130.51660.633621-0.30121.254260.11411.5906140.15190.6042220.43561.939970.09910.6862150.08671.028423-0.36071.960980.28171.0535161.17170.995624-0.37461.504290.16010.705017-0.22841.316525-0.20531.5010表5：原子的坐標(biāo)圖4：原子的位置示意圖為了評(píng)價(jià)該算法的精確性，在執(zhí)行除以上程序結(jié)果的基礎(chǔ)，繼續(xù)編寫如下程序：d=0.9607,0.4399,0.

51、8143,1.3765,1.2722,0.5294,0.6144,0.3766,0.6893,. 0.9488,0.8000,1.1090,1.1432,0.4758,1.3402,0.7006,0.4945,1.0559,. 0.6810,0.3587,0.3351,0.2878,1.1346,0.3870,0.7511,0.4439,0.8363,. 0.3208,0.1574,1.2736,0.5781,0.9254,0.6401,0.2467,0.4727,1.3840,. 0.4366,1.0307,1.3904,0.5725,0.7660,0.4394,1.0952,1.0422,

52、1.8255,.1.4325,1.0851,0.4995,1.2277,1.1271,0.7060,0.8052;f=distance(x,d)x=p;a=0;for i=1:52a=a+f(i);enda得到結(jié)果該模型的精度結(jié)果為a=0.8072，由此可以看出以上建立的無(wú)約束最優(yōu)化模型求解得到的任意兩個(gè)原子位置關(guān)系的結(jié)果是可靠的。五、模型的評(píng)價(jià)與改進(jìn)在本模型中運(yùn)用了無(wú)約束最優(yōu)化模型的求解，調(diào)用了lsqnonlin命令。并計(jì)算了結(jié)果的精確度，從精確度的角度驗(yàn)證結(jié)果的可靠性。該模型的還可以用調(diào)用基本命令fminunc來(lái)做，類似調(diào)用lsqnonlin命令首先編寫對(duì)應(yīng)條件第和第個(gè)原子間距離的M文件，

53、然后寫程序：clear all;d=0.9607,0.4399,0.8143,1.3765,1.2722,0.5294,0.6144,0.3766,0.6893,0.9488,0.8000,1.1090,1.1432,0.4758,1.3402,0.7006,0.4945,1.0559,0.6810,0.3587,0.3351,0.2878,1.1346,0.3870,0.7511,0.4439,0.8363,0.3208,0.1574,1.2736,0.5781,0.9254,0.6401,0.2467,0.4727,1.3840,0.4366,1.0307,1.3904,0.5725,0.

54、7660,0.4394,1.0952,1.0422,1.8255,1.4325,1.0851,0.4995,1.2277,1.1271,0.7060,0.8052;x0=zeros(1,26),ones(1,24);opt = optimset(tolx,1e-16,tolf,1e-16);x,z,exit1,out1 = fminunc(distance1,x0,opt,d);p=0,0;x(1:24),x(25:48);a=0,x(1:24);b=0,x(25:48);plot(a,b,*);可以得到原子的分布示意圖如圖5所示。也可以計(jì)算得到得到結(jié)果該模型的精度結(jié)果為f =0.0375圖4

55、：原子的位置示意圖 實(shí)驗(yàn)8 一、問(wèn)題的提出 給藥方案設(shè)計(jì)需要依據(jù)藥物吸收與排除過(guò)程的原理。藥物進(jìn)入機(jī)體后隨血液輸送到全身，不斷地被吸收、分布、代謝，最終排出體外。藥物在血液中的濃度，即單位體積血液中的藥物含量，稱血藥濃度。在最簡(jiǎn)單的一室模型中，將整個(gè)機(jī)體看作一個(gè)房室，稱中心室，室內(nèi)的血藥濃度是均勻的。這里我們用一室模型，討論在口服給藥方式下血藥濃度的變化規(guī)律，及根據(jù)實(shí)驗(yàn)數(shù)據(jù)擬合參數(shù)的方法?？诜o藥方式相當(dāng)于先有一個(gè)將藥從腸胃吸收入血液的過(guò)程，這個(gè)過(guò)程可簡(jiǎn)吸收室c1（t）中心室 c（t）k1k化為在藥物進(jìn)入中心室之前有一個(gè)吸收室如右圖所示，記中心室和吸收室的容積分別為，而t時(shí)刻的血藥濃度分別為，

56、；中心室的排除速率為k，吸收速率為（這里k和分別是中心室和吸收室血藥濃度變化率與濃度本身的比例系數(shù)）。設(shè)t=0時(shí)刻口服劑量為d的藥物，容易寫出吸收室的血藥濃度的微分方程為中心室血藥濃度的變化率由兩部分組成：與c成正比的排除（比例系數(shù)）；與成正比的吸收（比例系數(shù)）。再考慮到中心室和吸收室的容積分別為，得到的微分方程為由以上兩個(gè)微分方程不難解出中心室血藥濃度設(shè)t=0時(shí)刻口服一定劑量的藥物，實(shí)驗(yàn)數(shù)據(jù)如表6所示，請(qǐng)由此確定k，b。t0.0830.1670.250.50.7511.5c(t)10.921.127.336.435.538.434.8t2.2534681012c(t)24.223.615.78.28