Matlab在物理上的应用举例

1. 单列波%%单列波t=0:0.001:10;A=input('振幅A=');w=input('频率w=');a=input('a=');y=A.*sin(w.*t+a);plot(t,y);pause(1),sound(y);ylabel('y'),xlabel('t')2. %%光的单缝衍射现象Lambda=500e-9; %a=input('a='); % 可取0.2e-3,1e-3,2e-3三种情况z=1 %ymax=3*Lambda*z/a; %Ny=51; %ys=linspace(-ymax,ymax,Ny); %NPoints=51; %yPoint=linspace(-a/2,a/2,NPoints); %for j=1:Ny %L=sqrt((ys(j)-yPoint).^2+z^2); %Phi=2*pi.*(L-z)./Lambda; %SumCos=sum(cos(Phi)); %SumSin=sum(sin(Phi)); %B(j)=(SumCos^2+SumSin^2)/NPoints^2; %endclf,plot(ys,B,'*',ys,B);grid; %3. %%用毕奥－沙伐尔定律计算电流环产生的磁场mu0=4*pi*1e-7;I0=5.0;Rh=1;C0=mu0/(4*pi)*I0;NGx=21;NGy=21;x=linspace(-Rh,Rh,NGx);y=linspace(-3,3,20);y=x;Nh=20;theta0=linspace(0,2*pi,Nh+1);theta1=theta0(1:Nh);y1=Rh*cos(theta1);z1=Rh*sin(theta1);theta2=theta0(2:Nh+1);y2=Rh*cos(theta2);z2=Rh*sin(theta2);xc=0;yc=(y2+y1)/2;zc=(z2+z1)/2;for i=1:NGyfor j=1:NGxrx=x(j)-xc;ry=y(i)-yc;rz=0-zc;r3=sqrt(rx.^2+ry.^2+rz.^2).^3;dlXr_x=dly.*rz-dlz.*ry;dlXr_y=dly.*rx-dlx.*rz;Bx(i,j)=sum(C0*dlXr_x./r3);By(i,j)=sum(C0*dlXr_y./r3); endendclf;quiver(x,y,Bx,By);4. %%多普勒效应x0=500;v=50;y0=20;c=330;w=1000;t=0:0.001:30;r=sqrt((x0-v*t).^2+y0.^2);t1=t-r/c;u=sin(w*t)+sin(1.1*w*t);u1=sin(w*t1)+sin(1.1*w*t1); sound(u);pause(5);sound(u1);5.亥姆霍兹线圈clear allmu0=4*pi*1e-7;I0=5.0;Rh=1;C0=mu0/(4*pi)*I0;NGx=21;NGy=21;x=linspace(-Rh,Rh,NGx);y=linspace(-Rh,Rh,NGy);Nh=20;theta0=linspace(0,2*pi,Nh+1); theta1=theta0(1:Nh);y1=Rh*cos(theta1);z1=Rh*sin(theta1);theta2=theta0(2:Nh+1);y2=Rh*cos(theta2);z2=Rh*sin(theta2);xc=0;yc=(y2+y1)/2;zc=(z2+z1)/2;for i=1:NGyfor j=1:NGxrx=x(j)-xc;ry=y(i)-yc;rz=0-zc;r3=sqrt(rx.^2+ry.^2+rz.^2).^3;dlXr_x=dly.*rz-dlz.*ry;dlXr_y=dly.*rx-dlx.*rz;Bx(i,j)=sum(C0*dlXr_x./r3);By(i,j)=sum(C0*dlXr_y./r3);endendBax=Bx(:,11:21)+Bx(:,1:11);Bay=By(:,11:21)+By(:,1:11);subplot(1,2,1)mesh(x(11:21),y,Bax);xlabel('x');ylabel('y'); subplot(1,2,2),plot(y,Bax),grid,xlabel('y');ylabel('Bx'); ６．库仑引力clear all;N=input('电荷数目N：');for ic=1:N %fprintf('-----\n 对电荷#%g\n',ic);rc=input('电荷位置[x y]（米）：');x(ic)=rc(1); %y(ic)=rc(2); %q(ic)=input('输入电荷量（库仑）'); endE0=8.85e-12; %C0=1/(4*pi*E0); %for ic=1:NFx=0.0; Fy=0:0; %for jc=1:N %if(ic~=jc) %xij=x(ic)-x(jc); yij=y(ic)-y(jc);Rij=sqrt(xij^2+yij^2); %Fx=Fx+C0*q(ic)*xij/Rij^3;Fy=Fy+C0*q(ic)*yij/Rij^3;endendfprintf('其它电荷作用在电荷#%g上的合力为：\n',ic); %fprintf('x－分量：%gN\n',Fx);fprintf('y－分量：%gN\n',Fy)end７．李萨如图形% lisaru.msyms t a1 a2 w1 w2x=cos(w1.*t+a1);y=sin(w2.*t+a2);a1=input('a1=');a2=input('a2=');w1=input('w1=');w2=input('w2=');tf=10;Ns=1000;t=linspace(0,tf,Ns);dt=tf/(Ns-1); %分Ns个点,求出时间增量dtxplot=eval(x);yplot=eval(y); %计算Ns个点的位置x(t),y(t)figure(gcf);subplot(1,2,1),for i=1:750plot(yplot(1:i),xplot(1:i)); %画点的轨迹图axis('equal'); grid ; %使两轴比例相同pause(0.01)end８．耦合振子m1=2;m2=2;K1=16;K2=4; % x0=[1;0];xd0=[6;6];tf=10; %M=[m1,0;0,m2];K=[K1+K2,-K2;-K2,K1+K2] ; %u(:,s)�[u,L]=eig(K,M) ; %t=linspace(0,tf,101);x=zeros(2,101); % for s=1:2 %alfa=sqrt(u(:,s)'*M*u(:,s)) ; %u(:,s)= u(:,s)/alfa; %w(j)=sqrt(L(j,j)); % xt=u(:,j)*(u(:,j)'*M*x0*cos(w(j)*t)+u(:,j)'*M*xd0/w(j)*sin(w(j)*t));x=x+xt; %endfor r=1:2 %x1,x1?? subplot(2,1,r)plot(t.x(r,:)),grid;xlabel('xxx');ylabel(['yyy',num2str(r)]);end９．拍频%%t=0:0.001:10;a1=input('??1=');w1=input('??1=');a2=input('??2=');w2=input('??2=');y1=a1*sin(w1*t);y2=a2*sin(w2*t);y=y1+y2;subplot(3,1,1),plot(t,y1),ylabel('y1')subplot(3,1,2),plot(t,y2),ylabel('y2')subplot(3,1,3),plot(t,y),ylabel('y'),xlabel('t')pause,sound(y1);pause(5),sound(y2);pause(5),sound(y),pausesubplot(1,1,1)１０。

物块下滑m1=input('m1=');m2=input('m2=');theta=input('theta(度)=');theta=theta*pi/180; g=9.81;A=[m1*cos(theta),-m1,-sin(theta),0;...m1*sin(theta),0,cos(theta),0;...0,m2, -sin(theta),0;...0,0,-cos(theta),1];B=[0,m1*g,0,m2*g]'; X=A\B;a1=X(1), a2=X(2), N1=X(3),N2=X(4)１１。

循环%%xunhuanR=8.31;gama=1.4; %(注:由于在MATLAB中无拉丁文,所以用gama代替?)nMoles=0.5;P(1)=1e5;V(1)=0.012;WTotal=0;QTotal=0;iPoint=1;NCurve=100;PPlot=P(1);VPlot=V(1);%变量和图形初始化，输入气体的摩尔数nMoles，初始压力P（1），初始容积V（1），%气体常数R＝8.314；给定起始总功WTotal=0;点序号iPoint=1;画等温线用%的点数NCurve=100;P－V图第一点坐标PPlot=P(1);VPlot=V(1)T(1)=P(1)*V(1)/(nMoles*R); %算出初始温度%为了进入循环，先要设两个不相等的PathType和QuitType值QuitType=5;PathType=0;while(PathType~=QuitType) %在菜单上选择‘退出’之前不断循环，%选择路径类型或退出iPoint=iPoint+1; %下一点fprintf('对过程＃％g\n',iPoint-1);PathType=menu(sprintf('过程％g：选择下一路径',iPoint-1),...'等压','等容','等温','绝热','退出'); %图形界面菜单生成语句switch PathTypecase 1 %等压路径V(iPoint)=input('输入新容积：');P(iPoint)=P(iPoint-1); %压力不变T(iPoint)=P(iPoint)*V(iPoint)/(nMoles*R); %按新容积算出温度W=P(iPoint)*(V(iPoint)-V(iPoint-1)); %计算等压过程所做的功Q=(gama*nMoles*R/(gama-1))*(T(iPoint)-T(iPoint-1));VPlot=[VPlot,V(iPoint)]; %加上新的容积和压力点，用以绘图PPlot=[PPlot,P(iPoint)];case 2 %等容路径P(iPoint)=input('输入新压力：');V(iPoint)=V(iPoint-1); %容积不变T(iPoint)=P(iPoint)*V(iPoint)/(nMoles*R); %按新压力算出温度W=0; %等容路径上所做的功为零Q=(nMoles*R/(gama-1))*(T(iPoint)-T(iPoint-1));VPlot=[VPlot,V(iPoint)]; %加上绘图用的新容积和压力点PPlot=[PPlot,P(iPoint)];case 3 %等温路径V(iPoint)=input('输入新容积：');T(iPoint)=T(iPoint-1); %温度不变P(iPoint)=nMoles*R*T(iPoint)/V(iPoint); %按新容积求新压力W=nMoles*R*T(iPoint)*log(V(iPoint)/V(iPoint-1)); %求所做的功Q=W;%用元素群运算求等温路径上的P和V，加进绘图数据中VNew=linspace(V(iPoint-1),V(iPoint),NCurve);PNew=nMoles*R*T(iPoint)./VNew;VPlot=[VPlot,VNew]; %将新的V，P点加入绘图数据中PPlot=[PPlot,PNew];case 4 %绝热路径V(iPoint)=input('输入新容积:');P(iPoint)=P(iPoint-1)*V(iPoint-1)^gama/V(iPoint)^gama; %按新容积求新压力T(iPoint)=T(iPoint-1)*V(iPoint-1)^(gama-1)/V(iPoint)^(gama-1); %按新容积求新温度W=-(nMoles*R/(gama-1))*(T(iPoint)-T(iPoint-1)); %求所做的功Q=0;VNew=linspace(V(iPoint-1),V(iPoint),NCurve);PNew=P(iPoint-1)*V(iPoint-1)^gama./VNew.^gama;VPlot=[VPlot,VNew]; %将新的V，P点加入绘图数据中PPlot=[PPlot,PNew];otherwiseend%画出到目前为止的PV图if(PathType~=QuitType)WTotal=WTotal+W; %将新做的功加进总功if(Q>0)QTotal=QTotal+Q;endfigure(gcf);plot(V,P,'o',VPlot,PPlot,'-')%图形窗移前，绘图size=axis;axis([size(1)*0.9,size(2)*1.1,size(3)*0.9,size(4)*1.1]);%上两句用于美观图象，使曲线不要紧靠边框%标注语句略endendWTotalQTotaly=WTotal/QTotal１２。