source: LMDZ6/branches/Amaury_dev/libf/dyn3d_common/divergst.f90 @ 5209

! $Header$

SUBROUTINE divergst(klevel, x, y, div)
  USE lmdz_ssum_scopy, ONLY: ssum
  USE lmdz_comgeom

  USE lmdz_dimensions, ONLY: iim, jjm, llm, ndm
  USE lmdz_paramet
  IMPLICIT NONE

  ! P. Le Van

  !  ******************************************************************
  !  ... calcule la divergence a tous les niveaux d'1 vecteur de compos. x et y...
  !       x et y  etant des composantes contravariantes   ...
  !  ****************************************************************
  !  x  et  y  sont des arguments  d'entree pour le s-prog
  !    div      est  un argument  de sortie pour le s-prog


  !   -------------------------------------------------------------------
  !




  INTEGER :: klevel
  REAL :: x(ip1jmp1, klevel), y(ip1jm, klevel), div(ip1jmp1, klevel)
  INTEGER :: ij, l, i
  REAL :: aiy1(iip1), aiy2(iip1)
  REAL :: sumypn, sumyps


  DO l = 1, klevel

    DO ij = iip2, ip1jm - 1
      div(ij + 1, l) = x(ij + 1, l) - x(ij, l) + y(ij - iim, l) - y(ij + 1, l)
    END DO

    ! ....  correction pour  div( 1,j,l)  ......
    ! ....   div(1,j,l)= div(iip1,j,l) ....

    !DIR$ IVDEP
    DO ij = iip2, ip1jm, iip1
      div(ij, l) = div(ij + iim, l)
    END DO

    ! ....  calcul  aux poles  .....


    DO i = 1, iim
      aiy1(i) = y(i, l)
      aiy2(i) = y(i + ip1jmi1, l)
    END DO
    sumypn = SSUM (iim, aiy1, 1)
    sumyps = SSUM (iim, aiy2, 1)
    DO i = 1, iip1
      div(i, l) = - sumypn / iim
      div(i + ip1jm, l) = sumyps / iim
    END DO

  END DO
  RETURN
END SUBROUTINE divergst