南郵應(yīng)用物理計(jì)算物理實(shí)踐

....南京郵電大學(xué)實(shí)驗(yàn)報(bào)告課程名稱：計(jì)算物理實(shí)踐專業(yè)：應(yīng)用物理學(xué)學(xué)號(hào)：姓名：完成日期：年月

目錄一、簡(jiǎn)單物理實(shí)驗(yàn)的模擬及實(shí)驗(yàn)數(shù)據(jù)處理 11.1問(wèn)題描述： 11.2單擺運(yùn)動(dòng)原理 11.3模型的建立 11.4流程圖 21.5Matlab程序設(shè)計(jì)仿真 21.6Matlab程序 31.7單擺演示截圖 3二、方程組的數(shù)值解法 42.1問(wèn)題描述： 42.2原理分析 42.2.1二分法理論 42.2.2分析求解 42.3Matlab程序 52.4Matlab程序運(yùn)行結(jié)果： 5三、靜電場(chǎng)問(wèn)題的計(jì)算 73.1問(wèn)題描述： 73.2原理分析 73.2．1簡(jiǎn)單迭代法 73.2.2有限差分： 83.2.3解題過(guò)程 93.3Matlab程序設(shè)計(jì)仿真 93.4Matlab仿真結(jié)果 10四、熱傳導(dǎo)方程和波動(dòng)方程的差分解法 114.1問(wèn)題描述 114.2原理分析 114.3具體步驟 134.4MATLAB程序設(shè)計(jì)仿真 134.5MATLAB程序運(yùn)行結(jié)果 13結(jié)束語(yǔ) 15參考文獻(xiàn) 16附錄1： 17附錄2： 18附錄3： 19附錄4： 21..一、簡(jiǎn)單物理實(shí)驗(yàn)的模擬及實(shí)驗(yàn)數(shù)據(jù)處理1.1問(wèn)題描述：編寫單擺運(yùn)動(dòng)演示程序。在不考慮空氣阻力和很小的假設(shè)下.單位質(zhì)量小球做理想簡(jiǎn)諧運(yùn)動(dòng).此時(shí)。取g=9.8.L=1,0=pi/4.1.2單擺運(yùn)動(dòng)原理設(shè)在某一時(shí)刻.單擺的擺線偏離垂直線的角位移為.將重力mg分解為徑向力F和切向力T.則T的大小為mgsin,切向加速度為a=Ld2dt2.根據(jù)牛頓第二定律得方程ma=mLd從而單擺運(yùn)動(dòng)的微分方程為d2dt由于sin=—33!+55當(dāng)很小時(shí).sin≈所以單擺的微分方程可表示為d2dt2=上式表明.當(dāng)很小時(shí).單擺的角加速度與角位移成正比.但方向相反.且方程的解可表示為=Acos（ωt+φ）②1.3模型的建立建立物理模型.假設(shè)單擺運(yùn)動(dòng)過(guò)程中的小球中心點(diǎn)的坐標(biāo)為x,y。根據(jù)幾何關(guān)系x,y滿足如下關(guān)系：x在不考慮空氣阻力和很小的假設(shè)下=0costgL與②所以.單擺的運(yùn)動(dòng)方程為x=L將0=pi/4,g=9.8,L=1帶入得.x=1.4流程圖開始開始輸入0輸入0t=0,dt=0.005t=0,dt=0.005t=1?Yt=1?Y=-L*cosAngle*cosY=-L*cosAngle*X=L*sin(Angle*結(jié)束t=t+dt結(jié)束t=t+dt圖1.1程序流程圖1.5Matlab程序設(shè)計(jì)仿真通過(guò)set函數(shù)將變現(xiàn)和小球的圖像句柄.加入X、Y的動(dòng)態(tài)參量.它們得軌跡變化即為上面分析的軌跡方程。用line函數(shù)畫出初始的位置.并將句柄分別給sphere、lp.然后每隔dt時(shí)間刷新一次。嘗試對(duì)dt的設(shè)置發(fā)現(xiàn).值在0.0005的效果比較好。1.6Matlab程序程序見(jiàn)附錄1。1.7單擺演示截圖圖1.2單擺演示圖二、方程組的數(shù)值解法2.1問(wèn)題描述：二分法求解方程x3+4x2－10＝0在區(qū)間［1.2］內(nèi)的根.精度自設(shè)。2.2原理分析2.2.1二分法理論f(x)∈Ca,b,單調(diào).f(a)f(b)<0f(x)=0在（a,b）有唯一根。設(shè)f(x)在[a,b]上連續(xù).f(x)=0在[a,b]上存在唯一解.a第一步.計(jì)算f(a0)f(x0).若f(a0)f(x0)<0,則 x*∈(a0,b0)第k步.計(jì)算f(ak-1)f(xk-1)<0,則x*∈[ak-1,xk-1],記ak?k,x*∈[ak,bk]且xk=ak+bk2.所以數(shù)列｛即xk→x*,從而當(dāng)k充分大.x*≈xk且可由2.2.2分析求解令f(x)=-x3+4x2－10,f(x)在［1.2］上連續(xù).且f(1)f(2)<0,則f(x)=0在［1.2］上有唯一解.記a=1,b=2,x=1.5.然后計(jì)算f(a)f(b),若f(a)f(b)<0.則x∈(1,1.5),此時(shí)記a=1,b=1.5,否則x∈(1.5,2),記a=1.5,b=2.對(duì)兩種情形均有x∈[a,b],記x=a+b2,按照同樣的方法依次向下計(jì)算.直到求出的相鄰兩個(gè)x的值之差絕對(duì)誤差小于0.00005.2.3Matlab程序程序見(jiàn)附錄2.2.4Matlab程序運(yùn)行結(jié)果：x=1.2500x=1.3750x=1.3125x=1.3438x=1.3594x=1.3672x=1.3633x=1.3652x=1.3643x=1.3647x=1.3650x=1.3651x=1.3652x=1.3652x=1.3652x=1.3652x=1.3652x=1.36523,f(x)=0.00001三、靜電場(chǎng)問(wèn)題的計(jì)算3.1問(wèn)題描述：設(shè)兩個(gè)同軸矩形金屬槽如圖3-1所示.外金屬槽電位為0.內(nèi)金屬槽電位為100V.求內(nèi)電位分布.并繪出電位分布圖。圖3-13.2原理分析3.2．1簡(jiǎn)單迭代法對(duì)某一網(wǎng)格點(diǎn)設(shè)一初值.這個(gè)初值完全可以任意給定.稱為初值電位。雖然.問(wèn)題的最終結(jié)果與初值無(wú)關(guān).但若初值選擇得當(dāng).則計(jì)算步驟會(huì)得到簡(jiǎn)化（當(dāng)利用計(jì)算機(jī)來(lái)實(shí)現(xiàn)迭代計(jì)算時(shí).為了簡(jiǎn)化程序.初值點(diǎn)為一般可取值為零）。初值電位給定后.按一個(gè)固定順序（點(diǎn)的順序是從左到右.從下到上）依次計(jì)算每點(diǎn)的點(diǎn)位.即利用φi,j=14(φi,j+1+φi,j-1簡(jiǎn)單迭代法的特點(diǎn)是用之前一次迭代得到的網(wǎng)絡(luò)點(diǎn)電位作為下一次迭代的初值。如在（i,j）點(diǎn)在n+1次迭代時(shí)計(jì)算公式為：φ3.2.2有限差分：二維拉普拉斯方程（1）有限差分法的網(wǎng)格劃分.通常采用完全有規(guī)律的分布方式.這樣可使每個(gè)離散點(diǎn)上得到相同形式的差分方程.有效的提高解題速度.經(jīng)常采用的是正方形網(wǎng)格劃分。圖3-2迭代法網(wǎng)格劃分設(shè)網(wǎng)格節(jié)點(diǎn)（i,j）的電位為.其上下左右四個(gè)節(jié)點(diǎn)的電位分別為在h充分小的情況下.可以為基點(diǎn)進(jìn)行泰勒級(jí)數(shù)展開：把以上四式相加.在相加的過(guò)程中.h的所有奇次方項(xiàng)都抵消了。得到的結(jié)果的精度為h的二次項(xiàng)。（2）由于場(chǎng)中任意點(diǎn)都滿足泊松方程：式中為場(chǎng)源.則式（2）可變?yōu)椋海?）對(duì)于無(wú)源場(chǎng)..則二維拉普拉斯方程的有限差分形式為：（4）上式表示任一點(diǎn)的電位等于圍繞它的四個(gè)等間距點(diǎn)的電位的平均值.距離h越小則結(jié)果越精確.用式（4）可以近似的求解二維拉普拉斯方程。邊界條件：|內(nèi)槽=100；|外槽=03.2.3解題過(guò)程解：在直角坐標(biāo)系中.金屬槽中的電位函數(shù)φ滿足拉普拉斯方程：?其邊界條件滿足混合型編值問(wèn)題的邊界條件：φ?φ取步長(zhǎng)h=1,x、y方向的網(wǎng)格數(shù)為m=16,n=10,共有16*10=160個(gè)網(wǎng)孔17*11=187個(gè)節(jié)點(diǎn).其中槽內(nèi)節(jié)點(diǎn)（電位待求點(diǎn)）有15*9=135個(gè).界節(jié)點(diǎn)（電位已知點(diǎn)）有187—135=52個(gè)。設(shè)迭代精度為103.3Matlab程序設(shè)計(jì)仿真源程序見(jiàn)附錄三3.4Matlab仿真結(jié)果圖3-3運(yùn)行結(jié)果圖四、熱傳導(dǎo)方程和波動(dòng)方程的差分解法4.1問(wèn)題描述求熱傳導(dǎo)方程混合問(wèn)題：的數(shù)值解.其中N、ｈ、k值等參數(shù)自?。▽⒂?jì)算結(jié)果圖形化）。4.2原理分析二維熱傳導(dǎo)方程的初、邊值混合問(wèn)題與一維的相似.在確定差分方程格式并給出定解條件后.按時(shí)間序號(hào)分層計(jì)算.只是每一層是由二維點(diǎn)陣組成.通常稱為網(wǎng)格。各向同性介質(zhì)中無(wú)熱源的二維熱傳導(dǎo)方程為：.(,,,).①初始條件是：②設(shè)時(shí)間步長(zhǎng)為τ.空間步長(zhǎng)為h.二維平面xoy分為M×N的網(wǎng)格.并使..（）.則有（）.對(duì)節(jié)點(diǎn).在時(shí)刻（即時(shí)刻）有：③將差分格式③代入偏微分方程①中.可得：④（i=1,2,…,N-1,j=1,2,就所設(shè)定的具體問(wèn)題來(lái)討論邊界條件:在邊界的和區(qū)域以及整個(gè).邊界均為絕熱壁；而在邊界的區(qū)域?yàn)榕c恒溫?zé)嵩聪噙B的口。和兩邊界溫度始終為0.實(shí)際上也是與恒溫源相連的。也就是說(shuō).對(duì)于絕熱壁應(yīng)滿足：,（上述邊界條件的差分近似式為：即：（）對(duì)于與恒溫源相連的邊界.在熱傳導(dǎo)過(guò)程中始終有恒定的熱流.常可取歸一化值.例如高溫?zé)嵩纯扇　?”.而低溫?zé)嵩纯扇　?”。對(duì)于這種情況邊界條件還有：綜合上述初值、邊值混合問(wèn)題.并設(shè)初始時(shí)刻各點(diǎn)溫度均為零.則上述差分格式可歸納為：u(i,k+1)=a*u(i+1,k)+(1-2*a)*u(i,k)+a*u(i-1,k)u(0,k)=0u(1,k)=1可以證明.對(duì)于二維熱傳導(dǎo)方程.若滿足.則差分格式式④或式⑥就是穩(wěn)定的差分格式.一般的講.對(duì)于n維拋物線型微分方程差分格式穩(wěn)定的充分條件是：。4.3具體步驟1.給定、、和.取h=0.1.2.計(jì)算初值和邊值：；；；；；用差分格式計(jì)算；4.4MATLAB程序設(shè)計(jì)仿真源程序見(jiàn)附錄四MATLAB程序運(yùn)行結(jié)果圖4-1運(yùn)行結(jié)果圖結(jié)束語(yǔ)計(jì)算物理學(xué)是利用電子計(jì)算機(jī)進(jìn)行數(shù)據(jù)采集、數(shù)值計(jì)算和數(shù)字仿真來(lái)發(fā)現(xiàn)和研究物理現(xiàn)象與物理規(guī)律的一門現(xiàn)代交叉學(xué)科。經(jīng)過(guò)本次計(jì)算物理學(xué)實(shí)驗(yàn)周的學(xué)習(xí).我學(xué)會(huì)了怎樣利用計(jì)算機(jī)來(lái)進(jìn)行數(shù)據(jù)分析和用Matlab來(lái)進(jìn)行數(shù)學(xué)以及物理問(wèn)題的仿真而得到問(wèn)題的結(jié)果。我認(rèn)識(shí)到自己對(duì)于以前學(xué)習(xí)過(guò)的一些課程掌握得還不夠.Matlab編程語(yǔ)言的運(yùn)用也不夠熟練。通過(guò)這次實(shí)驗(yàn)也很好的鞏固了以前學(xué)習(xí)的一些知識(shí)點(diǎn).將理論運(yùn)用于實(shí)踐。這次實(shí)驗(yàn)讓我認(rèn)識(shí)到數(shù)理方程的實(shí)用性.掌握了利用差分代替微分來(lái)求解波動(dòng)方程、熱傳導(dǎo)方程、拉普拉斯方程等的基本原理和方法。本次實(shí)踐涉及到的二維拉普拉斯方程以及二維熱傳導(dǎo)方程的解題方法.都是先將連續(xù)的方程以及邊界條件離散化.再用計(jì)算機(jī)進(jìn)行計(jì)算.因?yàn)橛?jì)算機(jī)智能對(duì)離散的數(shù)值進(jìn)行計(jì)算。對(duì)于非線性方程的求解往往是采用迭代的方法求解.本次實(shí)踐主要涉及了Newton迭代法的重要思想.也是將連續(xù)的方程離散化后再進(jìn)行計(jì)算。通過(guò)本次實(shí)踐.我對(duì)Matlab的使用更加熟練了.解決問(wèn)題以及整理數(shù)據(jù)等能力有很大提高.希望我們?cè)趯W(xué)?？梢杂懈嗟臋C(jī)會(huì)實(shí)踐.有更多的收獲。參考文獻(xiàn)[1]陳鍾賢.計(jì)算物理學(xué).哈爾濱工業(yè)大學(xué)出版社.2001.3[2]楊振華.酈志新.數(shù)學(xué)實(shí)驗(yàn).科學(xué)出版社.2010.2[3]林亮.吳群英.數(shù)值分析方法與實(shí)驗(yàn)：基于MATLAB實(shí)現(xiàn).高等教育出版社.2012.9[4]李慶楊.王能超.易大義.數(shù)值分析.華中科技大學(xué)出版社.2006.7[5]鐘季康.鮑鴻吉.大學(xué)物理習(xí)題計(jì)算機(jī)解法—MATLAB編程應(yīng)用.機(jī)械工業(yè)出版社.2008.1[6]何紅雨.電磁場(chǎng)數(shù)值計(jì)算法與MATLAB實(shí)現(xiàn).華中科技大學(xué)出版社附錄1：plot([-0.2;0.2],[0;0],'color','y','linestyle','-','linewidth',10);g=9.8;l=1;theta0=pi/4;x0=l*sin(theta0);y0=(-1)*l*cos(theta0);axis([-0.75,0.75,-1.25,0]);axis('on');head=line(x0,y0,'color','r','linestyle','.','erasemode','xor','markersize',40);body=line([0;x0],[0;y0],'color','b','linestyle','-','erasemode','xor');t=0;dt=0.005;while1t=t+dttheta=theta0*cos(sqrt(g/l)*t);x=l*sin(theta);y=(-a)*l*cos(theta);set(head,'xdata',x,'ydata',y);set(body,'xdata',[0;x],'ydata',[0;y]);drawnow;end附錄2：a=1;b=2;f=@(x)x^3+4*x^2-10;c=(a+b)/2;whileabs(b-a)>1e-5iff(c)*f(b)<0a=c;elseb=c;endc=(a+b)/2;x=cendfprintf('\nx=%.5f,f(x)=%.5f\n',x,f(x));附錄3：FI=100;hx=23;hy=23;%定外框網(wǎng)絡(luò)；x軸1~23,1~23CX1=9;CX2=15;CY1=9;CY2=15;%定內(nèi)框位置：x軸9~15.y軸9~15DX=1+CX2-CX1;DY=1+CY2-CY1;v=ones(hy,hx)*50;%每個(gè)網(wǎng)格點(diǎn)的初值mesh1=ones(hy,hx)*2;%每個(gè)網(wǎng)格點(diǎn)的計(jì)算標(biāo)志值v(1,:)=zeros(1,hx);%第一行的初值（下邊界）mesh1(1,:)=zeros(1,hx);%第一行標(biāo)志值為零（下邊界）v(hy,:)=zeros(1,hx);%上邊界初值mesh1(hy,:)=zeros(1,hx);%上邊界標(biāo)志值為零v(:,1)=zeros(hy,1);%左邊界初值mesh1(:,1)=zeros(hy,1);%左邊界標(biāo)志值為零v(:,hx)=zeros(hy,1);%右邊界初值mesh1(:,hx)=zeros(hy,1);%右邊界標(biāo)志值為零v(CY1:CY2,CY1:CY2)=ones(DX,DY)*FI;%內(nèi)框區(qū)初值mesh1(CY1:CY2,CX1:CX2)=ones(DX,DY);%內(nèi)框區(qū)標(biāo)志值k=0;difmax=1.0;%最大誤差while(difmax>1.0e-6)k=k+1;%計(jì)算迭代次數(shù)difmax=0.0;fori=2:hy-1%從2到20行循環(huán)forj=2:hx-1%從2到20列循環(huán)m=mesh1(i,j);%?。╥,j）點(diǎn)標(biāo)志值if(m>=2)%標(biāo)志判斷vold=v(i,j);%取該點(diǎn)的原值v(i,j)=(1/4)*(v(i-1,j)+v(i,j-1)+v(i+1,j)+v(i,j+1));%拉普拉斯方程差分式dif=v(i,j)-vold;%前后兩次迭代值的差dif=abs(dif);%取絕對(duì)值if(dif>difmax)difmax=dif;e