source: LMDZ5/trunk/libf/phylmd/iniphysiq.F @ 1941

!
! $Id: iniphysiq.F 1907 2013-11-26 13:10:46Z jghattas $
!
c
c
      SUBROUTINE iniphysiq(ngrid,nlayer,
     $           punjours,
     $           pdayref,ptimestep,
     $           plat,plon,parea,pcu,pcv,
     $           prad,pg,pr,pcpp,iflag_phys)
      USE dimphy, only : klev
      USE mod_grid_phy_lmdz, only : klon_glo
      USE mod_phys_lmdz_para, only : klon_omp,klon_omp_begin,
     &                               klon_omp_end,klon_mpi_begin
      USE comgeomphy, only : airephy,cuphy,cvphy,rlond,rlatd
      USE phyaqua_mod, only: iniaqua
      IMPLICIT NONE
c
c=======================================================================
c
c   Initialisation of the physical constants and some positional and 
c   geometrical arrays for the physics
c
c
c    ngrid                 Size of the horizontal grid.
c                          All internal loops are performed on that grid.
c    nlayer                Number of vertical layers.
c    pdayref               Day of reference for the simulation
c
c=======================================================================
 
cym#include "dimensions.h"
cym#include "dimphy.h"
cym#include "comgeomphy.h"
#include "YOMCST.h"
#include "iniprint.h"

      REAL,INTENT(IN) :: prad ! radius of the planet (m)
      REAL,INTENT(IN) :: pg ! gravitational acceleration (m/s2)
      REAL,INTENT(IN) :: pr ! ! reduced gas constant R/mu
      REAL,INTENT(IN) :: pcpp ! specific heat Cp
      REAL,INTENT(IN) :: punjours ! length (in s) of a standard day
      INTEGER,INTENT(IN) :: ngrid ! number of horizontal grid points in the physics
      INTEGER,INTENT(IN) :: nlayer ! number of atmospheric layers
      REAL,INTENT(IN) :: plat(ngrid) ! latitudes of the physics grid
      REAL,INTENT(IN) :: plon(ngrid) ! longitudes of the physics grid
      REAL,INTENT(IN) :: parea(klon_glo) ! area (m2)
      REAL,INTENT(IN) :: pcu(klon_glo) ! cu coeff. (u_covariant = cu * u)
      REAL,INTENT(IN) :: pcv(klon_glo) ! cv coeff. (v_covariant = cv * v)
      INTEGER,INTENT(IN) :: pdayref ! reference day of for the simulation
      REAL,INTENT(IN) :: ptimestep !physics time step (s)
      INTEGER,INTENT(IN) :: iflag_phys ! type of physics to be called

      INTEGER :: ibegin,iend,offset 
      CHARACTER (LEN=20) :: modname='iniphysiq'
      CHARACTER (LEN=80) :: abort_message
 
      IF (nlayer.NE.klev) THEN
         write(lunout,*) 'STOP in ',trim(modname)
         write(lunout,*) 'Problem with dimensions :'
         write(lunout,*) 'nlayer     = ',nlayer
         write(lunout,*) 'klev   = ',klev
         abort_message = ''
         CALL abort_gcm (modname,abort_message,1)
      ENDIF

      IF (ngrid.NE.klon_glo) THEN
         write(lunout,*) 'STOP in ',trim(modname)
         write(lunout,*) 'Problem with dimensions :'
         write(lunout,*) 'ngrid     = ',ngrid
         write(lunout,*) 'klon   = ',klon_glo
         abort_message = ''
         CALL abort_gcm (modname,abort_message,1)
      ENDIF

!$OMP PARALLEL PRIVATE(ibegin,iend) 
!$OMP+         SHARED(parea,pcu,pcv,plon,plat)
      
      offset=klon_mpi_begin-1
      airephy(1:klon_omp)=parea(offset+klon_omp_begin:
     &                          offset+klon_omp_end)
      cuphy(1:klon_omp)=pcu(offset+klon_omp_begin:offset+klon_omp_end)
      cvphy(1:klon_omp)=pcv(offset+klon_omp_begin:offset+klon_omp_end)
      rlond(1:klon_omp)=plon(offset+klon_omp_begin:offset+klon_omp_end)
      rlatd(1:klon_omp)=plat(offset+klon_omp_begin:offset+klon_omp_end)

      ! suphel => initialize some physical constants (orbital parameters,
      !           geoid, gravity, thermodynamical constants, etc.) in the
      !           physics
      call suphel
      
!$OMP END PARALLEL

      ! check that physical constants set in 'suphel' are coherent
      ! with values set in the dynamics:
      if (RDAY.ne.punjours) then
        write(lunout,*) "iniphysiq: length of day discrepancy!!!"
        write(lunout,*) "  in the dynamics punjours=",punjours
        write(lunout,*) "   but in the physics RDAY=",RDAY
        if (abs(RDAY-punjours).gt.0.01) then
          ! stop here if the relative difference is more than 1%
          abort_message = 'length of day discrepancy'
          CALL abort_gcm (modname,abort_message,1)
        endif
      endif
      if (RG.ne.pg) then
        write(lunout,*) "iniphysiq: gravity discrepancy !!!"
        write(lunout,*) "     in the dynamics pg=",pg
        write(lunout,*) "  but in the physics RG=",RG
        if (abs(RG-pg).gt.0.01) then
          ! stop here if the relative difference is more than 1%
          abort_message = 'gravity discrepancy'
          CALL abort_gcm (modname,abort_message,1)
        endif
      endif
      if (RA.ne.prad) then
        write(lunout,*) "iniphysiq: planet radius discrepancy !!!"
        write(lunout,*) "   in the dynamics prad=",prad
        write(lunout,*) "  but in the physics RA=",RA
        if (abs(RA-prad).gt.0.01) then
          ! stop here if the relative difference is more than 1%
          abort_message = 'planet radius discrepancy'
          CALL abort_gcm (modname,abort_message,1)
        endif
      endif
      if (RD.ne.pr) then
        write(lunout,*)"iniphysiq: reduced gas constant discrepancy !!!"
        write(lunout,*)"     in the dynamics pr=",pr
        write(lunout,*)"  but in the physics RD=",RD
        if (abs(RD-pr).gt.0.01) then
          ! stop here if the relative difference is more than 1%
          abort_message = 'reduced gas constant discrepancy'
          CALL abort_gcm (modname,abort_message,1)
        endif
      endif
      if (RCPD.ne.pcpp) then
        write(lunout,*)"iniphysiq: specific heat discrepancy !!!"
        write(lunout,*)"     in the dynamics pcpp=",pcpp
        write(lunout,*)"  but in the physics RCPD=",RCPD
        if (abs(RCPD-pcpp).gt.0.01) then
          ! stop here if the relative difference is more than 1%
          abort_message = 'specific heat discrepancy'
          CALL abort_gcm (modname,abort_message,1)
        endif
      endif

! Additional initializations for aquaplanets
!$OMP PARALLEL
      if (iflag_phys>=100) then
        call iniaqua(klon_omp,rlatd,rlond,iflag_phys)
      endif
!$OMP END PARALLEL

!      RETURN
!9999  CONTINUE
!      abort_message ='Cette version demande les fichier rnatur.dat
!     & et surf.def'
!      CALL abort_gcm (modname,abort_message,1)

      END