Changeset 4482 for LMDZ6/branches/LMDZ_ECRad/libf/dyn3d_common
- Timestamp:
- Mar 29, 2023, 3:14:27 PM (2 years ago)
- Location:
- LMDZ6/branches/LMDZ_ECRad
- Files:
-
- 8 edited
Legend:
- Unmodified
- Added
- Removed
-
LMDZ6/branches/LMDZ_ECRad
- Property svn:mergeinfo changed
-
LMDZ6/branches/LMDZ_ECRad/libf/dyn3d_common/comvert_mod.F90
r2602 r4482 10 10 11 11 PUBLIC :: ap,bp,presnivs,dpres,sig,ds,pa,preff,nivsigs,nivsig, & 12 aps,bps,scaleheight,pseudoalt,disvert_type, pressure_exner 12 aps,bps,scaleheight,pseudoalt,disvert_type, pressure_exner, & 13 presinter 13 14 14 15 REAL ap(llm+1) ! hybrid pressure contribution at interlayers 15 16 REAL bp (llm+1) ! hybrid sigma contribution at interlayer 16 17 REAL presnivs(llm) ! (reference) pressure at mid-layers 18 REAL presinter(llm+1) ! (reference) pressure at interlayers 17 19 REAL dpres(llm) 18 20 REAL sig(llm+1) -
LMDZ6/branches/LMDZ_ECRad/libf/dyn3d_common/conf_planete.F90
r2600 r4482 30 30 31 31 !Reference surface pressure (Pa) 32 preff=101325. 32 ! 101080 : specific value for CMIP5 aqua/terra planets 33 ! "Specify the initial dry mass to be equivalent to 34 ! a global mean surface pressure (101325 minus 245) Pa." 35 preff=101080. 33 36 CALL getin('preff', preff) 37 34 38 ! Reference pressure at which hybrid coord. become purely pressure 35 39 ! pa=50000. 36 40 pa=preff/2. 37 41 CALL getin('pa', pa) 42 38 43 ! Gravity 39 44 g=9.80665 45 40 46 CALL getin('g',g) 41 47 ! Molar mass of the atmosphere 48 42 49 molmass = 28.9644 43 50 CALL getin('molmass',molmass) 44 51 ! kappa=R/Cp et Cp 52 45 53 kappa = 2./7. 46 54 CALL getin('kappa',kappa) 55 47 56 cpp=8.3145/molmass/kappa*1000. 48 57 CALL getin('cpp',cpp) 49 58 ! Radius of the planet 59 50 60 rad = 6371229. 51 61 CALL getin('radius',rad) 62 52 63 ! Length of a standard day (s) 53 64 daysec=86400. 54 65 CALL getin('daysec',daysec) 66 55 67 ! Rotation rate of the planet: 56 68 ! Length of a solar day, in standard days 57 69 daylen = 1. 70 58 71 CALL getin('daylen',daylen) 59 72 ! Number of days (standard) per year: 73 60 74 year_day = 365.25 61 75 CALL getin('year_day',year_day) 62 76 ! Omega 63 77 ! omeg=2.*pi/86400. 78 64 79 omeg=2.*pi/daysec*(1./daylen+1./year_day) 65 80 CALL getin('omeg',omeg) -
LMDZ6/branches/LMDZ_ECRad/libf/dyn3d_common/control_mod.F90
r4146 r4482 29 29 INTEGER,SAVE :: ip_ebil_dyn 30 30 LOGICAL,SAVE :: offline 31 CHARACTER(len=4),SAVE :: config_inca32 31 CHARACTER(len=10),SAVE :: planet_type ! planet type ('earth','mars',...) 33 32 LOGICAL,SAVE :: output_grads_dyn ! output dynamics diagnostics in -
LMDZ6/branches/LMDZ_ECRad/libf/dyn3d_common/disvert.F90
r2786 r4482 11 11 use assert_m, only: assert 12 12 USE comvert_mod, ONLY: ap, bp, aps, bps, nivsigs, nivsig, dpres, presnivs, & 13 pseudoalt, pa, preff, scaleheight 13 pseudoalt, pa, preff, scaleheight, presinter 14 14 USE logic_mod, ONLY: ok_strato 15 15 … … 35 35 ! dpres(llm) !--- PRESSURE DIFFERENCE FOR EACH LAYER 36 36 ! presnivs(llm) !--- PRESSURE AT EACH MID-LAYER 37 ! presinter(llm+1) !--- PRESSURE AT EACH INTERFACE 37 38 ! scaleheight !--- VERTICAL SCALE HEIGHT (Earth: 8kms) 38 39 ! nivsig(llm+1) !--- SIGMA INDEX OF EACH LAYER INTERFACE … … 355 356 max(ap(l+1)+bp(l+1)*preff, 1.e-10)) 356 357 ENDDO 358 DO l=1, llmp1 359 presinter(l)= ( ap(l)+bp(l)*preff) 360 write(lunout, *)'PRESINTER(', l, ')=', presinter(l) 361 ENDDO 357 362 358 363 write(lunout, *) trim(modname),': PRESNIVS ' … … 395 400 ELSE IF(SIGN(f1,f4)/=f1) THEN; x2=xx; f2=f4 396 401 ELSE IF(SIGN(f2,f4)/=f2) THEN; x1=xx; f1=f4 397 ELSE; CALL abort_gcm("ridders",'Algorithm failed (which is odd...' )402 ELSE; CALL abort_gcm("ridders",'Algorithm failed (which is odd...', 1) 398 403 END IF 399 404 IF(ABS(10.*LOG(ABS(x2-x1)))<=1E-5) EXIT !--- ERROR ON SIG <= 0.01m -
LMDZ6/branches/LMDZ_ECRad/libf/dyn3d_common/grilles_gcm_netcdf_sub.F90
r3803 r4482 2 2 ! $Id: $ 3 3 ! 4 ! This subroutine creates the file grilles_gcm.nc containg longitudes and 5 ! latitudes in degrees for grid u and v. This subroutine is called from 6 ! ce0l. This subroutine corresponds to the first 7 ! part in the program create_fausse_var. 4 ! This subroutine creates the grilles_gcm.nc file, containing: 5 ! -> longitudes and latitudes in degrees for dynamical grids u, v and scalaire, 6 ! and the following variables added for INCA (informative anyway) 7 ! -> vertical levels "presnivs" 8 ! -> mask (land/sea), area (grid), phis=surface geopotential height = phis/g 8 9 ! 10 ! The subroutine is called in dynphy_lonlat/phylmd/ce0l.F90. 11 9 12 SUBROUTINE grilles_gcm_netcdf_sub(masque,phis) 10 13 11 14 USE comconst_mod, ONLY: cpp, kappa, g, omeg, daysec, rad, pi 12 15 USE comvert_mod, ONLY: presnivs, preff, pa 16 use netcdf, only: nf90_def_var, nf90_int, nf90_float, nf90_put_var 13 17 14 18 IMPLICIT NONE … … 19 23 INCLUDE "netcdf.inc" 20 24 21 25 !======================== 22 26 REAL,DIMENSION(iip1,jjp1),INTENT(IN) :: masque ! masque terre/mer 23 27 REAL,DIMENSION(iip1,jjp1),INTENT(IN) :: phis ! geopotentiel au sol 24 28 25 REAL temp(iim+1,jjm+1) 26 ! Attributs netcdf sortie 29 INTEGER status,i,j 30 31 ! Attributs netcdf output 27 32 INTEGER ncid_out,rcode_out 28 INTEGER out_lonuid,out_lonvid,out_latuid,out_latvid,out_levid 29 INTEGER out_varid 33 34 INTEGER out_lonuid,out_lonvid,out_latuid,out_latvid 35 INTEGER out_uid,out_vid,out_tempid 30 36 INTEGER out_lonudim,out_lonvdim 31 INTEGER out_latudim,out_latvdim,out_dim( 3)37 INTEGER out_latudim,out_latvdim,out_dim(2) 32 38 INTEGER out_levdim 33 34 INTEGER start(4),COUNT(4) 35 36 INTEGER status,i,j 37 REAL rlatudeg(jjp1),rlatvdeg(jjm),rlevdeg(llm) 39 ! 40 INTEGER :: presnivs_id 41 INTEGER :: mask_id,area_id,phis_id 42 ! 43 INTEGER start(2),COUNT(2) 44 45 ! Variables 46 REAL rlatudeg(jjp1),rlatvdeg(jjm),rlev(llm) 38 47 REAL rlonudeg(iip1),rlonvdeg(iip1) 39 40 REAL dlon1(iip1),dlon2(iip1),dlat1(jjp1),dlat2(jjp1) 41 REAL acoslat,dxkm,dykm,resol(iip1,jjp1) 42 REAL,DIMENSION(iip1,jjp1) :: phis_loc 48 REAL uwnd(iip1,jjp1),vwnd(iip1,jjm),temp(iip1,jjp1) 49 ! 43 50 INTEGER masque_int(iip1,jjp1) 44 INTEGER :: phis_id45 INTEGER :: area_id46 INTEGER :: mask_id47 INTEGER :: presnivs_id48 51 REAL :: phis_loc(iip1,jjp1) 52 53 !======================== 54 ! CALCULATION of latu, latv, lonu, lonv in deg. 55 ! --------------------------------------------------- 49 56 rad = 6400000 50 57 omeg = 7.272205e-05 … … 64 71 rlatudeg(j)=rlatu(j)*180./pi 65 72 ENDDO 73 66 74 DO j=1,jjm 67 75 rlatvdeg(j)=rlatv(j)*180./pi … … 72 80 rlonvdeg(i)=rlonv(i)*180./pi + 360. 73 81 ENDDO 74 75 76 ! 2 ----- OUVERTURE DE LA SORTIE NETCDF 77 ! --------------------------------------------------- 78 ! CREATION OUTPUT 79 ! ouverture fichier netcdf de sortie out 82 83 ! CALCULATION of "false" variables on u, v, s grids 84 ! --------------------------------------------------- 85 DO i=1,iip1 86 DO j=1,jjp1 87 uwnd(i,j)=MOD(i,2)+MOD(j,2) 88 temp(i,j)=MOD(i,2)+MOD(j,2) 89 ENDDO 90 DO j=1,jjm 91 vwnd(i,j)=MOD(i,2)+MOD(j,2) 92 END DO 93 ENDDO 94 95 ! CALCULATION of local vars for presnivs, mask, sfc. geopot. height 96 ! --------------------------------------------------- 97 rlev(:) = presnivs(:) 98 phis_loc(:,:) = phis(:,:)/g 99 masque_int(:,:) = nINT(masque(:,:)) 100 101 102 ! OPEN output netcdf file 103 !------------------------- 80 104 status=NF_CREATE('grilles_gcm.nc',IOR(NF_CLOBBER,NF_64BIT_OFFSET),ncid_out) 81 105 CALL handle_err(status) 106 107 ! DEFINE output dimensions 108 !------------------------- 82 109 status=NF_DEF_DIM(ncid_out,'lonu',iim+1,out_lonudim) 83 110 CALL handle_err(status) … … 88 115 status=NF_DEF_DIM(ncid_out,'latv',jjm,out_latvdim) 89 116 CALL handle_err(status) 90 91 92 ! Longitudes en u 93 status=NF_DEF_VAR(ncid_out,'lonu',NF_FLOAT,1,out_lonudim, out_lonuid) 117 ! 118 status=NF_DEF_DIM(ncid_out,'lev',llm,out_levdim) 119 CALL handle_err(status) 120 121 ! DEFINE output variables 122 !------------------------- 123 ! Longitudes on "u" dynamical grid 124 status=NF90_DEF_VAR(ncid_out,'lonu',NF90_FLOAT,out_lonudim, out_lonuid) 94 125 CALL handle_err(status) 95 126 status=NF_PUT_ATT_TEXT(ncid_out,out_lonuid,'units', 12,'degrees_east') 96 status=NF_PUT_ATT_TEXT(ncid_out,out_lonuid,'long_name',9,'Longitude en u') 97 98 ! Longitudes en v 99 status=NF_DEF_VAR(ncid_out,'lonv',NF_FLOAT,1,out_lonvdim, out_lonvid) 127 status=NF_PUT_ATT_TEXT(ncid_out,out_lonuid,'long_name',19,'Longitude on u grid') 128 ! Longitudes on "v" dynamical grid 129 status=NF90_DEF_VAR(ncid_out,'lonv',NF90_FLOAT,out_lonvdim, out_lonvid) 100 130 CALL handle_err(status) 101 131 status=NF_PUT_ATT_TEXT(ncid_out,out_lonvid,'units', 12,'degrees_east') 102 status=NF_PUT_ATT_TEXT(ncid_out,out_lonvid,'long_name', 9,'Longitude en v') 103 104 ! Latitude en u 105 status=NF_DEF_VAR(ncid_out,'latu',NF_FLOAT,1,out_latudim, out_latuid) 132 status=NF_PUT_ATT_TEXT(ncid_out,out_lonvid,'long_name', 19,'Longitude on v grid') 133 ! Latitudes on "u" dynamical grid 134 status=NF90_DEF_VAR(ncid_out,'latu',NF90_FLOAT,out_latudim, out_latuid) 106 135 CALL handle_err(status) 107 136 status=NF_PUT_ATT_TEXT(ncid_out,out_latuid,'units', 13,'degrees_north') 108 status=NF_PUT_ATT_TEXT(ncid_out,out_latuid,'long_name', 8,'Latitude en u') 109 110 ! Latitude en v 111 status=NF_DEF_VAR(ncid_out,'latv',NF_FLOAT,1,out_latvdim, out_latvid) 137 status=NF_PUT_ATT_TEXT(ncid_out,out_latuid,'long_name', 18,'Latitude on u grid') 138 ! Latitudes on "v" dynamical grid 139 status=NF90_DEF_VAR(ncid_out,'latv',NF90_FLOAT,out_latvdim, out_latvid) 112 140 CALL handle_err(status) 113 141 status=NF_PUT_ATT_TEXT(ncid_out,out_latvid,'units', 13,'degrees_north') 114 status=NF_PUT_ATT_TEXT(ncid_out,out_latvid,'long_name', 8,'Latitude en v') 115 116 ! ecriture de la grille u 142 status=NF_PUT_ATT_TEXT(ncid_out,out_latvid,'long_name', 18,'Latitude on v grid') 143 ! "u" lat/lon dynamical grid 117 144 out_dim(1)=out_lonudim 118 145 out_dim(2)=out_latudim 119 status=NF_DEF_VAR(ncid_out,'grille_u',NF_FLOAT,2,out_dim, out_varid) 120 CALL handle_err(status) 121 status=NF_PUT_ATT_TEXT(ncid_out,out_varid,'units', 6,'Kelvin') 122 status=NF_PUT_ATT_TEXT(ncid_out,out_varid,'long_name', 16,'Grille aux point u') 123 124 ! ecriture de la grille v 146 status=NF90_DEF_VAR(ncid_out,'grille_u',NF90_FLOAT,out_dim, out_uid) 147 CALL handle_err(status) 148 status=NF_PUT_ATT_TEXT(ncid_out,out_uid,'units', 3,'m/s') 149 status=NF_PUT_ATT_TEXT(ncid_out,out_uid,'long_name', 21,'u-wind dynamical grid') 150 ! "v" lat/lon dynamical grid 125 151 out_dim(1)=out_lonvdim 126 152 out_dim(2)=out_latvdim 127 status=NF_DEF_VAR(ncid_out,'grille_v',NF_FLOAT,2,out_dim, out_varid) 128 CALL handle_err(status) 129 status=NF_PUT_ATT_TEXT(ncid_out,out_varid,'units', 6,'Kelvin') 130 status=NF_PUT_ATT_TEXT(ncid_out,out_varid,'long_name', 16,'Grille aux point v') 131 132 ! ecriture de la grille u 153 status=NF90_DEF_VAR(ncid_out,'grille_v',NF90_FLOAT,out_dim, out_vid) 154 CALL handle_err(status) 155 status=NF_PUT_ATT_TEXT(ncid_out,out_vid,'units', 3,'m/s') 156 status=NF_PUT_ATT_TEXT(ncid_out,out_vid,'long_name', 21,'v-wind dynamical grid') 157 ! "s" (scalar) lat/lon dynamical grid 133 158 out_dim(1)=out_lonvdim 134 159 out_dim(2)=out_latudim 135 status=NF_DEF_VAR(ncid_out,'grille_s',NF_FLOAT,2,out_dim, out_varid) 136 CALL handle_err(status) 137 status=NF_PUT_ATT_TEXT(ncid_out,out_varid,'units', 6,'Kelvin') 138 status=NF_PUT_ATT_TEXT(ncid_out,out_varid,'long_name',16,'Grille aux point u') 139 140 status=NF_ENDDEF(ncid_out) 141 write(*,*) "COUCOU 6" 142 CALL handle_err(status) 143 ! 5) ----- FERMETURE DES FICHIERS NETCDF------------------ 144 ! -------------------------------------------------------- 145 ! 3-b- Ecriture de la grille pour la sortie 146 ! rajoute l'ecriture de la grille 147 148 #ifdef NC_DOUBLE 149 status=NF_PUT_VARA_DOUBLE(ncid_out,out_lonuid,1,iim+1,rlonudeg) 150 CALL handle_err(status) 151 status=NF_PUT_VARA_DOUBLE(ncid_out,out_lonvid,1,iim+1,rlonvdeg) 152 CALL handle_err(status) 153 status=NF_PUT_VARA_DOUBLE(ncid_out,out_latuid,1,jjm+1,rlatudeg) 154 CALL handle_err(status) 155 status=NF_PUT_VARA_DOUBLE(ncid_out,out_latvid,1,jjm,rlatvdeg) 156 CALL handle_err(status) 157 #else 158 status=NF_PUT_VARA_REAL(ncid_out,out_lonuid,1,iim+1,rlonudeg) 159 CALL handle_err(status) 160 status=NF_PUT_VARA_REAL(ncid_out,out_lonvid,1,iim+1,rlonvdeg) 161 CALL handle_err(status) 162 status=NF_PUT_VARA_REAL(ncid_out,out_latuid,1,jjm+1,rlatudeg) 163 CALL handle_err(status) 164 status=NF_PUT_VARA_REAL(ncid_out,out_latvid,1,jjm,rlatvdeg) 165 CALL handle_err(status) 166 #endif 167 168 start(1)=1 169 start(2)=1 170 start(3)=1 171 start(4)=1 172 173 COUNT(1)=iim+1 174 COUNT(2)=jjm+1 175 COUNT(3)=1 176 COUNT(4)=1 177 178 DO j=1,jjm+1 179 DO i=1,iim+1 180 temp(i,j)=MOD(i,2)+MOD(j,2) 181 ENDDO 182 ENDDO 183 184 #ifdef NC_DOUBLE 185 status=NF_PUT_VARA_DOUBLE(ncid_out,out_varid,start, count,temp) 186 CALL handle_err(status) 187 #else 188 status=NF_PUT_VARA_REAL(ncid_out,out_varid,start, count,temp) 189 CALL handle_err(status) 190 #endif 191 192 ! On re-ouvre le fichier pour rajouter 4 nouvelles variables necessaire pour INCA 193 ! lev - phis - aire - mask 194 ! rlevdeg(:) = presnivs 195 rlevdeg(:) = presnivs(:) 196 phis_loc(:,:) = phis(:,:)/g 197 198 ! niveaux de pression verticaux 199 status = NF_REDEF (ncid_out) 200 CALL handle_err(status) 201 status=NF_DEF_DIM(ncid_out,'lev',llm,out_levdim) 202 CALL handle_err(status) 203 status=NF_DEF_VAR(ncid_out,'presnivs',NF_FLOAT,1,out_levdim,& 204 presnivs_id) 205 CALL handle_err(status) 206 207 ! fields 160 status=NF90_DEF_VAR(ncid_out,'grille_s',NF90_FLOAT,out_dim, out_tempid) 161 CALL handle_err(status) 162 status=NF_PUT_ATT_TEXT(ncid_out,out_tempid,'units', 6,'Kelvin') 163 status=NF_PUT_ATT_TEXT(ncid_out,out_tempid,'long_name',21,'scalar dynamical grid') 164 ! 165 ! for INCA : 166 ! vertical levels "presnivs" 167 status=NF90_DEF_VAR(ncid_out,'presnivs',NF90_FLOAT,out_levdim, presnivs_id) 168 CALL handle_err(status) 169 status=NF_PUT_ATT_TEXT(ncid_out,presnivs_id,'units', 2,'Pa') 170 status=NF_PUT_ATT_TEXT(ncid_out,presnivs_id,'long_name',15,'Vertical levels') 171 ! surface geopotential height: named "phis" as the sfc geopotential, is actually phis/g 208 172 out_dim(1)=out_lonvdim 209 173 out_dim(2)=out_latudim 210 211 status = nf_def_var(ncid_out,'phis',NF_FLOAT,2,out_dim,phis_id) 212 CALL handle_err(status) 213 status = nf_def_var(ncid_out,'aire',NF_FLOAT,2,out_dim,area_id) 214 CALL handle_err(status) 215 status = nf_def_var(ncid_out,'mask',NF_INT ,2,out_dim,mask_id) 216 CALL handle_err(status) 217 218 status=NF_ENDDEF(ncid_out) 219 CALL handle_err(status) 220 221 ! ecriture des variables 222 #ifdef NC_DOUBLE 223 status=NF_PUT_VARA_DOUBLE(ncid_out,presnivs_id,1,llm,rlevdeg) 224 CALL handle_err(status) 225 #else 226 status=NF_PUT_VARA_REAL(ncid_out,out_levid,1,llm,rlevdeg) 227 CALL handle_err(status) 228 #endif 229 230 start(1)=1 231 start(2)=1 232 start(3)=1 233 start(4)=0 174 status = nf90_def_var(ncid_out,'phis',NF90_FLOAT,out_dim,phis_id) 175 CALL handle_err(status) 176 status=NF_PUT_ATT_TEXT(ncid_out,phis_id,'units', 1,'m') 177 status=NF_PUT_ATT_TEXT(ncid_out,phis_id,'long_name',27,'surface geopotential height') 178 ! gridcell area 179 status = nf90_def_var(ncid_out,'aire',NF90_FLOAT,out_dim,area_id) 180 CALL handle_err(status) 181 status=NF_PUT_ATT_TEXT(ncid_out,area_id,'units', 2,'m2') 182 status=NF_PUT_ATT_TEXT(ncid_out,area_id,'long_name',13,'gridcell area') 183 ! land-sea mask (nearest integer approx) 184 status = nf90_def_var(ncid_out,'mask',NF90_INT,out_dim,mask_id) 185 CALL handle_err(status) 186 status=NF_PUT_ATT_TEXT(ncid_out,mask_id,'long_name',27,'land-sea mask (nINT approx)') 187 188 ! END the 'define' mode in netCDF file 189 status=NF_ENDDEF(ncid_out) 190 CALL handle_err(status) 191 192 ! WRITE the variables 193 !------------------------- 194 ! 1D : lonu, lonv,latu,latv ; INCA : presnivs 195 status=NF90_PUT_VAR(ncid_out,out_lonuid,rlonudeg,[1],[iip1]) 196 CALL handle_err(status) 197 status=NF90_PUT_VAR(ncid_out,out_lonvid,rlonvdeg,[1],[iip1]) 198 CALL handle_err(status) 199 status=NF90_PUT_VAR(ncid_out,out_latuid,rlatudeg,[1],[jjp1]) 200 CALL handle_err(status) 201 status=NF90_PUT_VAR(ncid_out,out_latvid,rlatvdeg,[1],[jjm]) 202 CALL handle_err(status) 203 status=NF90_PUT_VAR(ncid_out,presnivs_id,rlev,[1],[llm]) 204 CALL handle_err(status) 205 206 ! 2D : grille_u,grille_v,grille_s ; INCA: phis,aire,mask 207 start(:)=1 234 208 COUNT(1)=iip1 235 COUNT(2)=jjp1 236 COUNT(3)=1 237 COUNT(4)=0 238 239 status = nf_put_vara_double(ncid_out, phis_id,start,count, phis_loc) 240 CALL handle_err(status) 241 status = nf_put_vara_double(ncid_out, area_id,start,count, aire) 242 CALL handle_err(status) 243 masque_int(:,:) = nINT(masque(:,:)) 244 status = nf_put_vara_int(ncid_out, mask_id,start,count,masque_int) 245 CALL handle_err(status) 246 247 ! fermeture du fichier netcdf 209 210 COUNT(2)=jjp1 ! for "u" and "s" grids 211 status=NF90_PUT_VAR(ncid_out,out_uid,uwnd,start, count) 212 CALL handle_err(status) 213 COUNT(2)=jjm ! for "v" grid 214 status=NF90_PUT_VAR(ncid_out,out_vid,vwnd,start, count) 215 CALL handle_err(status) 216 COUNT(2)=jjp1 ! as "s" grid, for all the following vars 217 status=NF90_PUT_VAR(ncid_out,out_tempid,temp,start, count) 218 CALL handle_err(status) 219 status = nf90_put_var(ncid_out, phis_id, phis_loc,start,count) 220 CALL handle_err(status) 221 status = nf90_put_var(ncid_out, area_id, aire,start,count) 222 CALL handle_err(status) 223 status = nf90_put_var(ncid_out, mask_id,masque_int,start,count) 224 CALL handle_err(status) 225 226 ! CLOSE netcdf file 248 227 CALL ncclos(ncid_out,rcode_out) 228 write(*,*) "END grilles_gcm_netcdf_sub OK" 249 229 250 230 END SUBROUTINE grilles_gcm_netcdf_sub 251 252 231 253 232 -
LMDZ6/branches/LMDZ_ECRad/libf/dyn3d_common/infotrac.F90
r4203 r4482 3 3 MODULE infotrac 4 4 5 USE strings_mod, ONLY: msg, find, strIdx, strFind, strParse, dispTable, int2str, reduceExpr, & 6 cat, fmsg, test, strTail, strHead, strStack, strReduce, bool2str, maxlen, testFile 7 USE readTracFiles_mod, ONLY: trac_type, readTracersFiles, addPhase, indexUpdate, nphases, ancestor, & 8 isot_type, old2newName, delPhase, getKey_init, tran0, & 9 keys_type, initIsotopes, getPhase, known_phases, getKey, setGeneration, & 10 maxTableWidth 5 USE strings_mod, ONLY: msg, fmsg, maxlen, cat, dispTable, int2str, bool2str, strStack, strParse 6 USE readTracFiles_mod, ONLY: trac_type, readTracersFiles, tracers, setGeneration, itZonIso, nzone, tran0, isoZone, & 7 delPhase, niso, getKey, isot_type, readIsotopesFile, isotope, maxTableWidth, iqIsoPha, nphas, ixIso, isoPhas, & 8 addPhase, iH2O, nbIso, isoSelect, testTracersFiles, isoKeys, indexUpdate, isoCheck, nzone, ntiso, isoName, & 9 addKey 11 10 IMPLICIT NONE 12 11 … … 14 13 15 14 !=== FOR TRACERS: 16 PUBLIC :: in fotrac_init!--- Initialization of the tracers17 PUBLIC :: tracers, type_trac , types_trac!--- Full tracers database, tracers type keyword15 PUBLIC :: init_infotrac !--- Initialization of the tracers 16 PUBLIC :: tracers, type_trac !--- Full tracers database, tracers type keyword 18 17 PUBLIC :: nqtot, nbtr, nqo, nqCO2, nqtottr !--- Main dimensions 19 18 PUBLIC :: conv_flg, pbl_flg !--- Convection & boundary layer activation keys 20 19 21 20 !=== FOR ISOTOPES: General 22 PUBLIC :: isot opes,nbIso !--- Derived type, full isotopes families database + nb of families21 PUBLIC :: isot_type, nbIso !--- Derived type, full isotopes families database + nb of families 23 22 PUBLIC :: isoSelect, ixIso !--- Isotopes family selection tool + selected family index 24 23 !=== FOR ISOTOPES: Specific to water 25 PUBLIC :: iH2O , tnat, alpha_ideal !--- H2O isotopes index, natural abundance, fractionning coeff.24 PUBLIC :: iH2O !--- H2O isotopes class index 26 25 PUBLIC :: min_qParent, min_qMass, min_ratio !--- Min. values for various isotopic quantities 27 26 !=== FOR ISOTOPES: Depending on the selected isotopes family … … 34 33 !=== FOR BOTH TRACERS AND ISOTOPES 35 34 PUBLIC :: getKey !--- Get a key from "tracers" or "isotope" 36 37 INTERFACE isoSelect; MODULE PROCEDURE isoSelectByIndex, isoSelectByName; END INTERFACE isoSelect38 35 39 36 !=== CONVENTIONS FOR TRACERS NUMBERS: … … 71 68 ! | phase | Phases list ("g"as / "l"iquid / "s"olid) | / | [g][l][s] | 72 69 ! | component | Name(s) of the merged/cumulated section(s) | / | coma-separated names | 73 ! | iadv | Advection scheme number | iadv | 1-20,30 exc. 3-9,15,19 |74 70 ! | iGeneration | Generation (>=1) | / | | 75 ! | isAdvected | advected tracers flag (.TRUE. if iadv >= 0) | / | nqtrue .TRUE. values |76 ! | isInPhysics | tracers not extracted from the main table in physics | / | nqtottr .TRUE. values |77 71 ! | iqParent | Index of the parent tracer | iqpere | 1:nqtot | 78 72 ! | iqDescen | Indexes of the childs (all generations) | iqfils | 1:nqtot | 79 73 ! | nqDescen | Number of the descendants (all generations) | nqdesc | 1:nqtot | 80 ! | nqChilds | Number of childs (1st generation only) | nqfils | 1:nqtot | 74 ! | nqChildren | Number of childs (1st generation only) | nqfils | 1:nqtot | 75 ! | keys | key/val pairs accessible with "getKey" routine | / | | 76 ! | iadv | Advection scheme number | iadv | 1,2,10-20(exc.15,19),30| 77 ! | isAdvected | advected tracers flag (.TRUE. if iadv >= 0) | / | nqtrue .TRUE. values | 78 ! | isInPhysics | tracers not extracted from the main table in physics | / | nqtottr .TRUE. values | 81 79 ! | iso_iGroup | Isotopes group index in isotopes(:) | / | 1:nbIso | 82 80 ! | iso_iName | Isotope name index in isotopes(iso_iGroup)%trac(:) | iso_indnum | 1:niso | 83 81 ! | iso_iZone | Isotope zone index in isotopes(iso_iGroup)%zone(:) | zone_num | 1:nzone | 84 82 ! | iso_iPhas | Isotope phase index in isotopes(iso_iGroup)%phas(:) | phase_num | 1:nphas | 85 ! | keys | key/val pairs accessible with "getKey" routine | / | |86 83 ! +-------------+------------------------------------------------------+-------------+------------------------+ 87 84 ! … … 103 100 104 101 !=== DIMENSIONS OF THE TRACERS TABLES AND OTHER SCALAR VARIABLES 105 INTEGER, SAVE :: nqtot, & !--- Tracers nb in dynamics (incl. higher moments + H2O) 106 nbtr, & !--- Tracers nb in physics (excl. higher moments + H2O) 107 nqo, & !--- Number of water phases 108 nbIso, & !--- Number of available isotopes family 109 nqtottr, & !--- Number of tracers passed to phytrac (TO BE DELETED ?) 110 nqCO2 !--- Number of tracers of CO2 (ThL) 111 CHARACTER(LEN=maxlen), SAVE :: type_trac !--- Keyword for tracers type(s) 112 CHARACTER(LEN=maxlen), SAVE, ALLOCATABLE :: types_trac(:) !--- Keyword for tracers type(s), parsed version 113 114 !=== DERIVED TYPES EMBEDDING MOST INFORMATIONS ABOUT TRACERS AND ISOTOPES FAMILIES 115 TYPE(trac_type), TARGET, SAVE, ALLOCATABLE :: tracers(:) !=== TRACERS DESCRIPTORS VECTOR 116 TYPE(isot_type), TARGET, SAVE, ALLOCATABLE :: isotopes(:) !=== ISOTOPES PARAMETERS VECTOR 117 118 !=== ALIASES FOR CURRENTLY SELECTED ISOTOPES FAMILY OF VARIABLES EMBEDDED IN THE VECTOR "isotopes" 119 TYPE(isot_type), SAVE, POINTER :: isotope !--- CURRENTLY SELECTED ISOTOPES FAMILY DESCRIPTOR 120 INTEGER, SAVE :: ixIso, iH2O !--- Index of the selected isotopes family and H2O family 121 LOGICAL, SAVE :: isoCheck !--- Flag to trigger the checking routines 122 TYPE(keys_type), SAVE, POINTER :: isoKeys(:) !--- ONE SET OF KEYS FOR EACH ISOTOPE (LISTED IN isoName) 123 CHARACTER(LEN=maxlen), SAVE, POINTER :: isoName(:), & !--- ISOTOPES NAMES FOR THE CURRENTLY SELECTED FAMILY 124 isoZone(:), & !--- TAGGING ZONES FOR THE CURRENTLY SELECTED FAMILY 125 isoPhas !--- USED PHASES FOR THE CURRENTLY SELECTED FAMILY 126 INTEGER, SAVE :: niso, nzone, & !--- NUMBER OF ISOTOPES, TAGGING ZONES AND PHASES 127 nphas, ntiso !--- NUMBER OF PHASES AND ISOTOPES + ISOTOPIC TAGGING TRACERS 128 INTEGER, SAVE, POINTER ::itZonIso(:,:), & !--- INDEX IN "isoTrac" AS f(tagging zone idx, isotope idx) 129 iqIsoPha(:,:) !--- INDEX IN "qx" AS f(isotopic tracer idx, phase idx) 130 131 !=== VARIABLES FOR ISOTOPES INITIALIZATION AND FOR INCA 132 REAL, SAVE, ALLOCATABLE :: tnat(:), & !--- Natural relative abundance of water isotope (niso) 133 alpha_ideal(:) !--- Ideal fractionning coefficient (for initial state) (niso) 134 INTEGER, SAVE, ALLOCATABLE :: conv_flg(:), & !--- Convection activation ; needed for INCA (nbtr) 135 pbl_flg(:) !--- Boundary layer activation ; needed for INCA (nbtr) 102 INTEGER, SAVE :: nqtot, & !--- Tracers nb in dynamics (incl. higher moments + H2O) 103 nbtr, & !--- Tracers nb in physics (excl. higher moments + H2O) 104 nqo, & !--- Number of water phases 105 nqtottr, & !--- Number of tracers passed to phytrac (TO BE DELETED ?) 106 nqCO2 !--- Number of tracers of CO2 (ThL) 107 CHARACTER(LEN=maxlen), SAVE :: type_trac !--- Keyword for tracers type 108 109 !=== VARIABLES FOR INCA 110 INTEGER, SAVE, ALLOCATABLE :: conv_flg(:), & !--- Convection activation ; needed for INCA (nbtr) 111 pbl_flg(:) !--- Boundary layer activation ; needed for INCA (nbtr) 136 112 137 113 CONTAINS 138 114 139 SUBROUTINE in fotrac_init140 USE control_mod, ONLY: planet_type , config_inca115 SUBROUTINE init_infotrac 116 USE control_mod, ONLY: planet_type 141 117 #ifdef REPROBUS 142 118 USE CHEM_REP, ONLY: Init_chem_rep_trac … … 176 152 CHARACTER(LEN=2) :: suff(9) !--- Suffixes for schemes of order 3 or 4 (Prather) 177 153 CHARACTER(LEN=3) :: descrq(30) !--- Advection scheme description tags 178 CHARACTER(LEN=maxlen) :: msg1 !--- String for messages154 CHARACTER(LEN=maxlen) :: msg1, texp, ttp !--- Strings for messages and expanded tracers type 179 155 INTEGER :: fType !--- Tracers description file type ; 0: none 180 156 !--- 1/2/3: "traceur.def"/"tracer.def"/"tracer_*.def" 181 157 INTEGER :: nqtrue !--- Tracers nb from tracer.def (no higher order moments) 182 158 INTEGER :: iad !--- Advection scheme number 183 INTEGER :: ic, i p, np, iq, jq, it, nt, im, nm, ix, iz, nz, k!--- Indexes and temporary variables159 INTEGER :: ic, iq, jq, it, nt, im, nm, iz, k !--- Indexes and temporary variables 184 160 LOGICAL :: lerr, ll 185 161 CHARACTER(LEN=1) :: p 186 162 TYPE(trac_type), ALLOCATABLE, TARGET :: ttr(:) 187 163 TYPE(trac_type), POINTER :: t1, t(:) 188 TYPE(isot_type), POINTER :: iso189 190 CHARACTER(LEN=maxlen), ALLOCATABLE :: tnom_0(:), tnom_transp(:) !--- Tracer short name + transporting fluid name191 CHARACTER(LEN=maxlen) :: tchaine192 164 INTEGER :: ierr 193 165 194 CHARACTER(LEN=*), PARAMETER :: modname="in fotrac_init"166 CHARACTER(LEN=*), PARAMETER :: modname="init_infotrac" 195 167 !------------------------------------------------------------------------------------------------------------------------------ 196 168 ! Initialization : … … 202 174 203 175 CALL msg('type_trac = "'//TRIM(type_trac)//'"', modname) 204 IF(strParse(type_trac, '|', types_trac, n=nt)) CALL abort_gcm(modname,'can''t parse "type_trac = '//TRIM(type_trac)//'"',1) 205 206 !--------------------------------------------------------------------------------------------------------------------------- 207 DO it = 1, nt !--- nt>1=> "type_trac": coma-separated keywords list 208 !--------------------------------------------------------------------------------------------------------------------------- 209 !--- MESSAGE ABOUT THE CHOSEN CONFIGURATION 210 msg1 = 'For type_trac = "'//TRIM(types_trac(it))//'":' 211 SELECT CASE(types_trac(it)) 212 CASE('inca'); CALL msg(TRIM(msg1)//' coupling with INCA chemistry model, config_inca='//config_inca, modname) 213 CASE('inco'); CALL msg(TRIM(msg1)//' coupling jointly with INCA and CO2 cycle', modname) 214 CASE('repr'); CALL msg(TRIM(msg1)//' coupling with REPROBUS chemistry model', modname) 215 CASE('co2i'); CALL msg(TRIM(msg1)//' you have chosen to run with CO2 cycle', modname) 216 CASE('coag'); CALL msg(TRIM(msg1)//' tracers are treated for COAGULATION tests', modname) 217 CASE('lmdz'); CALL msg(TRIM(msg1)//' tracers are treated in LMDZ only', modname) 218 CASE DEFAULT; CALL abort_gcm(modname,'type_trac='//TRIM(types_trac(it))//' not possible yet.',1) 219 END SELECT 220 221 !--- COHERENCE TEST BETWEEN "type_trac" AND "config_inca" 222 IF(ANY(['inca', 'inco'] == types_trac(it)) .AND. ALL(['aero', 'aeNP', 'chem'] /= config_inca)) & 223 CALL abort_gcm(modname, 'Incoherence between type_trac and config_inca. Please modify "run.def"',1) 224 225 !--- COHERENCE TEST BETWEEN "type_trac" AND PREPROCESSING KEYS 226 SELECT CASE(types_trac(it)) 227 CASE('inca', 'inco') 176 177 !--- MESSAGE ABOUT THE CHOSEN CONFIGURATION 178 msg1 = 'For type_trac = "'//TRIM(type_trac)//'":' 179 SELECT CASE(type_trac) 180 CASE('inca'); CALL msg(TRIM(msg1)//' coupling with INCA chemistry model', modname) 181 CASE('inco'); CALL msg(TRIM(msg1)//' coupling jointly with INCA and CO2 cycle', modname) 182 CASE('repr'); CALL msg(TRIM(msg1)//' coupling with REPROBUS chemistry model', modname) 183 CASE('co2i'); CALL msg(TRIM(msg1)//' you have chosen to run with CO2 cycle', modname) 184 CASE('coag'); CALL msg(TRIM(msg1)//' tracers are treated for COAGULATION tests', modname) 185 CASE('lmdz'); CALL msg(TRIM(msg1)//' tracers are treated in LMDZ only', modname) 186 CASE DEFAULT; CALL abort_gcm(modname,'type_trac='//TRIM(type_trac)//' not possible yet.',1) 187 END SELECT 188 189 !--- COHERENCE TEST BETWEEN "type_trac" AND PREPROCESSING KEYS 190 SELECT CASE(type_trac) 191 CASE('inca', 'inco') 228 192 #ifndef INCA 229 230 #endif 231 193 CALL abort_gcm(modname, 'You must add cpp key INCA and compile with INCA code', 1) 194 #endif 195 CASE('repr') 232 196 #ifndef REPROBUS 233 234 #endif 235 197 CALL abort_gcm(modname, 'You must add cpp key REPROBUS and compile with REPROBUS code', 1) 198 #endif 199 CASE('coag') 236 200 #ifndef CPP_StratAer 237 CALL abort_gcm(modname, 'You must add cpp key StratAer and compile with StratAer code', 1) 238 #endif 239 END SELECT 240 241 !--------------------------------------------------------------------------------------------------------------------------- 242 END DO 243 !--------------------------------------------------------------------------------------------------------------------------- 244 245 !--- DISABLE "config_inca" OPTION FOR A RUN WITHOUT "INCA" IF IT DIFFERS FROM "none" 246 IF(fmsg('Setting config_inca="none" as you do not couple with INCA model', & 247 modname, ALL(types_trac /= 'inca') .AND. ALL(types_trac /= 'inco') .AND. config_inca /= 'none')) config_inca = 'none' 248 249 nqCO2 = 0; IF(ANY(types_trac == 'inco')) nqCO2 = 1 201 CALL abort_gcm(modname, 'You must add cpp key StratAer and compile with StratAer code', 1) 202 #endif 203 END SELECT 204 205 nqCO2 = COUNT( [type_trac == 'inco', type_trac == 'co2i'] ) 250 206 251 207 !============================================================================================================================== 252 208 ! 1) Get the numbers of: true (first order only) tracers "nqtrue", water tracers "nqo" (vapor/liquid/solid) 253 209 !============================================================================================================================== 254 IF(readTracersFiles(type_trac, fType, tracers)) CALL abort_gcm(modname, 'problem with tracers file(s)',1) 210 texp = type_trac !=== EXPANDED VERSION OF "type_trac", WITH "|" SEPARATOR 211 IF(texp == 'inco') texp = 'co2i|inca' 212 IF(texp /= 'lmdz') texp = 'lmdz|'//TRIM(texp) 213 214 !=== DETERMINE THE TYPE OF THE INPUT TRACERS DESCRIPTION FILE 215 IF(testTracersFiles(modname, texp, fType, .TRUE.)) CALL abort_gcm(modname, 'problem with tracers file(s)',1) 216 ttp = type_trac; IF(fType /= 1) ttp = texp 217 218 IF(readTracersFiles(ttp, type_trac == 'repr')) CALL abort_gcm(modname, 'problem with tracers file(s)',1) 219 !--------------------------------------------------------------------------------------------------------------------------- 255 220 IF(fType == 0) CALL abort_gcm(modname, 'Missing tracers file: "traceur.def", "tracer.def" or "tracer_<keyword>.def file.',1) 256 221 !--------------------------------------------------------------------------------------------------------------------------- 257 IF(fType == 1 .AND. ANY(['inca','inco'] == type_trac)) THEN !=== FOUND OLD STYLE INCA "traceur.def" (single type_trac)222 IF(fType == 1 .AND. ANY(['inca','inco']==type_trac)) THEN !=== FOUND OLD STYLE INCA "traceur.def" 258 223 !--------------------------------------------------------------------------------------------------------------------------- 259 224 #ifdef INCA 260 nqo = SIZE(tracers) 261 IF(nqCO2==1 .AND. nqo==4) nqo = 3 !--- Force nqo to 3 (ThL) 225 nqo = SIZE(tracers) - nqCO2 262 226 CALL Init_chem_inca_trac(nqINCA) !--- Get nqINCA from INCA 263 227 nbtr = nqINCA + nqCO2 !--- Number of tracers passed to phytrac … … 268 232 CALL init_transport(solsym_inca, conv_flg_inca, pbl_flg_inca, hadv_inca, vadv_inca) 269 233 ALLOCATE(ttr(nqtrue)) 270 ttr(1:nqo+nqCO2) 271 ttr(1 : nqo )%component 272 ttr(1+nqo:nqCO2+nqo )%component 273 ttr(1+nqo+nqCO2:nqtrue)%component 274 ttr(1+nqo :nqtrue)%name 275 ttr(1+nqo+nqCO2:nqtrue)%parent 276 ttr(1+nqo+nqCO2:nqtrue)%phase 234 ttr(1:nqo+nqCO2) = tracers 235 ttr(1 : nqo )%component = 'lmdz' 236 ttr(1+nqo:nqCO2+nqo )%component = 'co2i' 237 ttr(1+nqo+nqCO2:nqtrue)%component = 'inca' 238 ttr(1+nqo :nqtrue)%name = [('CO2 ', k=1, nqCO2), solsym_inca] 239 ttr(1+nqo+nqCO2:nqtrue)%parent = tran0 240 ttr(1+nqo+nqCO2:nqtrue)%phase = 'g' 277 241 lerr = getKey('hadv', had, ky=tracers(:)%keys) 278 242 lerr = getKey('vadv', vad, ky=tracers(:)%keys) 279 hadv(1:nqo ) = had(:); hadv(nqo+1:nqtrue) = hadv_inca280 vadv(1:nqo ) = vad(:); vadv(nqo+1:nqtrue) = vadv_inca243 hadv(1:nqo+nqCO2) = had(:); hadv(1+nqo+nqCO2:nqtrue) = hadv_inca 244 vadv(1:nqo+nqCO2) = vad(:); vadv(1+nqo+nqCO2:nqtrue) = vadv_inca 281 245 CALL MOVE_ALLOC(FROM=ttr, TO=tracers) 282 CALL setGeneration(tracers) !--- SET FIELDS %iGeneration, %gen0Name 246 DO iq = 1, nqtrue 247 t1 => tracers(iq) 248 CALL addKey('name', t1%name, t1%keys) 249 CALL addKey('component', t1%component, t1%keys) 250 CALL addKey('parent', t1%parent, t1%keys) 251 CALL addKey('phase', t1%phase, t1%keys) 252 END DO 253 IF(setGeneration(tracers)) CALL abort_gcm(modname,'See above',1) !- SET FIELDS %iGeneration, %gen0Name 283 254 DEALLOCATE(had, hadv_inca, vad, vadv_inca, conv_flg_inca, pbl_flg_inca, solsym_inca) 284 255 #endif 285 256 !--------------------------------------------------------------------------------------------------------------------------- 286 ELSE !=== FOUND NEW STYLE TRACERS CONFIGURATION FILE(S)257 ELSE !=== OTHER CASES (OLD OR NEW FORMAT, NO INCA MODULE) 287 258 !--------------------------------------------------------------------------------------------------------------------------- 288 259 nqo = COUNT(delPhase(tracers(:)%name) == 'H2O' & … … 300 271 !--------------------------------------------------------------------------------------------------------------------------- 301 272 302 CALL getKey_init(tracers) 303 273 #ifdef REPROBUS 304 274 !--- Transfert the number of tracers to Reprobus 305 #ifdef REPROBUS306 275 CALL Init_chem_rep_trac(nbtr, nqo, tracers(:)%name) 307 #endif 308 276 277 #endif 309 278 !============================================================================================================================== 310 279 ! 2) Calculate nqtot, number of tracers needed (greater if advection schemes 20 or 30 have been chosen). … … 380 349 CALL MOVE_ALLOC(FROM=ttr, TO=tracers) 381 350 382 !--- SET FIELDS %iqParent, %nqChild s, %iGeneration, %iqDescen, %nqDescen351 !--- SET FIELDS %iqParent, %nqChildren, %iGeneration, %iqDescen, %nqDescen 383 352 CALL indexUpdate(tracers) 384 353 … … 404 373 END DO 405 374 406 niso = 0; nzone=0; nphas=nqo; ntiso = 0; isoCheck=.FALSE. 407 IF(initIsotopes(tracers, isotopes)) CALL abort_gcm(modname, 'Problem when reading isotopes parameters', 1) 408 nbIso = SIZE(isotopes) 409 nqtottr = nqtot - COUNT(tracers%gen0Name == 'H2O' .AND. tracers%component == 'lmdz') 410 IF(nbIso/=0) THEN !--- ISOTOPES FOUND 411 412 !=== READ PHYSICAL PARAMETERS FROM "isotopes_params.def" FILE SPECIFIC TO WATER ISOTOPES 413 ! DONE HERE, AND NOT ONLY IN "infotrac_phy", BECAUSE SOME PHYSICAL PARAMS ARE NEEDED FOR RESTARTS (tnat, alpha_ideal) 414 CALL getKey_init(tracers, isotopes) 415 IF(isoSelect('H2O')) RETURN !--- Select water isotopes ; finished if no water isotopes 416 iH2O = ixIso !--- Keep track of water family index 417 IF(getKey('tnat' , tnat, isoName(1:niso))) CALL abort_gcm(modname, 'can''t read "tnat"', 1) 418 IF(getKey('alpha', alpha_ideal, isoName(1:niso))) CALL abort_gcm(modname, 'can''t read "alpha_ideal"', 1) 419 420 !=== MAKE SURE THE MEMBERS OF AN ISOTOPES FAMILY ARE PRESENT IN THE SAME PHASES 421 DO ix = 1, nbIso 422 iso => isotopes(ix) 423 !--- Check whether each isotope and tagging isotopic tracer is present in the same number of phases 424 DO it = 1, iso%ntiso 425 np = SUM([(COUNT(tracers(:)%name == addPhase(iso%trac(it), iso%phase(ip:ip))), ip=1, iso%nphas)]) 426 IF(np == iso%nphas) CYCLE 427 WRITE(msg1,'("Found ",i0," phases for ",a," instead of ",i0)')np, TRIM(iso%trac(it)), iso%nphas 428 CALL abort_gcm(modname, msg1, 1) 429 END DO 430 DO it = 1, iso%niso 431 nz = SUM([(COUNT(iso%trac == TRIM(iso%trac(it))//'_'//iso%zone(iz)), iz=1, iso%nzone)]) 432 IF(nz == iso%nzone) CYCLE 433 WRITE(msg1,'("Found ",i0," tagging zones for ",a," instead of ",i0)')nz, TRIM(iso%trac(it)), iso%nzone 434 CALL abort_gcm(modname, msg1, 1) 435 END DO 436 END DO 437 END IF 375 !=== READ PHYSICAL PARAMETERS FOR ISOTOPES ; DONE HERE BECAUSE dynetat0 AND iniacademic NEED "tnat" AND "alpha_ideal" 376 niso = 0; nzone = 0; nphas = nqo; ntiso = 0; isoCheck = .FALSE. 377 IF(readIsotopesFile()) CALL abort_gcm(modname, 'Problem when reading isotopes parameters', 1) 438 378 439 379 !--- Convection / boundary layer activation for all tracers … … 442 382 443 383 !--- Note: nqtottr can differ from nbtr when nmom/=0 444 ! IF(COUNT(tracers%iso_iName == 0) - COUNT(delPhase(tracers%name) == 'H2O' .AND. tracers%component == 'lmdz') /= nqtottr) & 445 ! CALL abort_gcm('infotrac_init', 'pb dans le calcul de nqtottr', 1) 384 nqtottr = nqtot - COUNT(delPhase(tracers%gen0Name) == 'H2O' .AND. tracers%component == 'lmdz') 385 IF(COUNT(tracers%iso_iName == 0) - COUNT(delPhase(tracers%name) == 'H2O' .AND. tracers%component == 'lmdz') /= nqtottr) & 386 CALL abort_gcm(modname, 'pb dans le calcul de nqtottr', 1) 446 387 447 388 !=== DISPLAY THE RESULTS … … 459 400 CALL msg('Information stored in infotrac :', modname) 460 401 IF(dispTable('isssssssssiiiiiiiii', & 461 ['iq ', 'name ', 'lName ', 'gen0N ', 'parent', 'type ', 'phase ', 'compon', 'is Adv ', 'isPhy', &402 ['iq ', 'name ', 'lName ', 'gen0N ', 'parent', 'type ', 'phase ', 'compon', 'isPhy ', 'isAdv ', & 462 403 'iadv ', 'iGen ', 'iqPar ', 'nqDes ', 'nqChld', 'iGroup', 'iName ', 'iZone ', 'iPhase'], & 463 cat(t%name, t%longName, t%gen0Name, t%parent, t%type, t%phase, t%component, bool2str(t%is Advected),&464 bool2str(t%is InPhysics)),&465 cat([(iq, iq=1, nqtot)], t%iadv, t%iGeneration, t%iqParent, t%nqDescen, t%nqChild s, t%iso_iGroup,&404 cat(t%name, t%longName, t%gen0Name, t%parent, t%type, t%phase, t%component, bool2str(t%isInPhysics), & 405 bool2str(t%isAdvected)), & 406 cat([(iq, iq=1, nqtot)], t%iadv, t%iGeneration, t%iqParent, t%nqDescen, t%nqChildren, t%iso_iGroup, & 466 407 t%iso_iName, t%iso_iZone, t%iso_iPhase), nColMax=maxTableWidth, nHead=2, sub=modname)) & 467 408 CALL abort_gcm(modname, "problem with the tracers table content", 1) 468 409 IF(niso > 0) THEN 469 410 CALL msg('Where, for isotopes family "'//TRIM(isotope%parent)//'":', modname) 470 CALL msg(' isoKeys = '//strStack(isoKeys%name), modname)411 CALL msg(' isoKeys%name = '//strStack(isoKeys%name), modname) 471 412 CALL msg(' isoName = '//strStack(isoName), modname) 472 413 CALL msg(' isoZone = '//strStack(isoZone), modname) … … 477 418 CALL msg('end', modname) 478 419 479 END SUBROUTINE infotrac_init 480 481 482 !============================================================================================================================== 483 !=== THE ROUTINE isoSelect IS USED TO SWITCH FROM AN ISOTOPE FAMILY TO ANOTHER: ISOTOPES DEPENDENT PARAMETERS ARE UPDATED 484 ! Single generic "isoSelect" routine, using the predefined index of the parent (fast version) or its name (first call). 485 !============================================================================================================================== 486 LOGICAL FUNCTION isoSelectByName(iName, lVerbose) RESULT(lerr) 487 IMPLICIT NONE 488 CHARACTER(LEN=*), INTENT(IN) :: iName 489 LOGICAL, OPTIONAL, INTENT(IN) :: lVerbose 490 INTEGER :: iIso 491 LOGICAL :: lV 492 lV = .FALSE.; IF(PRESENT(lVerbose)) lV = lVerbose 493 iIso = strIdx(isotopes(:)%parent, iName) 494 lerr = iIso == 0 495 IF(lerr) THEN 496 niso = 0; ntiso = 0; nzone=0; nphas=nqo; isoCheck=.FALSE. 497 CALL msg('no isotope family named "'//TRIM(iName)//'"', ll=lV) 498 RETURN 499 END IF 500 lerr = isoSelectByIndex(iIso, lV) 501 END FUNCTION isoSelectByName 502 !============================================================================================================================== 503 LOGICAL FUNCTION isoSelectByIndex(iIso, lVerbose) RESULT(lerr) 504 IMPLICIT NONE 505 INTEGER, INTENT(IN) :: iIso 506 LOGICAL, OPTIONAL, INTENT(IN) :: lVerbose 507 LOGICAL :: lv 508 lv = .FALSE.; IF(PRESENT(lVerbose)) lv = lVerbose 509 lerr = .FALSE. 510 IF(iIso == ixIso) RETURN !--- Nothing to do if the index is already OK 511 lerr = iIso<=0 .OR. iIso>nbIso 512 CALL msg('Inconsistent isotopes family index '//TRIM(int2str(iIso))//': should be > 0 and <= '//TRIM(int2str(nbIso))//'"',& 513 ll=lerr .AND. lV) 514 IF(lerr) RETURN 515 ixIso = iIso !--- Update currently selected family index 516 isotope => isotopes(ixIso) !--- Select corresponding component 517 isoKeys => isotope%keys; niso = isotope%niso 518 isoName => isotope%trac; ntiso = isotope%ntiso 519 isoZone => isotope%zone; nzone = isotope%nzone 520 isoPhas => isotope%phase; nphas = isotope%nphas 521 itZonIso => isotope%itZonIso; isoCheck = isotope%check 522 iqIsoPha => isotope%iqIsoPha 523 END FUNCTION isoSelectByIndex 524 !============================================================================================================================== 420 END SUBROUTINE init_infotrac 525 421 526 422 END MODULE infotrac -
LMDZ6/branches/LMDZ_ECRad/libf/dyn3d_common/iso_verif_dyn.F
r4050 r4482 64 64 function iso_verif_aberrant_nostop 65 65 : (x,iso,q,err_msg) 66 USE infotrac, ONLY: tnat66 USE infotrac, ONLY: isoName, getKey 67 67 implicit none 68 68 … … 74 74 ! locals 75 75 real qmin,deltaD 76 real deltaDmax,deltaDmin 76 real deltaDmax,deltaDmin,tnat 77 77 parameter (qmin=1e-11) 78 78 parameter (deltaDmax=200.0,deltaDmin=-999.9) … … 85 85 ! verifier que HDO est raisonable 86 86 if (q.gt.qmin) then 87 deltaD=(x/q/tnat(iso)-1)*1000 87 IF(getKey('tnat', tnat, isoName(iso))) THEN 88 err_msg = 'Missing isotopic parameter "tnat"' 89 iso_verif_aberrant_nostop=1 90 RETURN 91 END IF 92 deltaD=(x/q/tnat-1)*1000 88 93 if ((deltaD.gt.deltaDmax).or.(deltaD.lt.deltaDmin)) then 89 94 write(*,*) 'erreur detectee par iso_verif_aberrant:'
Note: See TracChangeset
for help on using the changeset viewer.