! $Id: iniphysiq.F90 2320 2015-07-01 13:57:32Z fairhead $ SUBROUTINE iniphysiq(ii,jj,nlayer,punjours, pdayref,ptimestep, & rlatu,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 vertical_layers_mod, ONLY : init_vertical_layers USE infotrac, ONLY: nqtot,nqo,nbtr,tname,ttext,type_trac,& niadv,conv_flg,pbl_flg,solsym,& nqfils,nqdesc,nqdesc_tot,iqfils,iqpere,& ok_isotopes,ok_iso_verif,ok_isotrac,& ok_init_iso,niso_possibles,tnat,& alpha_ideal,use_iso,iqiso,iso_num,& iso_indnum,zone_num,phase_num,& indnum_fn_num,index_trac,& niso,ntraceurs_zone,ntraciso 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 inifis_mod, ONLY: inifis USE infotrac_phy, ONLY: init_infotrac_phy USE phyaqua_mod, ONLY: iniaqua IMPLICIT NONE ! ======================================================================= ! Initialisation of the physical constants and some positional and ! geometrical arrays for the physics ! ======================================================================= include "dimensions.h" include "comvert.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) :: nlayer ! number of atmospheric layers INTEGER, INTENT (IN) :: ii ! number of atmospheric columns along longitudes INTEGER, INTENT (IN) :: jj ! number of atompsheric columns along latitudes REAL, INTENT (IN) :: rlatu(jj+1) ! latitudes of the physics grid REAL, INTENT (IN) :: rlonv(ii+1) ! longitudes of the physics grid REAL, INTENT (IN) :: aire(ii+1,jj+1) ! area of the dynamics grid (m2) REAL, INTENT (IN) :: cu((ii+1)*(jj+1)) ! cu coeff. (u_covariant = cu * u) REAL, INTENT (IN) :: cv((ii+1)*jj) ! 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 ! global array, on full physics grid: REAL,ALLOCATABLE :: latfi(:) REAL,ALLOCATABLE :: lonfi(:) REAL,ALLOCATABLE :: cufi(:) REAL,ALLOCATABLE :: cvfi(:) REAL,ALLOCATABLE :: airefi(:) IF (nlayer/=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) END IF !call init_phys_lmdz(ii,jj+1,llm,1,(/(jj-1)*ii+2/)) ! Generate global arrays on full physics grid ALLOCATE(latfi(klon_glo),lonfi(klon_glo),cufi(klon_glo),cvfi(klon_glo)) ALLOCATE(airefi(klon_glo)) IF (klon_glo>1) THEN ! general case ! North pole latfi(1)=rlatu(1) lonfi(1)=0. cufi(1) = cu(1) cvfi(1) = cv(1) DO j=2,jj DO i=1,ii latfi((j-2)*ii+1+i)= rlatu(j) lonfi((j-2)*ii+1+i)= rlonv(i) cufi((j-2)*ii+1+i) = cu((j-1)*(ii+1)+i) cvfi((j-2)*ii+1+i) = cv((j-1)*(ii+1)+i) ENDDO ENDDO ! South pole latfi(klon_glo)= rlatu(jj+1) lonfi(klon_glo)= 0. cufi(klon_glo) = cu((ii+1)*jj+1) cvfi(klon_glo) = cv((ii+1)*jj-ii) ! build airefi(), mesh area on physics grid CALL gr_dyn_fi(1,ii+1,jj+1,klon_glo,aire,airefi) ! Poles are single points on physics grid airefi(1)=sum(aire(1:ii,1)) airefi(klon_glo)=sum(aire(1:ii,jj+1)) ! Sanity check: do total planet area match between physics and dynamics? total_area_dyn=sum(aire(1:ii,1:jj+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 ELSE ! klon_glo==1, running the 1D model ! just copy over input values latfi(1)=rlatu(1) lonfi(1)=rlonv(1) cufi(1)=cu(1) cvfi(1)=cv(1) airefi(1)=aire(1,1) ENDIF ! of IF (klon_glo>1) !$OMP PARALLEL ! Initialize physical constants in physics: CALL inifis(punjours,prad,pg,pr,pcpp) ! copy over preff , ap(), bp(), etc CALL init_vertical_layers(nlayer,preff,scaleheight, & ap,bp,presnivs,pseudoalt) ! Initialize tracer names, numbers, etc. for physics CALL init_infotrac_phy(nqtot,nqo,nbtr,tname,ttext,type_trac,& niadv,conv_flg,pbl_flg,solsym,& nqfils,nqdesc,nqdesc_tot,iqfils,iqpere,& ok_isotopes,ok_iso_verif,ok_isotrac,& ok_init_iso,niso_possibles,tnat,& alpha_ideal,use_iso,iqiso,iso_num,& iso_indnum,zone_num,phase_num,& indnum_fn_num,index_trac,& niso,ntraceurs_zone,ntraciso) ! 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) END IF !$OMP END PARALLEL END SUBROUTINE iniphysiq