copyright (C) 2001 MSC-RPN COMM %%%MC2%%% *subroutine rhs3_p0 ( dimphy, wrk , lminx, lmaxx, lminy, lmaxy ) 1 implicit none * integer dimphy,lminx,lmaxx,lminy,lmaxy real*8 wrk(*) * *AUTHORs C. Girard & M. Desgagne * *OBJECT * * ******************************************************************* * * * * * FOR THE PURE LEAP-FROG SCHEME * * * * * * * * * CALCULATES * * * right-hand-side LINEAR Pm-terms: * * * Pu, Pv, Pw, Pb, Pp * * * FROM * * * time-levels t-Dt and t variables: * * * um, vm, wm, bm, pm * * * u0, v0, w0, b0, p0 * * * STORES * * * results in * * * up=Pu, vp=Pv, wp=Pw, bp=Pb, pp=Pp * * * * * * PERFORMS essentially the following calculations: * * * * * * ~ * * * up = um -dt0 * d(p0)/dX * * * * * * ~ * * * vp = vm -dt0 * d(p0)/dY * * * * * * ~ * * * wp = c03* wm+dt0*b0 -dt0 * d(p0)/dZ * * * * * * __Z * * * bp = bm - gama_star * pm - dt0*c06*n02g*w0 * * * * * * __Z * * * pp = [ c2r_star * pm + dt0*gc02*w0 ] - dt0*div(u0,v0,w0) * * * * * * * * ******************************************************************* * * #include "grd.cdk"
#include "consdyn_8.cdk"
#include "yomdyn.cdk"
#include "yomdyn1.cdk"
#include "dynpar.cdk"
#include "dynmem.cdk"
#include "topo.cdk"
#include "phymem.cdk"
#include "physnml.cdk"
#include "physcom.cdk"
#include "nbcpu.cdk"
#include "partopo.cdk"
#include "dtmdtp.cdk"
#include "alfmem.cdk"
* integer i, j, k, n, id, jd, iff, jf, km1, kp1 real*8 p25, pt5, one, dt0, c01k real s, d pointer (pas ,s(*)), (pad ,d(*)) real*8 q3,div pointer (paq3 , q3(minx:maxx,miny:maxy,*)), $ (padiv, div(minx:maxx,miny:maxy,*)) parameter (p25=0.25d0, pt5=0.5d0, one=1.0d0) *---------------------------------------------------------------------- * paq3 = loc(wrk( 1)) padiv = loc(wrk(dim3d+1)) * dt0 = dtm+dtp * id =1-hx*west jd =1-hy*south iff=ldni+(hx-1)*east jf =ldnj+(hy-1)*north * do k = gnk, 1, -1 if (k.eq.gnk) then do j = jd, jf do i = id, iff q3(i,j,gnk) = 0.0 div(i,j,gnk) = 0.0 end do end do else do j = jd, jf do i = id, iff q3(i,j,k) = sbxy(i,j) * pt5 * $ (-gg1(i+1,j,k) * (uu0(i+1,j,k) -uu0(i+1,j,km1) ) $ - gg1(i ,j,k) * (uu0(i ,j,k) -uu0(i ,j,km1) ) $ - gg2(i,j+1,k) * (vv0(i,j+1,k) -vv0(i,j+1,km1) ) $ - gg2(i,j ,k) * (vv0(i,j ,k) -vv0(i,j ,km1) ) ) div(i,j,k)= sbxy(i,j) * ((uu0(i+1,j,k) -uu0(i,j,k))*odx(1) + $ (vv0(i,j+1,k) -vv0(i,j,k))*ody(j) ) $ + gg0r(i,j,k)* (ww0(i,j,k+1) -ww0(i,j,k) + $ ( q3(i,j,k+1) + q3(i,j,k) ) * pt5 ) end do end do endif end do * do k = gnk, 1, -1 km1 = max( k-1 , 1 ) if (k.lt.gnk) then kp1 = min( k+1 , gnk-1 ) do j = jd, jf do i = id+west, iff uup(i,j,k) = uum(i,j,k) - dt0 * ( $ ( pp0(i,j,k ) - pp0(i-1,j,k) ) * odxu(1) + $ g0ur(i,j,k) * p25 * (-gg1(i,j,k+1) * $ ( pp0(i,j,kp1) + pp0(i-1,j,kp1) $ - pp0(i,j,k ) - pp0(i-1,j,k ) ) - $ gg1(i,j,k ) * $ ( pp0(i,j,k ) + pp0(i-1,j,k ) $ - pp0(i,j,km1) - pp0(i-1,j,km1) )) ) end do end do do j = jd+south, jf do i = id, iff vvp(i,j,k) = vvm(i,j,k) - dt0 * ( $ ( pp0(i,j,k ) - pp0(i,j-1,k) ) * odyv(j) + $ g0vr(i,j,k) * p25 * (-gg2(i,j,k+1) * $ ( pp0(i,j,kp1) + pp0(i,j-1,kp1) $ - pp0(i,j,k ) - pp0(i,j-1,k ) ) - $ gg2(i,j,k ) * $ ( pp0(i,j,k ) + pp0(i,j-1,k ) $ - pp0(i,j,km1) - pp0(i,j-1,km1) )) ) end do end do end if * if (k.eq.gnk) then c01k=one else c01k=c01 endif kp1 = min( k+1 , gnk ) do j = jd, jf do i = id, iff wwp(i,j,k) = c03*wwm(i,j,k) $ - dt0*(pp0(i,j,k)-pp0(i,j,k-1))*g0wr(i,j,k) $ + dt0 *bb0(i,j,k) bbp(i,j,k) = bbm(i,j,k) $ - gama_star*pt5*(ppm(i,j,k)+ppm(i,j,k-1)) $ - dt0*c06*n02g(i,j,k)*ww0(i,j,k) ppp(i,j,k) = c01k * c2r_star * ppm(i,j,k) $ + dt0*pt5*gc02(i,j,k)*(ww0(i,j,k)+ww0(i,j,kp1)) $ - dt0 * div(i,j,k) end do end do end do * *---------------------------------------------------------------------- return end