Index: trunk/LMDZ.MARS/libf/dynphy_lonlat/phymars/iniphysiq_mod.F90 =================================================================== --- trunk/LMDZ.MARS/libf/dynphy_lonlat/phymars/iniphysiq_mod.F90 (revision 1549) +++ trunk/LMDZ.MARS/libf/dynphy_lonlat/phymars/iniphysiq_mod.F90 (revision 1563) @@ -6,31 +6,18 @@ nbp, communicator, & punjours, pdayref,ptimestep, & - rlatu,rlatv,rlonu,rlonv,aire,cu,cv, & + rlatudyn,rlatvdyn,rlonudyn,rlonvdyn, & + airedyn,cudyn,cvdyn, & prad,pg,pr,pcpp,iflag_phys) -use dimphy, only : init_dimphy -use mod_grid_phy_lmdz, only : klon_glo, & ! number of atmospheric columns (on full grid) - regular_lonlat ! regular longitude-latitude grid type -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 geometry_mod, only: init_geometry -!use comgeomphy, only : initcomgeomphy, & -! cell_area, & ! physics grid area (m2) -! dx, & ! cu coeff. (u_covariant = cu * u) -! dy, & ! cv coeff. (v_covariant = cv * v) -! longitude, & ! longitudes (rad) -! latitude ! latitudes (rad) use infotrac, only : nqtot ! number of advected tracers use comgeomfi_h, only: ini_fillgeom use temps_mod, only: day_ini, hour_ini use phys_state_var_init_mod, only: phys_state_var_init -use physics_distribution_mod, only: init_physics_distribution -use regular_lonlat_mod, only: init_regular_lonlat, & - east, west, north, south, & - north_east, north_west, & - south_west, south_east -use mod_interface_dyn_phys, only: init_interface_dyn_phys + +use geometry_mod, only: cell_area, & ! physics grid area (m2) + longitude, & ! longitudes (rad) + latitude ! latitudes (rad) +! necessary to get klon_omp +USE mod_phys_lmdz_para, ONLY: klon_omp ! number of columns (on local omp grid) implicit none @@ -43,5 +30,4 @@ 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 (full grid) integer,intent(in) :: nlayer ! number of atmospheric layers integer,intent(in) :: ii ! number of atmospheric coulumns along longitudes @@ -49,177 +35,38 @@ integer,intent(in) :: nbp ! number of physics columns for this MPI process integer,intent(in) :: communicator ! MPI communicator -real,intent(in) :: rlatu(jj+1) ! latitudes of the physics grid -real,intent(in) :: rlatv(jj) ! latitude boundaries of the physics grid -real,intent(in) :: rlonv(ii+1) ! longitudes of the physics grid -real,intent(in) :: rlonu(ii+1) ! longitude boundaries 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) +real,intent(in) :: rlatudyn(jj+1) ! latitudes of the physics grid +real,intent(in) :: rlatvdyn(jj) ! latitude boundaries of the physics grid +real,intent(in) :: rlonvdyn(ii+1) ! longitudes of the physics grid +real,intent(in) :: rlonudyn(ii+1) ! longitude boundaries of the physics grid +real,intent(in) :: airedyn(ii+1,jj+1) ! area of the dynamics grid (m2) +real,intent(in) :: cudyn((ii+1)*(jj+1)) ! cu coeff. (u_covariant = cu * u) +real,intent(in) :: cvdyn((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,k -character(len=20) :: modname='iniphysiq' -character(len=80) :: abort_message -real :: total_area_phy, total_area_dyn -real :: pi + ! the common part for all planetary physics + !------------------------------------------ + ! --> initialize physics distribution, global fields and geometry + CALL inigeom(ii,jj,nlayer, & + nbp, communicator, & + rlatudyn,rlatvdyn, & + rlonudyn,rlonvdyn, & + airedyn,cudyn,cvdyn) -! boundaries, on global grid -real,allocatable :: boundslon_reg(:,:) -real,allocatable :: boundslat_reg(:,:) + ! the distinct part for all planetary physics + !------------------------------------------ -! global array, on full physics grid: -real,allocatable :: latfi_glo(:) -real,allocatable :: lonfi_glo(:) -real,allocatable :: cufi_glo(:) -real,allocatable :: cvfi_glo(:) -real,allocatable :: airefi_glo(:) -real,allocatable :: boundslonfi_glo(:,:) -real,allocatable :: boundslatfi_glo(:,:) - -! local arrays, on given MPI/OpenMP domain: -real,allocatable,save :: latfi(:) -real,allocatable,save :: lonfi(:) -real,allocatable,save :: cufi(:) -real,allocatable,save :: cvfi(:) -real,allocatable,save :: airefi(:) -real,allocatable,save :: boundslonfi(:,:) -real,allocatable,save :: boundslatfi(:,:) -!$OMP THREADPRIVATE (latfi,lonfi,cufi,cvfi,airefi,boundslonfi,boundslatfi) - - -pi=2.*asin(1.0) - -! Initialize Physics distibution and parameters and interface with dynamics -CALL init_physics_distribution(regular_lonlat,4, & - nbp,ii,jj+1,nlayer,communicator) -CALL init_interface_dyn_phys - -! init regular global longitude-latitude grid points and boundaries -ALLOCATE(boundslon_reg(ii,2)) -ALLOCATE(boundslat_reg(jj+1,2)) +!$OMP PARALLEL -DO i=1,ii - 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,jj - boundslat_reg(j,north)=rlatv(j-1) - boundslat_reg(j,south)=rlatv(j) -ENDDO -boundslat_reg(jj+1,north)= rlatv(jj) -boundslat_reg(jj+1,south)= -PI/2 - -! Write values in module regular_lonlat_mod -CALL init_regular_lonlat(ii,jj+1, rlonv(1:ii), rlatu, & - boundslon_reg, boundslat_reg) - -! Generate global arrays on full physics grid -allocate(latfi_glo(klon_glo),lonfi_glo(klon_glo)) -allocate(cufi_glo(klon_glo),cvfi_glo(klon_glo)) -allocate(airefi_glo(klon_glo)) -allocate(boundslonfi_glo(klon_glo,4)) -allocate(boundslatfi_glo(klon_glo,4)) - -! North pole -latfi_glo(1)=rlatu(1) -lonfi_glo(1)=0. -cufi_glo(1) = cu(1) -cvfi_glo(1) = cv(1) -boundslonfi_glo(1,north_east)=0 -boundslatfi_glo(1,north_east)=PI/2 -boundslonfi_glo(1,north_west)=2*PI -boundslatfi_glo(1,north_west)=PI/2 -boundslonfi_glo(1,south_west)=2*PI -boundslatfi_glo(1,south_west)=rlatv(1) -boundslonfi_glo(1,south_east)=0 -boundslatfi_glo(1,south_east)=rlatv(1) -DO j=2,jj - DO i=1,ii - k=(j-2)*ii+1+i - latfi_glo((j-2)*ii+1+i)= rlatu(j) - lonfi_glo((j-2)*ii+1+i)= rlonv(i) - cufi_glo((j-2)*ii+1+i) = cu((j-1)*(ii+1)+i) - cvfi_glo((j-2)*ii+1+i) = cv((j-1)*(ii+1)+i) - boundslonfi_glo(k,north_east)=rlonu(i) - boundslatfi_glo(k,north_east)=rlatv(j-1) - boundslonfi_glo(k,north_west)=rlonu(i+1) - boundslatfi_glo(k,north_west)=rlatv(j-1) - boundslonfi_glo(k,south_west)=rlonu(i+1) - boundslatfi_glo(k,south_west)=rlatv(j) - boundslonfi_glo(k,south_east)=rlonu(i) - boundslatfi_glo(k,south_east)=rlatv(j) - ENDDO -ENDDO -! South pole -latfi_glo(klon_glo)= rlatu(jj+1) -lonfi_glo(klon_glo)= 0. -cufi_glo(klon_glo) = cu((ii+1)*jj+1) -cvfi_glo(klon_glo) = cv((ii+1)*jj-ii) -boundslonfi_glo(klon_glo,north_east)= 0 -boundslatfi_glo(klon_glo,north_east)= rlatv(jj) -boundslonfi_glo(klon_glo,north_west)= 2*PI -boundslatfi_glo(klon_glo,north_west)= rlatv(jj) -boundslonfi_glo(klon_glo,south_west)= 2*PI -boundslatfi_glo(klon_glo,south_west)= -PI/2 -boundslonfi_glo(klon_glo,south_east)= 0 -boundslatfi_glo(klon_glo,south_east)= -Pi/2 - -! build airefi(), mesh area on physics grid -CALL gr_dyn_fi(1,ii+1,jj+1,klon_glo,aire,airefi_glo) -! Poles are single points on physics grid -airefi_glo(1)=sum(aire(1:ii,1)) -airefi_glo(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_glo(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 -! Now generate local lon/lat/cu/cv/area arrays -allocate(latfi(klon_omp),lonfi(klon_omp),cufi(klon_omp),cvfi(klon_omp)) -allocate(airefi(klon_omp)) -allocate(boundslonfi(klon_omp,4)) -allocate(boundslatfi(klon_omp,4)) -!call initcomgeomphy - -offset=klon_mpi_begin-1 -airefi(1:klon_omp)=airefi_glo(offset+klon_omp_begin:offset+klon_omp_end) -cufi(1:klon_omp)=cufi_glo(offset+klon_omp_begin:offset+klon_omp_end) -cvfi(1:klon_omp)=cvfi_glo(offset+klon_omp_begin:offset+klon_omp_end) -lonfi(1:klon_omp)=lonfi_glo(offset+klon_omp_begin:offset+klon_omp_end) -latfi(1:klon_omp)=latfi_glo(offset+klon_omp_begin:offset+klon_omp_end) -boundslonfi(1:klon_omp,:)=boundslonfi_glo(offset+klon_omp_begin:offset+klon_omp_end,:) -boundslatfi(1:klon_omp,:)=boundslatfi_glo(offset+klon_omp_begin:offset+klon_omp_end,:) - -! copy over local grid longitudes and latitudes -CALL init_geometry(klon_omp,lonfi,latfi,boundslonfi,boundslatfi, & - airefi,cufi,cvfi) - ! copy some fundamental parameters to physics ! and do some initializations -call init_dimphy(klon_omp,nlayer) ! Initialize dimphy module call phys_state_var_init(klon_omp,nlayer,nqtot, & day_ini,hour_ini,punjours,ptimestep, & prad,pg,pr,pcpp) -call ini_fillgeom(klon_omp,latfi,lonfi,airefi) +call ini_fillgeom(klon_omp,latitude,longitude,cell_area) +! work is needed to put what is in comgeomfi_h in geometry_mod? + call conf_phys(klon_omp,nlayer,nqtot) @@ -229,5 +76,4 @@ !$OMP END PARALLEL - end subroutine iniphysiq