copyright (C) 2001  MSC-RPN COMM  %%%MC2%%%
***s/r lhs3_wt

      subroutine lhs3_wt 2
      implicit none
*
*AUTHORs    C. Girard & M. Desgagne
*
*OBJECT
*
*
*    *******************************************************************
*    *                                                                 *
*    *    CALCULATES                                                   *
*    *                time-level t+Dt UPDATED variables                *
*    *                wp, bp                                           *
*    *       FROM                                                      *
*    *                right-hand sides                                 *
*    *                Qw=wp, Qb=bp                                     *
*    *       AND                                                       *
*    *                already updated variables                        *
*    *                up, vp, pp                                       *
*    *                                                                 *
*    *     PERFORMS   essentially the following calculations:          *
*    *                                                                 *
*    *                                    ~                       __Z  *
*    *          wp  =    Qw - dtp*nu0 * [ d(pp)/dZ - 2*c05*n02g * pp  ]*
*    *                                                                 *
*    *                                    __Z                          *
*    *          bp  =    Qb + gama_star * pp  - dtp*c06*n02g * wp      *
*    *                                                                 *
*    *                                                                 *
*    *                                       ~                         *
*    *      SEE  - rhs3_pm for definition of d/dZ                      *
*    *                                                                 *
*    *           - cxxpar for definitions of parameters c05, c06       *
*    *                               and gama_star                     *
*    *                                                                 *
*    *           - nacbar for calculation of  nu0b, gc02               *
*    *                                                                 *
*    *******************************************************************
*
*
#include "nbcpu.cdk"
#include "dynmem.cdk"
#include "consdyn_8.cdk"
#include "yomdyn.cdk"
#include "yomdyn1.cdk"
#include "dynpar.cdk"
#include "topo.cdk"
#include "dtmdtp.cdk"
#include "alfmem.cdk"
*
      integer i,j,k
      real*8 pt5, q3, wwp1, wwpgnk, nu0
      parameter (pt5 = 0.5d0)
*
*---------------------------------------------------------------
*
!$omp do
      do k = 1, gnk
*
         if (k.eq.1) then
            do j=1,ldnj-north
            do i=1,ldni-east
               nu0  =  dtp * nu0b(i,j,k)
               wwp1       = wwp(i,j,k)
               q3         = sbxy(i,j) * pt5 *
     $              (- gg1(i+1,j,k)*uup(i+1,j,k) - gg1(i,j,k)*uup(i,j,k)
     $               - gg2(i,j+1,k)*vvp(i,j+1,k) - gg2(i,j,k)*vvp(i,j,k) )
               wwp(i,j,k) = dhdt(i,j,k) - q3
               ppp(i,j,0) = ( wwp(i,j,k) - wwp1 + ( nu0*
     $                      (g0wr(i,j,k)-c05*n02g(i,j,k)))*ppp(i,j,k) )
     $                    / (nu0*(g0wr(i,j,k)+c05*n02g(i,j,k)))
               bbp(i,j,k) = bbp(i,j,k)
     $                   + gama_star*pt5*(ppp(i,j,k)+ppp(i,j,k-1))
     $                   - dtp*c06*n02g(i,j,k)*wwp(i,j,k)
            end do
            end do
*
         else if ((k.eq.gnk).and.(.not.flextop)) then
            do j=1,ldnj-north
            do i=1,ldni-east
               nu0  =  dtp * nu0b(i,j,k)
               wwpgnk = wwp(i,j,k)
               wwp(i,j,k) = 0.
               ppp(i,j,k) = ( - wwp(i,j,k) + wwpgnk + ( nu0*
     $                     (g0wr(i,j,k)+c05*n02g(i,j,k)))*ppp(i,j,k-1))
     $                    / (nu0*(g0wr(i,j,k)-c05*n02g(i,j,k)))
               bbp(i,j,k) = bbp(i,j,k)
     $                   + gama_star*pt5*(ppp(i,j,k)+ppp(i,j,k-1))
     $                   - dtp*c06*n02g(i,j,k)*wwp(i,j,k)
            end do
            end do
*
         else
            do j=1,ldnj-north
            do i=1,ldni-east
               nu0  =  dtp * nu0b(i,j,k)
               wwp(i,j,k) = wwp(i,j,k) - nu0 *
     $                     ( (ppp(i,j,k)-ppp(i,j,k-1))*g0wr(i,j,k)
     $                     - (ppp(i,j,k)+ppp(i,j,k-1))*c05*n02g(i,j,k))
               bbp(i,j,k) = bbp(i,j,k)
     $                   + gama_star*pt5*(ppp(i,j,k)+ppp(i,j,k-1))
     $                   - dtp*c06*n02g(i,j,k)*wwp(i,j,k)
            end do
            end do
*
         endif
*
      end do
!$omp enddo
*
*---------------------------------------------------------------
      return
      end