copyright (C) 2001  MSC-RPN COMM  %%%MC2%%%
*

      subroutine step8 (dimphy,trname,nt) 2,21
      implicit none
*
      integer dimphy,nt
      character*8 trname(nt)
*
#include "rec.cdk"
#include "consdyn_8.cdk"
#include "yomdyn.cdk"
#include "yomdyn1.cdk"
#include "dynmem.cdk"
#include "partopo.cdk"
#include "nbcpu.cdk"
#include "bcsdim.cdk"
#include "bcsmem.cdk"
#include "nesting.cdk"
#include "wrnmem.cdk"
#include "alfmem.cdk"
#include "dtmdtp.cdk"
#include "gcrk.cdk"
*
      integer lastep, stepno, err, nsol
      real, dimension (:), allocatable :: w1
      real*8, dimension (:), allocatable :: xyzd,wrk_sol1,wrk_sol2
      real wrn1,wrn2,dtf
      pointer (pawrn1, wrn1(*)), (pawrn2, wrn2(*))
*
      real*8 tss,pt5,one
      logical semi_imp,vert_imp,alfa_imp,damp_div,whenby1
      data whenby1 /.false./
      parameter (pt5 = 0.5d0, one = 1.0d0)
      save lastep,whenby1
*----------------------------------------------------------------------
      semi_imp=.true.
      vert_imp=.false.
      alfa_imp=.true.
      alfa_imp=alfa_imp.and.(semi_imp.or.vert_imp)
      damp_div=.false.
*
      tss = dble(grdt)
      if ( gnstepno .eq. 1 ) tss = pt5 * dble(grdt)
*
      dtm = (one-dble(grepsi)) * tss
      dtp = (one+dble(grepsi)) * tss
*
*     Switching T0 into T-   and T+ into T0 (flip pointers)
*     and re-apply lower boundary conditions W(i,j,1)=0 
*
*     Sink point
*
      call rpn_comm_xch_halo (pp0,minx,maxx,miny,maxy,ldni,ldnj,
     $              ndynvar*gnk+1,hx,hy,period_x,period_y,ldni,0)
      call flipnt ()
* 
      nsol = (5*dim_s*nz + 2*(maxx-minx+3)*(maxy-miny+3)*nz +
     $                     ny+2*(niterj-1))
      nsol = nsol + 7*dim3d + 3*nksp*dim3d + dim_s + dim_czz
*
      allocate (xyzd(10*dim3d),w1(max(10,ndynvar)*dim3d),wrk_sol1(nsol),
     $          wrk_sol2(dim3d*(im+1)+nthreads))
*
      call hpalloc (pawrn1, 5*dim3d, err,1)
      paalfa=  loc( wrn1(0*dim3d+1) )
      paap1r=  loc( wrn1(1*dim3d+1) )
      panu0b=  loc( wrn1(2*dim3d+1) )
      pan02g=  loc( wrn1(3*dim3d+1) )
      pagc02=  loc( wrn1(4*dim3d+1) )
*
      call hpalloc (pawrn2, 5*dim3d, err,1)
      pauur =  loc( wrn2(        1) )
      pavvr =  loc( wrn2(1*dim3d+1) )
      pawwr =  loc( wrn2(2*dim3d+1) )
      pabbr =  loc( wrn2(3*dim3d+1) )
      pappr =  loc( wrn2(4*dim3d+1) )
*
!$omp parallel
*
*     Update Metric Parameters if necessary
*
      call metric (w1, xyzd)
      call nacbar (alfa_imp)
*
*     Prepare Variables for Advection
*
!$omp single
      if (myproc.eq.0) print *,'P,R,Q-TERMS'
!$omp end single
*
      call uvadv3
*
      if (.not.nosolv) then
         if(semi_imp)then
            call rhs3_pm ( dimphy,w1,minx,maxx,miny,maxy )
         else if(vert_imp) then
            call rhs3_p0m( dimphy,w1,minx,maxx,miny,maxy )
         else
            call rhs3_p0 ( dimphy,w1,minx,maxx,miny,maxy )
         endif
         if(damp_div) then
            if(myproc.eq.0) print *,'damping 3d-divergence'
            call difdiv
         endif
         call rhs3_r0 ( xyzd )
      endif
*
      if (semi_lag) then
!$omp single
         if (myproc.eq.0) print *,'Semi-Lagrangian Advection'
!$omp end single
         call traject ( xyzd, w1, dtm, dtp )
         call fwrd3d  ( xyzd, w1, dtm, dtp )
      else
!$omp single
         if (myproc.eq.0) print *,'Eulerian Advection'
         call euler  ( w1, dtm, dtp )
!$omp end single
      endif
*
      call hpdeallc (pawrn2 ,err, 1)
*
      stepno = gnstepno
      if ( grdt .lt. 0.) then
         stepno = lastep - gnstepno
      else
         lastep = gnstepno
      endif
*
*     Update nesting data for the current validity time
*     (linear interpolation in time between unestt (t-dt)
*     and unestta).
*
      if (.not.ctebcs) call nest_intt (stepno, trname, nt, dtf, whenby1)
*
*     Specify Normal Wind Components on Model Boundaries
*
      call uvbdy (bcs_uu,bcs_uu(bcs_in),bcs_uu(bcs_iw),bcs_uu(bcs_ie),
     $            bcs_vv,bcs_vv(bcs_in),minxs,maxxs,minys,maxys,
     $                                  minxw,maxxw,minyw,maxyw)
*
      if (.not.nosolv) then
         if (semi_imp) then
*
*           Semi-Implicit Scheme
*
*           Pressure Variable: Solve the Elliptic Equation
*
!$omp single
            if (myproc.eq.0) print *,'Semi-Implicit Time Integration'
!$omp end single
            call solver2_p (wrk_sol1,wrk_sol2)
*
*           Other Variables: Back Substitution 
*
            call lhs3_uv
            call lhs3_wt
*
         else if(vert_imp) then
!$omp single
            if (myproc.eq.0) print *,'Semi-Implicit in vertical only'
            call sis_1d
            call rpn_comm_xch_halo(ppp(minx,miny,1),minx,maxx,miny,maxy,
     $                     ldni,ldnj,gnk,hx,hy,period_x,period_y,ldni,0)
            call lhs3_wt
!$omp end single
         else
!$omp single
            if (myproc.eq.0) print *,'Pure Leapfrog Scheme'
            call lhs3_qt
!$omp end single
         endif
*
      endif
*
!$omp end parallel
*
      deallocate    ( xyzd,w1,wrk_sol1,wrk_sol2 )
      call hpdeallc ( pawrn1 ,err, 1 )
*
*     Use Nesting Data or Boundary/Symmetry Conditions 
*                  to Update EXTERIOR Halo
*     Perform Gravity Wave Absorption on Global Boundaries and Model Top
*
      call glb_bound (wall,slab)
*
*----------------------------------------------------------------------
      return
      end