! ! $Id: iniphysiq.F 1403 2010-07-01 09:02:53Z fairhead $ ! MODULE iniphysiq_mod CONTAINS SUBROUTINE iniphysiq(iim,jjm,nlayer,punjours, pdayref,ptimestep, & rlatu,rlatv,rlonu,rlonv,aire,cu,cv, & prad,pg,pr,pcpp,iflag_phys) USE dimphy, ONLY: klev ! number of atmospheric levels USE mod_grid_phy_lmdz, ONLY: klon_glo ! number of atmospheric columns ! (on full grid) USE mod_phys_lmdz_para, ONLY: klon_omp, & ! number of columns (on local omp grid) klon_omp_begin, & ! start index of local omp subgrid klon_omp_end, & ! end index of local omp subgrid klon_mpi_begin ! start indes of columns (on local mpi grid) USE comgeomphy, ONLY: initcomgeomphy, & airephy, & ! physics grid area (m2) cuphy, & ! cu coeff. (u_covariant = cu * u) cvphy, & ! cv coeff. (v_covariant = cv * v) rlond, & ! longitudes rlatd ! latitudes USE infotrac, ONLY: nqtot, type_trac USE infotrac_phy, ONLY: init_infotrac_phy ! USE comcstphy, ONLY: rradius, & ! planet radius (m) ! rr, & ! recuced gas constant: R/molar mass of atm ! rg, & ! gravity ! rcpp ! specific heat of the atmosphere USE inifis_mod, ONLY: inifis USE phyaqua_mod, ONLY: iniaqua USE regular_lonlat_mod, ONLY : init_regular_lonlat, & east, west, north, south, & north_east, north_west, & south_west, south_east USE nrtype, ONLY: pi IMPLICIT NONE ! !======================================================================= ! Initialisation of the physical constants and some positional and ! geometrical arrays for the physics !======================================================================= 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) :: nlayer ! number of atmospheric layers INTEGER, INTENT (IN) :: iim ! number of atmospheric coulumns along longitudes INTEGER, INTENT (IN) :: jjm ! number of atompsheric columns along latitudes REAL, INTENT (IN) :: rlatu(jjm+1) ! latitudes of the physics grid REAL, INTENT (IN) :: rlatv(jjm) ! latitude boundaries of the physics grid REAL, INTENT (IN) :: rlonv(iim+1) ! longitudes of the physics grid REAL, INTENT (IN) :: rlonu(iim+1) ! longitude boundaries of the physics grid REAL, INTENT (IN) :: aire(iim+1,jjm+1) ! area of the dynamics grid (m2) REAL, INTENT (IN) :: cu((iim+1)*(jjm+1)) ! cu coeff. (u_covariant = cu * u) REAL, INTENT (IN) :: cv((iim+1)*jjm) ! 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 INTEGER :: i,j CHARACTER (LEN=20) :: modname='iniphysiq' CHARACTER (LEN=80) :: abort_message REAL :: total_area_phy, total_area_dyn ! boundaries, on global grid REAL,ALLOCATABLE :: boundslon_reg(:,:) REAL,ALLOCATABLE :: boundslat_reg(:,:) ! global array, on full physics grid: REAL,ALLOCATABLE :: latfi(:) REAL,ALLOCATABLE :: lonfi(:) REAL,ALLOCATABLE :: cufi(:) REAL,ALLOCATABLE :: cvfi(:) REAL,ALLOCATABLE :: airefi(:) 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 !call init_phys_lmdz(iim,jjm+1,llm,1,(/(jjm-1)*iim+2/)) ! init regular global longitude-latitude grid points and boundaries ALLOCATE(boundslon_reg(iim,2)) ALLOCATE(boundslat_reg(jjm+1,2)) DO i=1,iim boundslon_reg(i,east)=rlonu(i) boundslon_reg(i,west)=rlonu(i+1) ENDDO boundslat_reg(1,north)= PI/2 boundslat_reg(1,south)= rlatv(1) DO j=2,jjm boundslat_reg(j,north)=rlatv(j-1) boundslat_reg(j,south)=rlatv(j) ENDDO boundslat_reg(jjm+1,north)= rlatv(jjm) boundslat_reg(jjm+1,south)= -PI/2 ! Write values in module regular_lonlat_mod CALL init_regular_lonlat(iim,jjm+1, rlonv(1:iim), rlatu, & boundslon_reg, boundslat_reg) ! Generate global arrays on full physics grid ALLOCATE(latfi(klon_glo),lonfi(klon_glo),cufi(klon_glo),cvfi(klon_glo)) ALLOCATE(airefi(klon_glo)) ! North pole latfi(1)=rlatu(1) lonfi(1)=0. cufi(1) = cu(1) cvfi(1) = cv(1) DO j=2,jjm DO i=1,iim latfi((j-2)*iim+1+i)= rlatu(j) lonfi((j-2)*iim+1+i)= rlonv(i) cufi((j-2)*iim+1+i) = cu((j-1)*iim+1+i) cvfi((j-2)*iim+1+i) = cv((j-1)*iim+1+i) ENDDO ENDDO ! South pole latfi(klon_glo)= rlatu(jjm+1) lonfi(klon_glo)= 0. cufi(klon_glo) = cu((iim+1)*jjm+1) cvfi(klon_glo) = cv((iim+1)*jjm-iim) ! build airefi(), mesh area on physics grid CALL gr_dyn_fi(1,iim+1,jjm+1,klon_glo,aire,airefi) ! Poles are single points on physics grid airefi(1)=sum(aire(1:iim,1)) airefi(klon_glo)=sum(aire(1:iim,jjm+1)) ! Sanity check: do total planet area match between physics and dynamics? total_area_dyn=sum(aire(1:iim,1:jjm+1)) total_area_phy=sum(airefi(1:klon_glo)) IF (total_area_dyn/=total_area_phy) THEN WRITE (lunout, *) 'iniphysiq: planet total surface discrepancy !!!' WRITE (lunout, *) ' in the dynamics total_area_dyn=', total_area_dyn WRITE (lunout, *) ' but in the physics total_area_phy=', total_area_phy IF (abs(total_area_dyn-total_area_phy)>0.00001*total_area_dyn) THEN ! stop here if the relative difference is more than 0.001% abort_message = 'planet total surface discrepancy' CALL abort_gcm(modname, abort_message, 1) ENDIF ENDIF !$OMP PARALLEL ! Initialize physical constants in physics: CALL inifis(prad,pg,pr,pcpp) ! Initialize tracer names, numbers, etc. for physics CALL init_infotrac_phy(nqtot,type_trac) ! Now generate local lon/lat/cu/cv/area arrays CALL initcomgeomphy offset = klon_mpi_begin - 1 airephy(1:klon_omp) = airefi(offset+klon_omp_begin:offset+klon_omp_end) cuphy(1:klon_omp) = cufi(offset+klon_omp_begin:offset+klon_omp_end) cvphy(1:klon_omp) = cvfi(offset+klon_omp_begin:offset+klon_omp_end) rlond(1:klon_omp) = lonfi(offset+klon_omp_begin:offset+klon_omp_end) rlatd(1:klon_omp) = latfi(offset+klon_omp_begin:offset+klon_omp_end) ! Additional initializations for aquaplanets IF (iflag_phys>=100) THEN CALL iniaqua(klon_omp,rlatd,rlond,iflag_phys) ENDIF !$OMP END PARALLEL END SUBROUTINE iniphysiq END MODULE iniphysiq_mod