编辑代码

	parameter(li=83,lj=35)
	real z(li,lj),u(li,lj),v(li,lj),hh(li,lj),u1(li,lj),v1(li,lj)
	g=9.8;h=10.;am=2.;dt=30.;dx=250.;rb2=0.0025
	k1=nint(3600/dt)
	kt4=4*12*k1
	km=5*12*k1
	wt=2*3.1415927/12./3600.*dt
	z=0.;u=0.;v=0.;hh=h
	open(1,file='t00a.txt')
	open(2,file='t00s.txt')
	do k=1,km
	do i=3,li-2,2
	do j=3,lj-2,2
    	uv=dt*rb2*sqrt(u(i,j)*u(i,j)+v(i,j)*v(i,j))/(hh(i,j)*hh(i,j))
	u(i,j)=u(i,j)-dt*g*hh(i,j)*(z(i+1,j-1)+z(i+1,j+1)-z(i-1,j-1)-z(i-1,j+1))/(4*dx)-uv*u(i,j)
	v(i,j)=v(i,j)-dt*g*hh(i,j)*(z(i-1,j+1)+z(i+1,j+1)-z(i-1,j-1)-z(i+1,j-1))/(4*dx)-uv*v(i,j)
	enddo
	enddo
	do i=2,li-3,2
	do j=2,lj-1,2
	z(i,j)=z(i,j)-dt*(u(i+1,j-1)+u(i+1,j+1)-u(i-1,j-1)-u(i-1,j+1)+v(i+1,j+1)+v(i-1,j+1)-v(i+1,j-1)-v(i-1,j-1))/(4*dx)
	enddo
	enddo
	do j=2,lj-1,2
	z(li-1,j)=am*sin(wt*k)
	enddo
	do i=3,li-2,2
	do j=3,lj-2,2
	hh(i,j)=h+(z(i-1,j-1)+z(i+1,j-1)+z(i+1,j+1)+z(i-1,j+1))/4
	enddo
	enddo
    	u1=u/hh
    	v1=v/hh
	if(mod(k,k1).eq.0) print *,'k=',k
	write(2,'(1x,i6,9f6.2)') k,z(82,18),z(42,18),z(2,18),u1(81,17),v1(81,17),u1(41,17),v1(41,17),u1(3,17),v1(3,17)
	if(k.ge.kt4.and.mod(k,k1).eq.0) then
	write(1,'(1x,2hK=,i5)') k
	write(1,'(1x,2hZ:,17i6)')(j,j=2,lj-1,2)
	do i=2,li-1,2
	write(1,'(1x,i2,17f6.2)') i,(z(i,j),j=2,lj-1,2)
	enddo
	write(1,'(1x,2hU:,16i6)')(j,j=3,lj-2,2)
	do i=3,li-2,2
	write(1,'(1x,i2,16f6.2)') i,(u1(i,j),j=3,lj-2,2)
	enddo
	write(1,'(1x,2hV:,16i6)')(j,j=3,lj-2,2)
	do i=3,li-2,2
	write(1,'(1x,i2,16f6.2)') i,(v1(i,j),j=3,lj-2,2)
	enddo
	endif
	enddo
	close(1)
	close(2)
	end