卢瑟福散射演示动画

卢瑟福散射演示动画动力学方程：22222212222222212211yxyyxmeZZdtydyxxyxmeZZdtxd，初始条件：;|,0|;0|,|00000dyxdtdyvdtdxtttt用matlab求解微分方程前，进行转变：22222222222211yxyyxCdtydyxxyxCdtxd；令;42;31;2,1xyyxydtdyxxxdtdx；则动力学方程写为：;34;4)42(3;12;2)42(12/3222/322xdtdxxxxCdtdxxdtdxxxxCdtdx初始条件为：;)0(4;10)0(2;0)0(3;)0(10dxxxvx积分时间区段取：;25,00tftfigure('name','散射演示');%设置标题名字pausetime=.002;%设置暂停时间set(gca,'xlim',[-1520],'ylim',[-1818]);%设置图片显示坐标范围set(gcf,'doublebuffer','on')%消除抖动%axisequalholdonplot(0,0,'color','k','marker','.','markersize',40);%画靶原子line([-10-10],[-1010],'LineStyle','-');%入射粒子初始位置范围u=0:0.01*pi:2*pi;X=14*cos(u);Y=14*sin(u);plot(X,Y,'-k');%环形探测器len=0;forn=1:100%入射粒子数量v0=1;%入射粒子的初速度dr=rand;d=20*dr-10;%入射粒子的位置pp=plot(0,d,'r','marker','o','markersize',5);%入射粒子的初始位置t0=0;tf=25;x0=[v0,-10,0,d]';[t,x]=ode45('xdot',t0,tf,x0);len=length(t);form=1:lenset(pp,'xdata',x(m,2),'ydata',x(m,4));%入射粒子的即时位置plot(x(m,2),x(m,4),'o','markersize',1);ifm==lenplot(x(m,2),x(m,4),'ok','markersize',5);end%把入射粒子最后停留的位置记录下来pause(pausetime);%暂停一会drawnowendend%subplot(1,3,1);%plot(t,x(:,1),':b',t,x(:,2),'-r');%subplot(1,3,2);%plot(t,x(:,3),':b',t,x(:,4),'-r');%subplot(1,3,3);%plot(x(:,2),x(:,4),'-b');求解动力学方程组所用的函数。xdotfunctionxd=xdot(t,x)C=1;xd=zeros(4,1);xd(1)=C*((x(2))^2+(x(4))^2)^(-1.5)*x(2);xd(2)=x(1);xd(3)=C*((x(2))^2+(x(4))^2)^(-1.5)*x(4);xd(4)=x(3);可调整参数为：1入射粒子数量n2入射粒子初速度v03入射粒子位置d4入射粒子动力学参数C散射分布结果如图：取[nv0dC]=[100120*dr-101]可以看出散射角大致分布区域，以及存在大角度散射。甚至存在沿原路返回的入（出）射粒子。对于两个靶核的散射：增大两个靶核间距，使得入射粒子可以从中间穿过。多靶体系（两个靶原子）动力学方程：2222202222202222222022222022)()(1)()(1)()(1)()(1dyxdydyxmeZZdyxdydyxmeZZdtyddyxxdyxmeZZdyxxdyxmeZZdtxdBABA初始条件不变，变换后的动力学方程组为：;34);4())4(2(2)4())4(2(13;12;2))4(2(22))4(2(112/3222/3222/3222/322xdtdxdxdxxCdxdxxCdtdxxdtdxxdxxCxdxxCdtdx两个靶原子的情况：figure('name','2靶核散射演示');%设置标题名字pausetime=.002;%设置暂停时间set(gca,'xlim',[-1520],'ylim',[-1818]);set(gcf,'doublebuffer','on')%消除抖动%axisequalholdonplot(0,4,'color','k','marker','.','markersize',40);%画靶原子Aplot(0,-4,'color','k','marker','.','markersize',40);%画靶原子Bline([-10-10],[-1010],'LineStyle','-');u=0:0.01*pi:2*pi;X=14*cos(u);Y=14*sin(u);plot(X,Y,'-k');len=0;forn=1:100v0=1;%入射粒子的初速度dr=rand;d=20*dr-10;%入射粒子的位置pp=plot(0,d,'r','marker','o','markersize',5);%入射粒子的初始位置t0=0;tf=25;x0=[v0,-10,0,d]';[t,x]=ode45('xdo',t0,tf,x0);len=length(t);form=1:lenset(pp,'xdata',x(m,2),'ydata',x(m,4));%入射粒子的即时位置plot(x(m,2),x(m,4),'o','markersize',1);ifm==lenplot(x(m,2),x(m,4),'ok','markersize',5);end%把入射粒子最后停留的位置记录下来pause(pausetime);%暂停一会drawnowendend%subplot(1,3,1);%plot(t,x(:,1),':b',t,x(:,2),'-r');%subplot(1,3,2);%plot(t,x(:,3),':b',t,x(:,4),'-r');%subplot(1,3,3);%plot(x(:,2),x(:,4),'-b');用于求解微分方程组的函数：functionxdd=xdo(t,x)C1=1;C2=1;dd=4;xdd=zeros(4,1);xdd(1)=C1*((x(2))^2+(x(4)-dd)^2)^(-1.5)*x(2)+C2*((x(2))^2+(x(4)+dd)^2)^(-1.5)*x(2);xdd(2)=x(1);xdd(3)=C1*((x(2))^2+(x(4)-dd)^2)^(-1.5)*(x(4)-dd)+C2*((x(2))^2+(x(4)+dd)^2)^(-1.5)*(x(4)+dd);xdd(4)=x(3);