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