Changeset 2351 for LMDZ5/trunk/libf/dynlonlat_phylonlat/phylmd
- Timestamp:
- Aug 25, 2015, 5:14:59 PM (9 years ago)
- Location:
- LMDZ5/trunk/libf/dynlonlat_phylonlat/phylmd
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
LMDZ5/trunk/libf/dynlonlat_phylonlat/phylmd/ce0l.F90
r2349 r2351 26 26 USE filtreg_mod, ONLY: inifilr 27 27 USE iniphysiq_mod, ONLY: iniphysiq 28 USE mod_const_mpi, ONLY: comm_lmdz 28 29 #ifdef inca 29 30 USE indice_sol_mod, ONLY: nbsrf, is_oce, is_sic, is_ter, is_lic … … 106 107 CALL read_distrib() 107 108 CALL init_mod_hallo() 108 CALL Init_Phys_lmdz(iim,jjp1,llm,mpi_size,distrib_phys)109 CALL init_interface_dyn_phys()110 #else111 CALL Init_Phys_lmdz(iim,jjp1,llm,1,(/(jjm-1)*iim+2/))112 109 #endif 113 110 WRITE(lunout,*)'---> klon=',klon … … 126 123 127 124 CALL inifilr() 128 CALL iniphysiq(iim,jjm,llm,daysec,day_ini,dtphys/nsplit_phys, & 125 CALL iniphysiq(iim,jjm,llm, & 126 distrib_phys(mpi_rank),comm_lmdz, & 127 daysec,day_ini,dtphys/nsplit_phys, & 129 128 rlatu,rlatv,rlonu,rlonv,aire,cu,cv,rad,g,r,cpp,iflag_phys) 130 129 IF(pressure_exner) CALL test_disvert -
LMDZ5/trunk/libf/dynlonlat_phylonlat/phylmd/iniphysiq_mod.F90
r2347 r2351 6 6 CONTAINS 7 7 8 SUBROUTINE iniphysiq(ii,jj,nlayer,punjours, pdayref,ptimestep, & 8 SUBROUTINE iniphysiq(ii,jj,nlayer, & 9 nbp, communicator, & 10 punjours, pdayref,ptimestep, & 9 11 rlatu,rlatv,rlonu,rlonv,aire,cu,cv, & 10 12 prad,pg,pr,pcpp,iflag_phys) 11 USE dimphy, ONLY: klev ! number of atmospheric levels12 USE mod_grid_phy_lmdz, ONLY: klon_glo ! number of atmospheric columns13 ! (on full grid)13 USE dimphy, ONLY: init_dimphy 14 USE mod_grid_phy_lmdz, ONLY: klon_glo, & ! number of atmospheric columns (on full grid) 15 regular_lonlat ! regular longitude-latitude grid type 14 16 USE mod_phys_lmdz_para, ONLY: klon_omp, & ! number of columns (on local omp grid) 15 17 klon_omp_begin, & ! start index of local omp subgrid 16 18 klon_omp_end, & ! end index of local omp subgrid 17 19 klon_mpi_begin ! start indes of columns (on local mpi grid) 20 USE geometry_mod, ONLY : init_geometry 18 21 USE vertical_layers_mod, ONLY : init_vertical_layers 19 22 USE infotrac, ONLY: nqtot,nqo,nbtr,tname,ttext,type_trac,& … … 26 29 indnum_fn_num,index_trac,& 27 30 niso,ntraceurs_zone,ntraciso 31 #ifdef REPROBUS 32 USE CHEM_REP, ONLY : Init_chem_rep_phys 33 #endif 28 34 USE control_mod, ONLY: dayref,anneeref,day_step,nday,offline 29 USE comgeomphy, ONLY: initcomgeomphy, &30 airephy, & ! physics grid area (m2)31 cuphy, & ! cu coeff. (u_covariant = cu * u)32 cvphy, & ! cv coeff. (v_covariant = cv * v)33 rlond, & ! longitudes34 rlatd ! latitudes35 35 USE inifis_mod, ONLY: inifis 36 36 USE time_phylmdz_mod, ONLY: init_time … … 38 38 USE phystokenc_mod, ONLY: init_phystokenc 39 39 USE phyaqua_mod, ONLY: iniaqua 40 USE physics_distribution_mod, ONLY : init_physics_distribution 40 41 USE regular_lonlat_mod, ONLY : init_regular_lonlat, & 41 42 east, west, north, south, & 42 43 north_east, north_west, & 43 44 south_west, south_east 45 USE mod_interface_dyn_phys, ONLY : init_interface_dyn_phys 44 46 IMPLICIT NONE 45 47 … … 64 66 INTEGER, INTENT (IN) :: ii ! number of atmospheric columns along longitudes 65 67 INTEGER, INTENT (IN) :: jj ! number of atompsheric columns along latitudes 68 INTEGER, INTENT(IN) :: nbp ! number of physics columns for this MPI process 69 INTEGER, INTENT(IN) :: communicator ! MPI communicator 66 70 REAL, INTENT (IN) :: rlatu(jj+1) ! latitudes of the physics grid 67 71 REAL, INTENT (IN) :: rlatv(jj) ! latitude boundaries of the physics grid … … 76 80 77 81 INTEGER :: ibegin, iend, offset 78 INTEGER :: i,j 82 INTEGER :: i,j,k 79 83 CHARACTER (LEN=20) :: modname = 'iniphysiq' 80 84 CHARACTER (LEN=80) :: abort_message … … 86 90 87 91 ! global array, on full physics grid: 88 REAL,ALLOCATABLE :: latfi(:) 89 REAL,ALLOCATABLE :: lonfi(:) 90 REAL,ALLOCATABLE :: cufi(:) 91 REAL,ALLOCATABLE :: cvfi(:) 92 REAL,ALLOCATABLE :: airefi(:) 93 94 IF (nlayer/=klev) THEN 95 WRITE (lunout, *) 'nlayer = ', nlayer 96 WRITE (lunout, *) 'klev = ', klev 97 CALL abort_gcm(modname, 'Problem with dimensions', 1) 98 END IF 99 100 !call init_phys_lmdz(ii,jj+1,llm,1,(/(jj-1)*ii+2/)) 92 REAL,ALLOCATABLE :: latfi_glo(:) 93 REAL,ALLOCATABLE :: lonfi_glo(:) 94 REAL,ALLOCATABLE :: cufi_glo(:) 95 REAL,ALLOCATABLE :: cvfi_glo(:) 96 REAL,ALLOCATABLE :: airefi_glo(:) 97 REAL,ALLOCATABLE :: boundslonfi_glo(:,:) 98 REAL,ALLOCATABLE :: boundslatfi_glo(:,:) 99 100 ! local arrays, on given MPI/OpenMP domain: 101 REAL,ALLOCATABLE,SAVE :: latfi(:) 102 REAL,ALLOCATABLE,SAVE :: lonfi(:) 103 REAL,ALLOCATABLE,SAVE :: cufi(:) 104 REAL,ALLOCATABLE,SAVE :: cvfi(:) 105 REAL,ALLOCATABLE,SAVE :: airefi(:) 106 REAL,ALLOCATABLE,SAVE :: boundslonfi(:,:) 107 REAL,ALLOCATABLE,SAVE :: boundslatfi(:,:) 108 !$OMP THREADPRIVATE (latfi,lonfi,cufi,cvfi,airefi,boundslonfi,boundslatfi) 109 110 ! Initialize Physics distibution and parameters and interface with dynamics 111 CALL init_physics_distribution(regular_lonlat,4, & 112 nbp,ii,jj+1,nlayer,communicator) 113 CALL init_interface_dyn_phys 101 114 102 115 ! init regular global longitude-latitude grid points and boundaries … … 123 136 124 137 ! Generate global arrays on full physics grid 125 ALLOCATE(latfi(klon_glo),lonfi(klon_glo),cufi(klon_glo),cvfi(klon_glo)) 126 ALLOCATE(airefi(klon_glo)) 138 ALLOCATE(latfi_glo(klon_glo),lonfi_glo(klon_glo)) 139 ALLOCATE(cufi_glo(klon_glo),cvfi_glo(klon_glo)) 140 ALLOCATE(airefi_glo(klon_glo)) 141 ALLOCATE(boundslonfi_glo(klon_glo,4)) 142 ALLOCATE(boundslatfi_glo(klon_glo,4)) 143 127 144 128 145 IF (klon_glo>1) THEN ! general case 129 146 ! North pole 130 latfi(1)=rlatu(1) 131 lonfi(1)=0. 132 cufi(1) = cu(1) 133 cvfi(1) = cv(1) 147 latfi_glo(1)=rlatu(1) 148 lonfi_glo(1)=0. 149 cufi_glo(1) = cu(1) 150 cvfi_glo(1) = cv(1) 151 boundslonfi_glo(1,north_east)=0 152 boundslatfi_glo(1,north_east)=PI/2 153 boundslonfi_glo(1,north_west)=2*PI 154 boundslatfi_glo(1,north_west)=PI/2 155 boundslonfi_glo(1,south_west)=2*PI 156 boundslatfi_glo(1,south_west)=rlatv(1) 157 boundslonfi_glo(1,south_east)=0 158 boundslatfi_glo(1,south_east)=rlatv(1) 134 159 DO j=2,jj 135 160 DO i=1,ii 136 latfi((j-2)*ii+1+i)= rlatu(j) 137 lonfi((j-2)*ii+1+i)= rlonv(i) 138 cufi((j-2)*ii+1+i) = cu((j-1)*(ii+1)+i) 139 cvfi((j-2)*ii+1+i) = cv((j-1)*(ii+1)+i) 161 k=(j-2)*ii+1+i 162 latfi_glo(k)= rlatu(j) 163 lonfi_glo(k)= rlonv(i) 164 cufi_glo(k) = cu((j-1)*ii+1+i) 165 cvfi_glo(k) = cv((j-1)*ii+1+i) 166 boundslonfi_glo(k,north_east)=rlonu(i) 167 boundslatfi_glo(k,north_east)=rlatv(j-1) 168 boundslonfi_glo(k,north_west)=rlonu(i+1) 169 boundslatfi_glo(k,north_west)=rlatv(j-1) 170 boundslonfi_glo(k,south_west)=rlonu(i+1) 171 boundslatfi_glo(k,south_west)=rlatv(j) 172 boundslonfi_glo(k,south_east)=rlonu(i) 173 boundslatfi_glo(k,south_east)=rlatv(j) 140 174 ENDDO 141 175 ENDDO 142 176 ! South pole 143 latfi(klon_glo)= rlatu(jj+1) 144 lonfi(klon_glo)= 0. 145 cufi(klon_glo) = cu((ii+1)*jj+1) 146 cvfi(klon_glo) = cv((ii+1)*jj-ii) 147 148 ! build airefi(), mesh area on physics grid 149 CALL gr_dyn_fi(1,ii+1,jj+1,klon_glo,aire,airefi) 177 latfi_glo(klon_glo)= rlatu(jj+1) 178 lonfi_glo(klon_glo)= 0. 179 cufi_glo(klon_glo) = cu((ii+1)*jj+1) 180 cvfi_glo(klon_glo) = cv((ii+1)*jj-ii) 181 boundslonfi_glo(klon_glo,north_east)= 0 182 boundslatfi_glo(klon_glo,north_east)= rlatv(jj) 183 boundslonfi_glo(klon_glo,north_west)= 2*PI 184 boundslatfi_glo(klon_glo,north_west)= rlatv(jj) 185 boundslonfi_glo(klon_glo,south_west)= 2*PI 186 boundslatfi_glo(klon_glo,south_west)= -PI/2 187 boundslonfi_glo(klon_glo,south_east)= 0 188 boundslatfi_glo(klon_glo,south_east)= -Pi/2 189 190 ! build airefi_glo(), mesh area on physics grid 191 CALL gr_dyn_fi(1,ii+1,jj+1,klon_glo,aire,airefi_glo) 150 192 ! Poles are single points on physics grid 151 airefi (1)=sum(aire(1:ii,1))152 airefi (klon_glo)=sum(aire(1:ii,jj+1))193 airefi_glo(1)=sum(aire(1:ii,1)) 194 airefi_glo(klon_glo)=sum(aire(1:ii,jj+1)) 153 195 154 196 ! Sanity check: do total planet area match between physics and dynamics? 155 197 total_area_dyn=sum(aire(1:ii,1:jj+1)) 156 total_area_phy=sum(airefi (1:klon_glo))198 total_area_phy=sum(airefi_glo(1:klon_glo)) 157 199 IF (total_area_dyn/=total_area_phy) THEN 158 200 WRITE (lunout, *) 'iniphysiq: planet total surface discrepancy !!!' … … 167 209 ELSE ! klon_glo==1, running the 1D model 168 210 ! just copy over input values 169 latfi(1)=rlatu(1) 170 lonfi(1)=rlonv(1) 171 cufi(1)=cu(1) 172 cvfi(1)=cv(1) 173 airefi(1)=aire(1,1) 211 latfi_glo(1)=rlatu(1) 212 lonfi_glo(1)=rlonv(1) 213 cufi_glo(1)=cu(1) 214 cvfi_glo(1)=cv(1) 215 airefi_glo(1)=aire(1,1) 216 boundslonfi_glo(1,north_east)=rlonu(1) 217 boundslatfi_glo(1,north_east)=PI/2 218 boundslonfi_glo(1,north_west)=rlonu(2) 219 boundslatfi_glo(1,north_west)=PI/2 220 boundslonfi_glo(1,south_west)=rlonu(2) 221 boundslatfi_glo(1,south_west)=rlatv(1) 222 boundslonfi_glo(1,south_east)=rlonu(1) 223 boundslatfi_glo(1,south_east)=rlatv(1) 174 224 ENDIF ! of IF (klon_glo>1) 175 225 176 226 !$OMP PARALLEL DEFAULT(SHARED) COPYIN(/temps/) 177 ! Initialize physical constants in physics: 178 CALL inifis(punjours,prad,pg,pr,pcpp) 179 CALL init_time(annee_ref,day_ref,day_ini,start_time,nday,ptimestep) 180 181 ! Copy over "offline" settings 182 CALL init_phystokenc(offline,istphy) 227 ! Now generate local lon/lat/cu/cv/area/bounds arrays 228 ALLOCATE(latfi(klon_omp),lonfi(klon_omp),cufi(klon_omp),cvfi(klon_omp)) 229 ALLOCATE(airefi(klon_omp)) 230 ALLOCATE(boundslonfi(klon_omp,4)) 231 ALLOCATE(boundslatfi(klon_omp,4)) 232 233 234 offset = klon_mpi_begin - 1 235 airefi(1:klon_omp) = airefi_glo(offset+klon_omp_begin:offset+klon_omp_end) 236 cufi(1:klon_omp) = cufi_glo(offset+klon_omp_begin:offset+klon_omp_end) 237 cvfi(1:klon_omp) = cvfi_glo(offset+klon_omp_begin:offset+klon_omp_end) 238 lonfi(1:klon_omp) = lonfi_glo(offset+klon_omp_begin:offset+klon_omp_end) 239 latfi(1:klon_omp) = latfi_glo(offset+klon_omp_begin:offset+klon_omp_end) 240 boundslonfi(1:klon_omp,:) = boundslonfi_glo(offset+klon_omp_begin:offset+klon_omp_end,:) 241 boundslatfi(1:klon_omp,:) = boundslatfi_glo(offset+klon_omp_begin:offset+klon_omp_end,:) 242 243 ! copy over local grid longitudes and latitudes 244 CALL init_geometry(klon_omp,lonfi,latfi,boundslonfi,boundslatfi, & 245 airefi,cufi,cvfi) 183 246 184 247 ! copy over preff , ap(), bp(), etc 185 248 CALL init_vertical_layers(nlayer,preff,scaleheight, & 186 249 ap,bp,presnivs,pseudoalt) 250 251 ! Initialize physical constants in physics: 252 CALL inifis(punjours,prad,pg,pr,pcpp) 253 254 CALL init_time(annee_ref,day_ref,day_ini,start_time,nday,ptimestep) 255 256 ! Initialize dimphy module 257 CALL Init_dimphy(klon_omp,nlayer) 258 259 ! Copy over "offline" settings 260 CALL init_phystokenc(offline,istphy) 187 261 188 262 ! Initialize tracer names, numbers, etc. for physics … … 197 271 niso,ntraceurs_zone,ntraciso) 198 272 199 ! Now generate local lon/lat/cu/cv/area arrays 200 CALL initcomgeomphy 201 202 offset = klon_mpi_begin - 1 203 airephy(1:klon_omp) = airefi(offset+klon_omp_begin:offset+klon_omp_end) 204 cuphy(1:klon_omp) = cufi(offset+klon_omp_begin:offset+klon_omp_end) 205 cvphy(1:klon_omp) = cvfi(offset+klon_omp_begin:offset+klon_omp_end) 206 rlond(1:klon_omp) = lonfi(offset+klon_omp_begin:offset+klon_omp_end) 207 rlatd(1:klon_omp) = latfi(offset+klon_omp_begin:offset+klon_omp_end) 273 ! Initializations for Reprobus 274 IF (type_trac == 'repr') THEN 275 #ifdef REPROBUS 276 CALL Init_chem_rep_phys(klon_omp,nlayer) 277 #endif 278 ENDIF 208 279 209 280 ! Additional initializations for aquaplanets 210 281 IF (iflag_phys>=100) THEN 211 CALL iniaqua(klon_omp, rlatd, rlond,iflag_phys)282 CALL iniaqua(klon_omp,iflag_phys) 212 283 END IF 213 284 !$OMP END PARALLEL
Note: See TracChangeset
for help on using the changeset viewer.