Ignore:
Timestamp:
Mar 29, 2023, 3:14:27 PM (2 years ago)
Author:
lguez
Message:

Sync latest trunk changes to branch LMDZ_ECRad

Location:
LMDZ6/branches/LMDZ_ECRad
Files:
8 edited

Legend:

Unmodified
Added
Removed
  • LMDZ6/branches/LMDZ_ECRad

  • LMDZ6/branches/LMDZ_ECRad/libf/dyn3d_common/comvert_mod.F90

    r2602 r4482  
    1010
    1111PUBLIC :: 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
    1314
    1415REAL ap(llm+1) ! hybrid pressure contribution at interlayers
    1516REAL bp (llm+1) ! hybrid sigma contribution at interlayer
    1617REAL presnivs(llm) ! (reference) pressure at mid-layers
     18REAL presinter(llm+1) ! (reference) pressure at interlayers
    1719REAL dpres(llm)
    1820REAL sig(llm+1)
  • LMDZ6/branches/LMDZ_ECRad/libf/dyn3d_common/conf_planete.F90

    r2600 r4482  
    3030
    3131!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."
     35preff=101080.
    3336CALL getin('preff', preff)
     37
    3438! Reference pressure at which hybrid coord. become purely pressure
    3539! pa=50000.
    3640pa=preff/2.
    3741CALL getin('pa', pa)
     42
    3843! Gravity
    3944g=9.80665
     45
    4046CALL getin('g',g)
    4147! Molar mass of the atmosphere
     48
    4249molmass = 28.9644
    4350CALL getin('molmass',molmass)
    4451! kappa=R/Cp et Cp     
     52
    4553kappa = 2./7.
    4654CALL getin('kappa',kappa)
     55
    4756cpp=8.3145/molmass/kappa*1000.
    4857CALL getin('cpp',cpp)
    4958! Radius of the planet
     59
    5060rad = 6371229.
    5161CALL getin('radius',rad)
     62
    5263! Length of a standard day (s)
    5364daysec=86400.
    5465CALL getin('daysec',daysec)
     66
    5567! Rotation rate of the planet:
    5668! Length of a solar day, in standard days
    5769daylen = 1.
     70
    5871CALL getin('daylen',daylen)
    5972! Number of days (standard) per year:
     73
    6074year_day = 365.25
    6175CALL getin('year_day',year_day)
    6276! Omega
    6377! omeg=2.*pi/86400.
     78
    6479omeg=2.*pi/daysec*(1./daylen+1./year_day)
    6580CALL getin('omeg',omeg)
  • LMDZ6/branches/LMDZ_ECRad/libf/dyn3d_common/control_mod.F90

    r4146 r4482  
    2929  INTEGER,SAVE :: ip_ebil_dyn
    3030  LOGICAL,SAVE :: offline
    31   CHARACTER(len=4),SAVE :: config_inca
    3231  CHARACTER(len=10),SAVE :: planet_type ! planet type ('earth','mars',...)
    3332  LOGICAL,SAVE :: output_grads_dyn ! output dynamics diagnostics in
  • LMDZ6/branches/LMDZ_ECRad/libf/dyn3d_common/disvert.F90

    r2786 r4482  
    1111  use assert_m, only: assert
    1212  USE comvert_mod, ONLY: ap, bp, aps, bps, nivsigs, nivsig, dpres, presnivs, &
    13                          pseudoalt, pa, preff, scaleheight
     13                         pseudoalt, pa, preff, scaleheight, presinter
    1414  USE logic_mod, ONLY: ok_strato
    1515
     
    3535! dpres(llm)                 !--- PRESSURE DIFFERENCE FOR EACH LAYER
    3636! presnivs(llm)              !--- PRESSURE AT EACH MID-LAYER
     37! presinter(llm+1)           !--- PRESSURE AT EACH INTERFACE
    3738! scaleheight                !--- VERTICAL SCALE HEIGHT            (Earth: 8kms)
    3839! nivsig(llm+1)              !--- SIGMA INDEX OF EACH LAYER INTERFACE
     
    355356          max(ap(l+1)+bp(l+1)*preff, 1.e-10))
    356357  ENDDO
     358  DO l=1, llmp1
     359     presinter(l)= ( ap(l)+bp(l)*preff)
     360     write(lunout, *)'PRESINTER(', l, ')=', presinter(l)
     361  ENDDO
    357362
    358363  write(lunout, *) trim(modname),': PRESNIVS '
     
    395400      ELSE IF(SIGN(f1,f4)/=f1) THEN; x2=xx; f2=f4
    396401      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)
    398403      END IF
    399404      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  
    22! $Id: $
    33!
    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
    89!
     10! The subroutine is called in dynphy_lonlat/phylmd/ce0l.F90.
     11
    912SUBROUTINE grilles_gcm_netcdf_sub(masque,phis)
    1013
    1114  USE comconst_mod, ONLY: cpp, kappa, g, omeg, daysec, rad, pi
    1215  USE comvert_mod, ONLY: presnivs, preff, pa
     16  use netcdf, only: nf90_def_var, nf90_int, nf90_float, nf90_put_var
    1317 
    1418  IMPLICIT NONE
     
    1923  INCLUDE "netcdf.inc"
    2024
    21 
     25!========================
    2226  REAL,DIMENSION(iip1,jjp1),INTENT(IN)  :: masque ! masque terre/mer
    2327  REAL,DIMENSION(iip1,jjp1),INTENT(IN)  :: phis   ! geopotentiel au sol
    2428
    25   REAL temp(iim+1,jjm+1)
    26   ! Attributs netcdf sortie
     29  INTEGER status,i,j
     30 
     31  ! Attributs netcdf output
    2732  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
    3036  INTEGER out_lonudim,out_lonvdim
    31   INTEGER out_latudim,out_latvdim,out_dim(3)
     37  INTEGER out_latudim,out_latvdim,out_dim(2)
    3238  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)
    3847  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  !
    4350  INTEGER masque_int(iip1,jjp1)
    44   INTEGER :: phis_id
    45   INTEGER :: area_id
    46   INTEGER :: mask_id
    47   INTEGER :: presnivs_id
    48  
     51  REAL :: phis_loc(iip1,jjp1)
     52 
     53  !========================
     54  ! CALCULATION of latu, latv, lonu, lonv in deg.
     55  ! ---------------------------------------------------
    4956  rad = 6400000
    5057  omeg = 7.272205e-05
     
    6471     rlatudeg(j)=rlatu(j)*180./pi
    6572  ENDDO
     73 
    6674  DO j=1,jjm
    6775     rlatvdeg(j)=rlatv(j)*180./pi
     
    7280     rlonvdeg(i)=rlonv(i)*180./pi + 360.
    7381  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  !-------------------------
    80104  status=NF_CREATE('grilles_gcm.nc',IOR(NF_CLOBBER,NF_64BIT_OFFSET),ncid_out)
    81105  CALL handle_err(status)
     106 
     107  ! DEFINE output dimensions
     108  !-------------------------
    82109  status=NF_DEF_DIM(ncid_out,'lonu',iim+1,out_lonudim)
    83110  CALL handle_err(status)
     
    88115  status=NF_DEF_DIM(ncid_out,'latv',jjm,out_latvdim)
    89116  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)
    94125  CALL handle_err(status)
    95126  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)
    100130  CALL handle_err(status)
    101131  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)
    106135  CALL handle_err(status)
    107136  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)
    112140  CALL handle_err(status)
    113141  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
    117144  out_dim(1)=out_lonudim
    118145  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
    125151  out_dim(1)=out_lonvdim
    126152  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
    133158  out_dim(1)=out_lonvdim
    134159  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
    208172  out_dim(1)=out_lonvdim
    209173  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
    234208  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
    248227  CALL ncclos(ncid_out,rcode_out)
     228  write(*,*) "END grilles_gcm_netcdf_sub OK"
    249229
    250230END SUBROUTINE grilles_gcm_netcdf_sub
    251 
    252231
    253232
  • LMDZ6/branches/LMDZ_ECRad/libf/dyn3d_common/infotrac.F90

    r4203 r4482  
    33MODULE infotrac
    44
    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
    1110   IMPLICIT NONE
    1211
     
    1413
    1514   !=== FOR TRACERS:
    16    PUBLIC :: infotrac_init                                 !--- Initialization of the tracers
    17    PUBLIC :: tracers, type_trac, types_trac                !--- Full tracers database, tracers type keyword
     15   PUBLIC :: init_infotrac                                 !--- Initialization of the tracers
     16   PUBLIC :: tracers, type_trac                            !--- Full tracers database, tracers type keyword
    1817   PUBLIC :: nqtot,   nbtr,   nqo,   nqCO2,   nqtottr      !--- Main dimensions
    1918   PUBLIC :: conv_flg, pbl_flg                             !--- Convection & boundary layer activation keys
    2019
    2120   !=== FOR ISOTOPES: General
    22    PUBLIC :: isotopes, nbIso                              !--- Derived type, full isotopes families database + nb of families
     21   PUBLIC :: isot_type, nbIso                              !--- Derived type, full isotopes families database + nb of families
    2322   PUBLIC :: isoSelect, ixIso                              !--- Isotopes family selection tool + selected family index
    2423   !=== 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
    2625   PUBLIC :: min_qParent, min_qMass, min_ratio             !--- Min. values for various isotopic quantities
    2726   !=== FOR ISOTOPES: Depending on the selected isotopes family
     
    3433   !=== FOR BOTH TRACERS AND ISOTOPES
    3534   PUBLIC :: getKey                                        !--- Get a key from "tracers" or "isotope"
    36 
    37    INTERFACE isoSelect; MODULE PROCEDURE isoSelectByIndex, isoSelectByName; END INTERFACE isoSelect
    3835
    3936!=== CONVENTIONS FOR TRACERS NUMBERS:
     
    7168!  | phase       | Phases list ("g"as / "l"iquid / "s"olid)             | /           | [g][l][s]              |
    7269!  | 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 |
    7470!  | 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  |
    7771!  | iqParent    | Index of the parent tracer                           | iqpere      | 1:nqtot                |
    7872!  | iqDescen    | Indexes of the childs       (all generations)        | iqfils      | 1:nqtot                |
    7973!  | 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  |
    8179!  | iso_iGroup  | Isotopes group index in isotopes(:)                  | /           | 1:nbIso                |
    8280!  | iso_iName   | Isotope  name  index in isotopes(iso_iGroup)%trac(:) | iso_indnum  | 1:niso                 |
    8381!  | iso_iZone   | Isotope  zone  index in isotopes(iso_iGroup)%zone(:) | zone_num    | 1:nzone                |
    8482!  | iso_iPhas   | Isotope  phase index in isotopes(iso_iGroup)%phas(:) | phase_num   | 1:nphas                |
    85 !  | keys        | key/val pairs accessible with "getKey" routine       | /           |                        |
    8683!  +-------------+------------------------------------------------------+-------------+------------------------+
    8784!
     
    103100
    104101   !=== 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)
    136112
    137113CONTAINS
    138114
    139 SUBROUTINE infotrac_init
    140    USE control_mod, ONLY: planet_type, config_inca
     115SUBROUTINE init_infotrac
     116   USE control_mod, ONLY: planet_type
    141117#ifdef REPROBUS
    142118   USE CHEM_REP,    ONLY: Init_chem_rep_trac
     
    176152   CHARACTER(LEN=2)      ::   suff(9)                                !--- Suffixes for schemes of order 3 or 4 (Prather)
    177153   CHARACTER(LEN=3)      :: descrq(30)                               !--- Advection scheme description tags
    178    CHARACTER(LEN=maxlen) :: msg1                                     !--- String for messages
     154   CHARACTER(LEN=maxlen) :: msg1, texp, ttp                          !--- Strings for messages and expanded tracers type
    179155   INTEGER :: fType                                                  !--- Tracers description file type ; 0: none
    180156                                                                     !--- 1/2/3: "traceur.def"/"tracer.def"/"tracer_*.def"
    181157   INTEGER :: nqtrue                                                 !--- Tracers nb from tracer.def (no higher order moments)
    182158   INTEGER :: iad                                                    !--- Advection scheme number
    183    INTEGER :: ic, ip, np, iq, jq, it, nt, im, nm, ix, iz, nz, k      !--- Indexes and temporary variables
     159   INTEGER :: ic, iq, jq, it, nt, im, nm, iz, k                      !--- Indexes and temporary variables
    184160   LOGICAL :: lerr, ll
    185161   CHARACTER(LEN=1) :: p
    186162   TYPE(trac_type), ALLOCATABLE, TARGET :: ttr(:)
    187163   TYPE(trac_type), POINTER             :: t1, t(:)
    188    TYPE(isot_type), POINTER             :: iso
    189 
    190    CHARACTER(LEN=maxlen), ALLOCATABLE :: tnom_0(:), tnom_transp(:)   !--- Tracer short name + transporting fluid name
    191    CHARACTER(LEN=maxlen)              :: tchaine
    192164   INTEGER :: ierr
    193165
    194    CHARACTER(LEN=*), PARAMETER :: modname="infotrac_init"
     166   CHARACTER(LEN=*), PARAMETER :: modname="init_infotrac"
    195167!------------------------------------------------------------------------------------------------------------------------------
    196168! Initialization :
     
    202174   
    203175   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')
    228192#ifndef INCA
    229             CALL abort_gcm(modname, 'You must add cpp key INCA and compile with INCA code', 1)
    230 #endif
    231          CASE('repr')
     193         CALL abort_gcm(modname, 'You must add cpp key INCA and compile with INCA code', 1)
     194#endif
     195      CASE('repr')
    232196#ifndef REPROBUS
    233             CALL abort_gcm(modname, 'You must add cpp key REPROBUS and compile with REPROBUS code', 1)
    234 #endif
    235          CASE('coag')
     197         CALL abort_gcm(modname, 'You must add cpp key REPROBUS and compile with REPROBUS code', 1)
     198#endif
     199      CASE('coag')
    236200#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'] )
    250206
    251207!==============================================================================================================================
    252208! 1) Get the numbers of: true (first order only) tracers "nqtrue", water tracers "nqo" (vapor/liquid/solid)
    253209!==============================================================================================================================
    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   !---------------------------------------------------------------------------------------------------------------------------
    255220   IF(fType == 0) CALL abort_gcm(modname, 'Missing tracers file: "traceur.def", "tracer.def" or "tracer_<keyword>.def file.',1)
    256221   !---------------------------------------------------------------------------------------------------------------------------
    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"
    258223   !---------------------------------------------------------------------------------------------------------------------------
    259224#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
    262226      CALL Init_chem_inca_trac(nqINCA)                               !--- Get nqINCA from INCA
    263227      nbtr = nqINCA + nqCO2                                          !--- Number of tracers passed to phytrac
     
    268232      CALL init_transport(solsym_inca, conv_flg_inca, pbl_flg_inca, hadv_inca, vadv_inca)
    269233      ALLOCATE(ttr(nqtrue))
    270       ttr(1:nqo+nqCO2)                    = tracers
    271       ttr(1    :      nqo   )%component   = 'lmdz'
    272       ttr(1+nqo:nqCO2+nqo   )%component   = 'co2i'
    273       ttr(1+nqo+nqCO2:nqtrue)%component   = 'inca'
    274       ttr(1+nqo      :nqtrue)%name        = [('CO2     ', k=1, nqCO2), solsym_inca]
    275       ttr(1+nqo+nqCO2:nqtrue)%parent      = tran0
    276       ttr(1+nqo+nqCO2:nqtrue)%phase       = 'g'
     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'
    277241      lerr = getKey('hadv', had, ky=tracers(:)%keys)
    278242      lerr = getKey('vadv', vad, ky=tracers(:)%keys)
    279       hadv(1:nqo) = had(:); hadv(nqo+1:nqtrue) = hadv_inca
    280       vadv(1:nqo) = vad(:); vadv(nqo+1:nqtrue) = vadv_inca
     243      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
    281245      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
    283254      DEALLOCATE(had, hadv_inca, vad, vadv_inca, conv_flg_inca, pbl_flg_inca, solsym_inca)
    284255#endif
    285256   !---------------------------------------------------------------------------------------------------------------------------
    286    ELSE                                                              !=== FOUND NEW STYLE TRACERS CONFIGURATION FILE(S)
     257   ELSE                                                              !=== OTHER CASES (OLD OR NEW FORMAT, NO INCA MODULE)
    287258   !---------------------------------------------------------------------------------------------------------------------------
    288259      nqo    =        COUNT(delPhase(tracers(:)%name)     == 'H2O' &
     
    300271   !---------------------------------------------------------------------------------------------------------------------------
    301272
    302    CALL getKey_init(tracers)
    303 
     273#ifdef REPROBUS
    304274   !--- Transfert the number of tracers to Reprobus
    305 #ifdef REPROBUS
    306275   CALL Init_chem_rep_trac(nbtr, nqo, tracers(:)%name)
    307 #endif
    308 
     276
     277#endif
    309278!==============================================================================================================================
    310279! 2) Calculate nqtot, number of tracers needed (greater if advection schemes 20 or 30 have been chosen).
     
    380349   CALL MOVE_ALLOC(FROM=ttr, TO=tracers)
    381350
    382    !--- SET FIELDS %iqParent, %nqChilds, %iGeneration, %iqDescen, %nqDescen
     351   !--- SET FIELDS %iqParent, %nqChildren, %iGeneration, %iqDescen, %nqDescen
    383352   CALL indexUpdate(tracers)
    384353
     
    404373   END DO
    405374
    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)
    438378
    439379   !--- Convection / boundary layer activation for all tracers
     
    442382
    443383   !--- 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)
    446387
    447388   !=== DISPLAY THE RESULTS
     
    459400   CALL msg('Information stored in infotrac :', modname)
    460401   IF(dispTable('isssssssssiiiiiiiii', &
    461       ['iq    ', 'name  ', 'lName ', 'gen0N ', 'parent', 'type  ', 'phase ', 'compon', 'isAdv ', 'isPhy ', &
     402      ['iq    ', 'name  ', 'lName ', 'gen0N ', 'parent', 'type  ', 'phase ', 'compon', 'isPhy ', 'isAdv ', &
    462403       '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%isAdvected), &
    464                                                                                   bool2str(t%isInPhysics)),&
    465       cat([(iq, iq=1, nqtot)], t%iadv, t%iGeneration, t%iqParent, t%nqDescen, t%nqChilds, 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,  &
    466407                  t%iso_iName, t%iso_iZone, t%iso_iPhase), nColMax=maxTableWidth, nHead=2, sub=modname))   &
    467408      CALL abort_gcm(modname, "problem with the tracers table content", 1)
    468409   IF(niso > 0) THEN
    469410      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)
    471412      CALL msg('  isoName = '//strStack(isoName),      modname)
    472413      CALL msg('  isoZone = '//strStack(isoZone),      modname)
     
    477418   CALL msg('end', modname)
    478419
    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 !==============================================================================================================================
     420END SUBROUTINE init_infotrac
    525421
    526422END MODULE infotrac
  • LMDZ6/branches/LMDZ_ECRad/libf/dyn3d_common/iso_verif_dyn.F

    r4050 r4482  
    6464        function iso_verif_aberrant_nostop
    6565     :           (x,iso,q,err_msg)
    66         USE infotrac, ONLY: tnat
     66        USE infotrac, ONLY: isoName, getKey
    6767        implicit none
    6868       
     
    7474        ! locals
    7575        real qmin,deltaD
    76         real deltaDmax,deltaDmin
     76        real deltaDmax,deltaDmin,tnat
    7777        parameter (qmin=1e-11)
    7878        parameter (deltaDmax=200.0,deltaDmin=-999.9)
     
    8585        ! verifier que HDO est raisonable
    8686         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
    8893             if ((deltaD.gt.deltaDmax).or.(deltaD.lt.deltaDmin)) then
    8994                  write(*,*) 'erreur detectee par iso_verif_aberrant:'
Note: See TracChangeset for help on using the changeset viewer.