资源描述:
《MATLAB水波模拟源代码.doc》由会员上传分享,免费在线阅读,更多相关内容在教育资源-天天文库。
1、functionwaterwave%WATERWAVE%2DShallowWaterModel%%Lax-Wendrofffinitedifferencemethod.%Reflectiveboundaryconditions.%Randomwaterdropsinitiategravitywaves.%Surfaceplotdisplaysheightcoloredbymomentum.%Plottitleshowst=simulatedtimeandtv=ameasureoftotalvariation.%An
2、exactsolutiontotheconservationlawwouldhaveconstanttv.%Lax-Wendroffproducesnonphysicaloscillationsandincreasingtv.%%Information:%Author:cuixing%qq:%email:@qq.com%Parametersn=64;%gridsizeg=9.8;%gravitationalconstantdt=0.01;%hardwiredtimestepdx=1.0;dy=1.0;nplotst
3、ep=8;%plotintervalndrops=5;%maximumnumberofdropsdropstep=500;%dropintervalD=droplet(1.5,21);%simulateawaterdrop%Initializegraphics[surfplot,top,start,stop]=initgraphics(n);%Outerloop,restarts.whileget(stop,'value')==0set(start,'value',0)H=ones(n+2,n+2);U=zeros
4、(n+2,n+2);V=zeros(n+2,n+2);Hx=zeros(n+1,n+1);Ux=zeros(n+1,n+1);Vx=zeros(n+1,n+1);Hy=zeros(n+1,n+1);Uy=zeros(n+1,n+1);Vy=zeros(n+1,n+1);ndrop=ceil(rand*ndrops);nstep=0;%Innerloop,timesteps.whileget(start,'value')==0&&get(stop,'value')==0nstep=nstep+1;%Randomwat
5、erdropsifmod(nstep,dropstep)==0&&nstep<=ndrop*dropstepw=size(D,1);i=ceil(rand*(n-w))+(1:w);j=ceil(rand*(n-w))+(1:w);H(i,j)=H(i,j)+rand*D;end%ReflectiveboundaryconditionsH(:,1)=H(:,2);U(:,1)=U(:,2);V(:,1)=-V(:,2);H(:,n+2)=H(:,n+1);U(:,n+2)=U(:,n+1);V(:,n+2)=-V(
6、:,n+1);H(1,:)=H(2,:);U(1,:)=-U(2,:);V(1,:)=V(2,:);H(n+2,:)=H(n+1,:);U(n+2,:)=-U(n+1,:);V(n+2,:)=V(n+1,:);%Firsthalfstep%xdirectioni=1:n+1;j=1:n;%heightHx(i,j)=(H(i+1,j+1)+H(i,j+1))/2-dt/(2*dx)*(U(i+1,j+1)-U(i,j+1));%xmomentumUx(i,j)=(U(i+1,j+1)+U(i,j+1))/2-...
7、dt/(2*dx)*((U(i+1,j+1).^2./H(i+1,j+1)+g/2*H(i+1,j+1).^2)-...(U(i,j+1).^2./H(i,j+1)+g/2*H(i,j+1).^2));%ymomentumVx(i,j)=(V(i+1,j+1)+V(i,j+1))/2-...dt/(2*dx)*((U(i+1,j+1).*V(i+1,j+1)./H(i+1,j+1))-...(U(i,j+1).*V(i,j+1)./H(i,j+1)));%ydirectioni=1:n;j=1:n+1;%heigh
8、tHy(i,j)=(H(i+1,j+1)+H(i+1,j))/2-dt/(2*dy)*(V(i+1,j+1)-V(i+1,j));%xmomentumUy(i,j)=(U(i+1,j+1)+U(i+1,j))/2-...dt/(2*dy)*((V(i+1,j+1).*U(i+1,j+1)./H(i+1,j+1))-...(V(i+1,j).*U(i+1,j)