补充fortran和grads

HFTq

数值预报生产实习

HFTq

c     请补充计算地图投影放大系数mm 和科氏参数 ff子程序
c     mm:the amplification factor of map projection
c     ff:coriolis parameter
        subroutine mf(mm,ff,m0,n0)
                real mm(m0,n0),ff(m0,n0),vl(m0,n0)
                vk=0.7156
                vle=11423370
                va=6371000
                omege=7.292e-5
                vd=300000
                in=-4
                jn=25
                do j=1,n0
                        do i=1,m0
                                vl(i,j)=sqrt(((in-1+i)*vd)**2+((jn+1-j)*vd)**2)
                                mm(i,j)=(vk*vl(i,j))/(va*sqrt(1-((vle**(2/vk)-vl(i,j)**(2/vk))/(vle**(2/vk)+vl(i,j)**(2/vk)))**2))
                                ff(i,j)=2*omege*((vle**(2/vk)-vl(i,j)**(2/vk))/(vle**(2/vk)+vl(i,j)**(2/vk)))
                        end do
          end do

c-----------------------------------------------------------------c
c-----------------------------------------------------------------c
c     请补充计算绝对涡度平流子程序
        subroutine abv(av,af,h,mm,ff,m0,n0)
                real av(m0,n0),af(m0,n0),hh(m0,n0),mm(m0,n0),ff(m0,n0)
                real,parameter::d=300000,f=1.0e-4,g=9.8
                integer i,j
                do j=2,n0-1
                        do i=2,m0-1
                                av(i,j)=(g*mm(i,j))/(d**2.0*f)*(hh(i+1,j)+hh(i,j+1)+hh(i-1,j)+hh(i,j-1)-4*hh(i,j))+ff(i,j)
                        end do
        end do
                do j=3,n0-2
                        do i=3,m0-3
                                af(i,j)=-0.25*((hh(i+1,j)-hh(i-1,j))*(av(i,j+1)-av(i,j-1))-(hh(i,j+1)-hh(i,j-1))*(av(i+1,j)-av(i-1,j)))
                        end do
                end do
                end

c-----------------------------------------------------------------c
c-----------------------------------------------------------------c
c     请补充计算位势高度倾向子程序
        subroutine pzpt(zt,af,m0,n0)
                parameter(rer=1.5,rps=1.0e-5)
                dimension
                rn(m0,n0),zt(m0,n0),f(m0,n0),hh(m0,n0),af(m0,n0)
                do i=1,m0
                        do j=1,n0
                                hh(i,j)=0.0
                                zt(i,j)=hh(i,j)
                        end do
                end do
                do i=4,m0-3
                        do j=4,n0-3
                                rn(i,j)=zt(i+1,j)+zt(i,j+1)+zt(i-1,j)+zt(i,j-1)-4*zt(i,j)-af(i,j)
                                f(i,j)=zt(i,j)+rer*rn(i,j)/4
                                zt(i,j)=f(i,j)
                        end do
                end do
                m=1
                n=1
                do i=2,m0-1
                        do j=2,n0-1
                                if(abs(f(i,j)-f(m,n))>eps)then
                                        m=i
                                        n=j
                                end if
                        end do
                end do
                if(abs(f(i,j)-f(m,n))<eps)then
                        exit
                end if
                end do
                return
                end subroutine