南邮应用物理计算物理实践.doc

《南邮应用物理计算物理实践.doc》由会员分享，可在线阅读，更多相关《南邮应用物理计算物理实践.doc（30页珍藏版）》请在三一办公上搜索。

1、标准文档实 验 报 告课程名称： 计算物理实践 专 业： 应用物理学 学 号： 姓 名： 完成日期： 年 月 目 录一、简单物理实验的模拟及实验数据处理1问题描述：1单摆运动原理1模型的建立1流程图21.5 Matlab程序设计仿真2程序3单摆演示截图3二、方程组的数值解法4问题描述：4原理分析4二分法理论42.2.2 分析求解42.3 Matlab程序52.4 Matlab程序运行结果：5三、静电场问题的计算7问题描述：7原理分析71简单迭代法7有限差分：83.2.3 解题过程9程序设计仿真9仿真结果10四 、 热传导方程和波动方程的差分解法11问题描述11原理分析11具体步骤134.4 M

2、ATLAB程序设计仿真13程序运行结果13结束语15参考文献16附录1：17附录2：18附录3：19附录4：21实用文案标准文档一、简单物理实验的模拟及实验数据处理：编写单摆运动演示程序。在不考虑空气阻力和很小的假设下，单位质量小球做理想简谐运动，此时。取g=9.8，L=1, =pi/4.设在某一时刻，单摆的摆线偏离垂直线的角位移为，将重力mg 分解为径向力 F 和切向力T ，则T的大小为mg sin,切向加速度为=L.根据牛顿第二定律得方程ma=m=mg从而单摆运动的微分方程为 =由于=+.当很小时， 所以单摆的微分方程可表示为= 上式表明，当很小时，单摆的角加速度与角位移成正比，但方向相反

3、，且方程的解可表示为=A 建立物理模型，假设单摆运动过程中的小球中心点的坐标为(x, y)。根据几何关系x, y满足如下关系：在不考虑空气阻力和很小的假设下= 与式比较知 A= 所以，单摆的运动方程为将=pi/4,g=9.8,L=1带入得，开 始输入q5= t=1?q YY=-L*;X=L*;= N结束t=t+dt图1.1 程序流程图 Matlab程序设计仿真通过set函数将变现和小球的图像句柄，加入X、Y的动态参量，它们得轨迹变化即为上面分析的轨迹方程。用line函数画出初始的位置，并将句柄分别给sphere、lp，然后每隔dt时间刷新一次。尝试对dt的设置发现，值在的效果比较好。Matla

4、b程序程序见附录1。图1.2 单摆演示图二、方程组的数值解法2.1问题描述：二分法求解方程x3+4x2100在区间1，2内的根，精度自设。二分法理论f(x)单调， f(a)f(b)0f(x)=0在（a,b）有唯一根。设f(x)在a,b上连续，f(x)=0在a,b上存在唯一解，且f(a)f(b)0,记第一步，计算f()f().若f()f()0,则x*,记，,否则x*,记，,对两种情形均有x*,记第k步，计算f()f()0,则x*,记，,否则x*,对这两种情形均有x* ,记 ,k=1,2,.x*且所以数列满足x*- |()=)=.=)=即从而当k充分大，且可由x*- 控制精度。2.2.2 分析求解

5、令f(x)=-x3+4x210,f(x)在1，2上连续，且f(1)f(2) 0,则f(x)=0在1，2上有唯一解，记a=1,b=2,x=1.5.然后计算f(a)f(b),若f(a)f(b) 1e-5if f(c)*f(b)1.0e-6)k=k+1; %计算迭代次数difmax=0.0;for i=2:hy-1 %从2到20行循环for j=2:hx-1 %从2到20列循环m=mesh1(i,j); %取（ i, j）点标志值if(m=2) %标志判断vold=v(i,j); %取该点的原值v(i,j)=(1/4)*(v(i-1,j)+v(i,j-1)+v(i+1,j)+v(i,j+1); %拉

6、普拉斯方程差分式dif=v(i,j)-vold; %前后两次迭代值的差dif=abs(dif); %取绝对值if(difdifmax) difmax=dif; end %所有网格中取最大误值endend endendsubplot(1,2,1),mesh(v)%画三维曲面图axis(-2,hx+3,-2,hy+3,0,100)subplot(1,2,2),contour(v,13)%画等电位线图hold onx=1:1:hx;y=1:1:hy;xx,yy=meshgrid(x,y);%形成栅格Gx,Gy=gradient(v,0.6,0.6);%计算梯度%quiver(xx,yy,Gx,Gy,

7、-0.5,r)%根据梯度数画箭头axis(-2,hx+3,-2,hy+3)plot(1,1,hx,hx,1,1,hy,hy,1,1,k)plot(CX1,CX1,CX2,CX2,CX1,CY1,CY2,CY2,CY1,CY1,k)%画外框边线text(CX1+0.6,CY1+(CY2-CY1)/2,U=100,fontsize,10);%画内框边界线text(hx/2,hy+1,U=0,fontsize,10);%外框上边界标注text(hx/2,0,U=0,fontsize,10);%外框下边界标注text(-1.7,hy/2,U=0,fontsize,10);%外框左边界标注text(hx

8、+0.4,hy/2,U=0,fontsize,10);%外框右边界标注hold off附录4：h=0.1;for k=1:37 %从1到37循环37次u(1,k)=0;N=10;a=1/6 u(11,k)=1;for i=2:10 %从2到10共9步 u(i,1)=i*h*i*h; for k=1:36 for i=2:10 u(i,k+1)=a*u(i+1,k)+(1-2*a)*u(i,k)+a*u(i-1,k); end end end endmesh(u)%画图xlabel(X Axis),ylabel(Y Axis),zlabel(Temperature),title(Thermal Field Distribution)实用文案