copyright (C) 2001  MSC-RPN COMM  %%%MC2%%%
***s/r cnt3_s Evalue les coefficients czz(k,l)
*

      subroutine cnt3_s (czz,nx,ny,nz,niter)  2
      implicit none
* 
      integer nx,ny,nz,niter
      real*8 czz(1-(niter-1):nx+(niter-1),1-(niter-1):ny+(niter-1),nz,3)
*
*AUTHOR
*
*REVISION
*
*OBJECT
*     Calcule les coefficients de la derivee seconde en z: czz(i,j,k,3) 
*     servant dans l'equation de Helmoltz.
*     le terme d'Helmoltz y est ajoute.
*
*IMPLICIT
#include "dynmem.cdk"
#include "yomdyn.cdk"
#include "yomdyn1.cdk"
#include "dynpar.cdk"
#include "dtmdtp.cdk"
#include "alfmem.cdk"
*
**
      integer i,j,k,dim
      real*8 nu0,one,pt5
      parameter(one=1.0d0,pt5=0.5d0)
*----------------------------------------------------------------------
*
      dim=(niter-1)
*
!$omp do
      do j = 1-dim, ny+dim
         do i = 1-dim, nx+dim
            czz(i,j,1,1)     = 0.0
            czz(i,j,1,2)     = - c01*c2r_star
            czz(i,j,gnk-1,3) = 0.0
         end do
*
         do k=1,gnk-2 
*
         do i = 1-dim, nx+dim
*
         nu0 = dtp*dtp*nu0b(i,j,k+1)
*
         czz(i,j,k+1,1) =            (gg0r (i,j,k+1)+pt5*gc02(i,j,k+1))*
     $                               (g0wr (i,j,k+1)+c05*n02g(i,j,k+1))*
     $                                                     nu0
         czz(i,j,k+1,2) = - c01 * c2r_star  -
     $                               (gg0r (i,j,k+1)+pt5*gc02(i,j,k+1))*
     $                               (g0wr (i,j,k+1)-c05*n02g(i,j,k+1))*
     $                                                     nu0
*
         czz(i,j,k,2)   = czz(i,j,k,2) -
     $                               (gg0r (i,j,k  )-pt5*gc02(i,j,k  ))*
     $                               (g0wr (i,j,k+1)+c05*n02g(i,j,k+1))*
     $                                                     nu0
         czz(i,j,k,3)   =            (gg0r (i,j,k  )-pt5*gc02(i,j,k  ))*
     $                               (g0wr (i,j,k+1)-c05*n02g(i,j,k+1))*
     $                                                     nu0
         end do
         end do
      end do
!$omp enddo
*
      if(flextop) then
!$omp do
         do j = 1-dim, ny+dim
         do i = 1-dim, nx+dim
*
            nu0 = dtp*dtp*nu0b(i,j,gnk)
*
            czz(i,j,gnk  ,1) =                        gc02(i,j,gnk  ) *
     $                          (g0wr (i,j,gnk  )+c05*n02g(i,j,gnk  ))*
     $                                                     nu0
            czz(i,j,gnk  ,2) = - c2r_star   -
     $                                                gc02(i,j,gnk  ) *
     $                          (g0wr (i,j,gnk  )-c05*n02g(i,j,gnk  ))*
     $                                                     nu0
            czz(i,j,gnk  ,3) = 0.0
*
            czz(i,j,gnk-1,2) = czz(i,j,gnk-1,2) -
     $                          (gg0r (i,j,gnk-1)-pt5*gc02(i,j,gnk-1))*
     $                          (g0wr (i,j,gnk  )+c05*n02g(i,j,gnk  ))*
     $                                                     nu0
            czz(i,j,gnk-1,3) =  (gg0r (i,j,gnk-1)-pt5*gc02(i,j,gnk-1))*
     $                          (g0wr (i,j,gnk  )-c05*n02g(i,j,gnk  ))*
     $                                                     nu0
         enddo
         enddo
!$omp enddo
      endif
*
*----------------------------------------------------------------------
      return
      end