FORTRAN数值方法及其在物理学中应用2

1第二章物理图形、图象与计算机模拟§2.1简谐振动及其合成曲线模拟§2.2阻尼运动和阻尼振动的模拟§2.3驻波的模拟§2.5波的干涉和衍射图形模拟§2.4点电荷与点电荷系的等势线和电场模拟2§2.1简谐振动及其合成曲线模拟一、简谐振动的曲线和曲线txtv例1：画出曲线及对应的曲线，其中)cos(tAxtv0,,0.2A周期）(0)2cos()2cos()sin(TttvtAtAdtdxvm等分）秒，取NNTtT(/2/23Implicitreal*8(a-h,o-z)open(1,file=’x-t.dat’)open(2,file=’v-t.dat’)write(*,*)’inputA,w,phi,N’read(*,*)A,w,phi,Npi=3.1415926do10I=1,Nt=2.*pi/wt=t*float(I)/Nx=A*cos(w*t+phi)v=A*w*cos(w*t+phi+pi/2.)write(1,*)t,x10write(2,*)t,vend计算程序cos()xAtcos()2vAt40.00.20.40.60.81.01.21.41.61.82.0-8-6-4-202468tx=Acos(*t+)v=Acos(*t++/2)图形模拟Origin简介5二、简谐振动的合成1.同方向简谐振动的合成同频率情况)cos(111tAx)cos(222tAx)cos(21tAxxx其中2/112212221)]cos(2[AAAAA11221122sinsintgcoscosAAAA减弱（反相））（加强（同相）,12|,|,2,12211221KAAKAAA6例2：试给出两个同方向同频率简谐振动的合成程序open(1,file=’x1.dat’)open(2,file=’x2.dat’)open(3,file=’x.dat’)write(*,*)’inputA1,A2,w,phi1,phi2=?’read(*,*)A1,A2,w,phi1,phi2pi=3.1415926do10I=1,1000t=2.*pi/wt=t*float(I)/1000x1=A1*cos(w*t+phi1)x2=A2*cos(w*t+phi2)x=x1+x2write(1,*)t,x1write(2,*)t,x210write(3,*)t,xend)cos(111tAx)cos(222tAx70.000.050.100.150.20-0.20-0.15-0.10-0.050.000.050.100.150.20tx1=0.05cos(10t+0.6)x2=0.1cos(10t+0.6)x=x1+x2(同相)0.000.050.100.150.20-0.15-0.10-0.050.000.050.100.15tx1=0.05cos(10t+0.6)x2=0.1cos(10t+1.6)x=x1+x2(反相)图形模拟两个同方向同频率简谐振动的合成80.000.050.100.150.20-0.15-0.10-0.050.000.050.100.15tx1=0.05cos(10t+0.6)x2=0.1cos(10t+0.2)x=x1+x2EX2-1:编程完成例2。两个同方向同频率简谐振动的合成9不同频率情况若，会出现拍的现象。21,)cos(1111tAx)cos(2222tAx合振动不再是简谐振动，利用旋转矢量法可以求得合振动的振幅为2/11212212221)]}()cos[(2{tAAAAA||212中间经历的时间称为周期，显然，频率：||2||121221AAA||21AAA振幅在和间周期性地变化，属振动调制。合振动振幅从一次极大到相邻的另一次极大。10假设两个分振动振幅都为，圆频率相差较小，取它们的初相位1A21,tAx111costAx212cos此时合成运动的位移可写成：])(21cos[])(21cos[21212121ttAxxx])(21cos[12tAx变化主要取决于，振幅按变化。])(21cos[12t])(21cos[2121tA都是零，则可以分别表示为：])(21cos[2121tAA其中。)(2112)(2112由于圆频率远大于圆频率，110510152025303540-2.5-2.0-1.5-1.0-0.50.00.51.01.52.02.5txx=x1+x2=cos10t+cos11t图形模拟两个同方向频率近似的简谐振动的合成122.两个相互垂直方向简谐振动的合成111222cos()cos()xAtyAt若，则有合振动方程：21)(sin)cos(21221221222212AAxyAyAxxAAyAyAx1221120:0①（一、三象限直线方程）xAAyAyAx1221120:（二、四象限直线方程）②（椭圆方程）1:222221212AyAx13质点轨迹曲线下图所示为两个频率相同、振幅相等、相互垂直而相位差-1.0-0.50.00.51.0-1.0-0.50.00.51.02-1=xy-1.0-0.50.00.51.0-1.0-0.50.00.51.02-1=xy-1.0-0.50.00.51.0-1.0-0.50.00.51.02-1=xy两个频率相同、振幅相等、相互垂直简谐振动的合成4/3,2/,4/为下的质点轨迹曲线。14若，但满足一定整数倍数比关系时，则会21下图所示为两个频率不同（满足）、振幅相等、2/3/21-1.0-0.50.00.51.0-1.0-0.50.00.51.02-1=xy-1.0-0.50.00.51.0-1.0-0.50.00.51.02-1=xy-1.0-0.50.00.51.0-1.5-1.0-0.50.00.51.01.52-1=0xy利萨如图形示意图出现利萨如图形：4/,8/相互垂直而相位差为0，下的质点轨迹曲线。15§2.2阻尼运动和阻尼振动的模拟一、阻尼情况下物体运动的曲线tv例3：质量为的摩托快艇以速度行驶，它受到的摩擦阻力与速度成正比，设比例系数为，则，试求关闭发动机后，对的变化规律。（取m0vKKvFvt)/004.0250/1,/100skgK/msmv解：物理分析与数学模型)(or00mKAevvevvAttmKKvdtdvmmaFvmKdtdv可用函数作图法方法1:16方法2：用（差商法）dtdvtv替代AvvmKtvdtdvtAvv)(1tAvvvii函数近似值作图法17open(1,file=’vt.dat’)write(*,*)’inputa,v0,t=?’read(*,*)a,v0,tv1=v0t0=0.v10=v0write(1,*)t0,v0,v10dt=t/1000.do10j=1,1000tt=t*float(j)/1000.v=v0*exp((-1.)*a*tt)v1=v1-a*v1*dt10write(1,*)tt,v,v1end05001000150020000246810K/m=1/250,v0=10m/s阻尼运动曲线tv模拟程序1()iivvAvt0()AtKvveAm18二、阻尼振动问题：弹簧振子阻尼振动的方程为—阻尼因子—弹簧振子的角频率)(0202022xdtdxdtxd0试用函数近似法作出位移与时间的函数变化曲线。19问题分析:解：将二阶微分方程化为一阶微分方程22dtxddtdvdtdxv202(,)dvvxfxvdttvxfv*),(tvxfvviii*),(1即tvx*而tvxxiii*1∴tv20实例说明：例4：画出当，，时，，，，秒下的曲线。1.0100t0.1x0.0v250t1000/25ttx计算程序：open(1,file='v-t.dat')open(2,file='x-t.dat')write(*,*)'inputB,w0,x0,v0,t=?'read(*,*)B,w0,v0,x0,tdt=t/1000.v=v0x=x0tt0=0.0write(1,*)tt0,v0write(2,*)tt0,x0do10j=1,1000tt=float(j)*dtf=-2.*B*v-w0**2*xv=v+f*dtx=x+v*dtwrite(1,*)tt,v10write(2,*)tt,xend210510152025-1.0-0.50.00.51.0t(s)Yx-t曲线v-t曲线阻尼振动曲线示意图图形模拟220510152025-1.0-0.50.00.51.0tx=0.1,=1,小阻尼情况=1,=1，临界阻尼情况=1.5,=1，大阻尼情况不同阻尼情况振动曲线示意图图形模拟23EX2-2:编程完成例3。BvAdtxd22EX2-3:一石子从空中静止下落，已知BA,tx式中为常数，试绘制石子下落的,8.9A,5.0B,10|0tx曲线。其中,0|0tv20t。作业24§2.3驻波的模拟定义：两列振幅、振动方向和频率都相同而传播方向相反的两列同类波相干叠加形成驻波。设有两列振动方向相同、振幅相同、频率相同的平面余弦波，,)(2cos1xtAy)(2cos2xtAy按叠加原理，合成的驻波的波函数为:)](2cos)(2[cos21xtxtAyyytxAy2cos2cos2轴的正、负方向传播。X分别沿25在值满足下式的各点，振幅为零－驻波波节处x,2,1,0,4)12(2)12(2kkxkx相邻两波节的距离为半波长，即：24)12(41)1(21kkxxkk讨论:tt2cos因子是时间的余弦函数，说明形成驻波后，各质点都在作同频率的谐振动。另一因子xA2cos2是坐标的余弦函数，说明各质点的振幅按余弦函数规律分布。txAy2cos2cos226x在值满足下式的各点，振幅最大－驻波波腹处,2,1,0,22kkxkx相邻两波腹间的距离也为半波长，即：222)1(1kkxxkkl波节处的质点振动的振幅为零，始终处于静止；l波腹处的质点振动的振幅最大，等于。l其他各处质点振动的振幅则在零与最大之间。l两相邻波节或两相邻波腹之间相距半波长。l波腹和相邻波节间的距离为，波腹和波节交替作等距离排列。A24/特征:讨论:27开始输入A,T和t=0,1,2,...Nx=0,1,2,...M计算y1=Acos2(t/T-x/)和y2=Acos2(t/T+x/)计算y=y1+y2画点（t,x）结束驻波模拟程序流图28implicitreal*8(a-h,o-z)open(1,file='zhubo.dat')open(2,file='zhubo1.dat')open(3,file='zhubo2.dat')write(*,*)'inputA,x,wavelength'read(*,*)A,x,wapi=3.1415926T=2.*pitime=0.5*Tdx=x/1000do10I=1,1000x=dx*floa