source: LMDZ6/trunk/libf/dyn3d_common/infotrac.F90 @ 4301

Last change on this file since 4301 was 4301, checked in by dcugnet, 20 months ago
  • for REPROBUS: simplify (and fix) the handling of exceptions for old HNO3 convention.
  • cleaning + comments added in readTracFiles_mod.
  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
  • Property svn:keywords set to Id
File size: 34.6 KB
Line 
1!$Id: infotrac.F90 4301 2022-10-20 11:57:21Z dcugnet $
2!
3MODULE infotrac
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, indexUpdate, addPhase, getKey, maxTableWidth, keys_type, &
8                                isot_type, setGeneration,   initIsotopes, delPhase, getKey_init, ancestor, tran0
9                               
10   IMPLICIT NONE
11
12   PRIVATE
13
14   !=== FOR TRACERS:
15   PUBLIC :: infotrac_init                                 !--- Initialization of the tracers
16   PUBLIC :: tracers, type_trac, types_trac                !--- Full tracers database, tracers type keyword
17   PUBLIC :: nqtot,   nbtr,   nqo,   nqCO2,   nqtottr      !--- Main dimensions
18   PUBLIC :: conv_flg, pbl_flg                             !--- Convection & boundary layer activation keys
19
20   !=== FOR ISOTOPES: General
21   PUBLIC :: isotopes,  nbIso                              !--- Derived type, full isotopes families database + nb of families
22   PUBLIC :: isoSelect, ixIso                              !--- Isotopes family selection tool + selected family index
23   !=== FOR ISOTOPES: Specific to water
24   PUBLIC :: iH2O, tnat, alpha_ideal                       !--- H2O isotopes index, natural abundance, fractionning coeff.
25   PUBLIC :: min_qParent, min_qMass, min_ratio             !--- Min. values for various isotopic quantities
26   !=== FOR ISOTOPES: Depending on the selected isotopes family
27   PUBLIC :: isotope, isoKeys                              !--- Selected isotopes database + associated keys (cf. getKey)
28   PUBLIC :: isoName, isoZone, isoPhas                     !--- Isotopes and tagging zones names, phases
29   PUBLIC :: niso,    nzone,   nphas,   ntiso              !---  " " numbers + isotopes & tagging tracers number
30   PUBLIC :: itZonIso                                      !--- idx "it" (in "isoName(1:niso)") = function(tagging idx, isotope idx)
31   PUBLIC :: iqIsoPha                                      !--- idx "iq" (in "qx") = function(isotope idx, phase idx) + aliases
32   PUBLIC :: isoCheck                                      !--- Run isotopes checking routines
33   !=== FOR BOTH TRACERS AND ISOTOPES
34   PUBLIC :: getKey                                        !--- Get a key from "tracers" or "isotope"
35
36   INTERFACE isoSelect; MODULE PROCEDURE isoSelectByIndex, isoSelectByName; END INTERFACE isoSelect
37
38!=== CONVENTIONS FOR TRACERS NUMBERS:
39!  |--------------------+-----------------------+-----------------+---------------+----------------------------|
40!  | water in different |    water tagging      |  water isotopes | other tracers | additional tracers moments |
41!  | phases: H2O_[gls]  |      isotopes         |                 |               |  for higher order schemes  |
42!  |--------------------+-----------------------+-----------------+---------------+----------------------------|
43!  |                    |                       |                 |               |                            |
44!  |<--     nqo      -->|<-- nqo*niso* nzone -->|<-- nqo*niso  -->|<--  nbtr   -->|<--        (nmom)        -->|         
45!  |                    |                                         |                                            |
46!  |                    |<-- nqo*niso*(nzone+1)  =   nqo*ntiso -->|<--    nqtottr = nbtr + nmom             -->|
47!  |                                                                              = nqtot - nqo*(ntiso+1)      |
48!  |                                                                                                           |
49!  |<--                        nqtrue  =  nbtr + nqo*(ntiso+1)                 -->|                            |
50!  |                                                                                                           |
51!  |<--                        nqtot   =  nqtrue + nmom                                                     -->|
52!  |                                                                                                           |
53!  |-----------------------------------------------------------------------------------------------------------|
54!  NOTES FOR THIS TABLE:
55!  * Used "niso", "nzone" and "ntiso" are components of "isotopes(ip)" for water (isotopes(ip)%parent == 'H2O'),
56!    since water is so far the sole tracers family, except passive CO2, removed from the main tracers table.
57!  * For water, "nqo" is equal to the more general field "isotopes(ip)%nphas".
58!  * "niso", "nzone", "ntiso", "nphas" are defined for other isotopic tracers families, if any.
59!
60!=== DERIVED TYPE EMBEDDING MOST OF THE TRACERS-RELATED QUANTITIES (LENGTH: nqtot)
61!    Each entry is accessible using "%" sign.
62!  |-------------+------------------------------------------------------+-------------+------------------------+
63!  |  entry      | Meaning                                              | Former name | Possible values        |
64!  |-------------+------------------------------------------------------+-------------+------------------------+
65!  | name        | Name (short)                                         | tname       |                        |
66!  | gen0Name    | Name of the 1st generation ancestor                  | /           |                        |
67!  | parent      | Name of the parent                                   | /           |                        |
68!  | longName    | Long name (with adv. scheme suffix) for outputs      | ttext       |                        |
69!  | type        | Type (so far: tracer or tag)                         | /           | tracer,tag             |
70!  | phase       | Phases list ("g"as / "l"iquid / "s"olid)             | /           | [g][l][s]              |
71!  | component   | Name(s) of the merged/cumulated section(s)           | /           | coma-separated names   |
72!  | iadv        | Advection scheme number                              | iadv        | 1-20,30 exc. 3-9,15,19 |
73!  | iGeneration | Generation (>=1)                                     | /           |                        |
74!  | isAdvected  | advected tracers flag (.TRUE. if iadv >= 0)          | /           | nqtrue  .TRUE. values  |
75!  | isInPhysics | tracers not extracted from the main table in physics | /           | nqtottr .TRUE. values  |
76!  | iqParent    | Index of the parent tracer                           | iqpere      | 1:nqtot                |
77!  | iqDescen    | Indexes of the childs       (all generations)        | iqfils      | 1:nqtot                |
78!  | nqDescen    | Number of the descendants   (all generations)        | nqdesc      | 1:nqtot                |
79!  | nqChilds    | Number of childs            (1st generation only)    | nqfils      | 1:nqtot                |
80!  | iso_iGroup  | Isotopes group index in isotopes(:)                  | /           | 1:nbIso                |
81!  | iso_iName   | Isotope  name  index in isotopes(iso_iGroup)%trac(:) | iso_indnum  | 1:niso                 |
82!  | iso_iZone   | Isotope  zone  index in isotopes(iso_iGroup)%zone(:) | zone_num    | 1:nzone                |
83!  | iso_iPhas   | Isotope  phase index in isotopes(iso_iGroup)%phas(:) | phase_num   | 1:nphas                |
84!  | keys        | key/val pairs accessible with "getKey" routine       | /           |                        |
85!  +-------------+------------------------------------------------------+-------------+------------------------+
86!
87!=== DERIVED TYPE EMBEDDING MOST OF THE ISOTOPES-RELATED QUANTITIES (LENGTH: nbIso, NUMBER OF ISOTOPES FAMILIES)
88!    Each entry is accessible using "%" sign.
89!  |-----------------+--------------------------------------------------+--------------------+-----------------+
90!  |  entry | length | Meaning                                          |    Former name     | Possible values |
91!  |-----------------+--------------------------------------------------+--------------------+-----------------+
92!  | parent          | Parent tracer (isotopes family name)             |                    |                 |
93!  | keys   | niso   | Isotopes keys/values pairs list + number         |                    |                 |
94!  | trac   | ntiso  | Isotopes + tagging tracers list + number         | / | ntraciso       |                 |
95!  | zone   | nzone  | Geographic tagging zones   list + number         | / | ntraceurs_zone |                 |
96!  | phase  | nphas  | Phases                     list + number         |                    | [g][l][s], 1:3  |
97!  | iqIsoPha        | Index in "qx"           = f(name(1:ntiso)),phas) | iqiso              | 1:nqtot         |
98!  | itZonIso        | Index in "trac(1:ntiso)"= f(zone, name(1:niso))  | index_trac         | 1:ntiso         |
99!  +-----------------+--------------------------------------------------+--------------------+-----------------+
100
101   REAL, PARAMETER :: min_qParent = 1.e-30, min_qMass = 1.e-18, min_ratio = 1.e-16 ! MVals et CRisi
102
103   !=== DIMENSIONS OF THE TRACERS TABLES AND OTHER SCALAR VARIABLES
104   INTEGER,                 SAVE :: nqtot,  &                   !--- Tracers nb in dynamics (incl. higher moments + H2O)
105                                    nbtr,   &                   !--- Tracers nb in physics  (excl. higher moments + H2O)
106                                    nqo,    &                   !--- Number of water phases
107                                    nbIso,  &                   !--- Number of available isotopes family
108                                    nqtottr, &                  !--- Number of tracers passed to phytrac (TO BE DELETED ?)
109                                    nqCO2                       !--- Number of tracers of CO2  (ThL)
110   CHARACTER(LEN=maxlen),   SAVE :: type_trac                   !--- Keyword for tracers type(s)
111   CHARACTER(LEN=maxlen),   SAVE, ALLOCATABLE :: types_trac(:)  !--- Keyword for tracers type(s), parsed version
112
113   !=== DERIVED TYPES EMBEDDING MOST INFORMATIONS ABOUT TRACERS AND ISOTOPES FAMILIES
114   TYPE(trac_type), TARGET, SAVE, ALLOCATABLE ::  tracers(:)    !=== TRACERS DESCRIPTORS VECTOR
115   TYPE(isot_type), TARGET, SAVE, ALLOCATABLE :: isotopes(:)    !=== ISOTOPES PARAMETERS VECTOR
116
117   !=== ALIASES FOR CURRENTLY SELECTED ISOTOPES FAMILY OF VARIABLES EMBEDDED IN THE VECTOR "isotopes"
118   TYPE(isot_type),         SAVE, POINTER :: isotope            !--- CURRENTLY SELECTED ISOTOPES FAMILY DESCRIPTOR
119   INTEGER,                 SAVE          :: ixIso, iH2O        !--- Index of the selected isotopes family and H2O family
120   LOGICAL,                 SAVE          :: isoCheck           !--- Flag to trigger the checking routines
121   TYPE(keys_type),         SAVE, POINTER :: isoKeys(:)         !--- ONE SET OF KEYS FOR EACH ISOTOPE (LISTED IN isoName)
122   CHARACTER(LEN=maxlen),   SAVE, POINTER :: isoName(:),   &    !--- ISOTOPES NAMES FOR THE CURRENTLY SELECTED FAMILY
123                                             isoZone(:),   &    !--- TAGGING ZONES  FOR THE CURRENTLY SELECTED FAMILY
124                                             isoPhas            !--- USED PHASES    FOR THE CURRENTLY SELECTED FAMILY
125   INTEGER,                 SAVE          ::  niso, nzone, &    !--- NUMBER OF ISOTOPES, TAGGING ZONES AND PHASES
126                                             nphas, ntiso       !--- NUMBER OF PHASES AND ISOTOPES + ISOTOPIC TAGGING TRACERS
127   INTEGER,                 SAVE, POINTER ::itZonIso(:,:), &    !--- INDEX IN "isoTrac" AS f(tagging zone idx,  isotope idx)
128                                            iqIsoPha(:,:)       !--- INDEX IN "qx"      AS f(isotopic tracer idx, phase idx)
129
130   !=== VARIABLES FOR ISOTOPES INITIALIZATION AND FOR INCA
131   REAL,                SAVE, ALLOCATABLE ::     tnat(:), &     !--- Natural relative abundance of water isotope        (niso)
132                                          alpha_ideal(:)        !--- Ideal fractionning coefficient (for initial state) (niso)
133   INTEGER,             SAVE, ALLOCATABLE :: conv_flg(:), &     !--- Convection     activation ; needed for INCA        (nbtr)
134                                              pbl_flg(:)        !--- Boundary layer activation ; needed for INCA        (nbtr)
135
136CONTAINS
137
138SUBROUTINE infotrac_init
139   USE control_mod, ONLY: planet_type, config_inca
140#ifdef REPROBUS
141   USE CHEM_REP,    ONLY: Init_chem_rep_trac
142#endif
143   IMPLICIT NONE
144!==============================================================================================================================
145!
146!   Auteur:  P. Le Van /L. Fairhead/F.Hourdin
147!   -------
148!
149!   Modifications:
150!   --------------
151!   05/94: F.Forget      Modif special traceur
152!   02/02: M-A Filiberti Lecture de traceur.def
153!   01/22: D. Cugnet     Nouveaux tracer.def et tracer_*.def + encapsulation (types trac_type et isot_type)
154!
155!   Objet:
156!   ------
157!   GCM LMD nouvelle grille
158!
159!==============================================================================================================================
160!   ... modification de l'integration de q ( 26/04/94 ) ....
161!------------------------------------------------------------------------------------------------------------------------------
162! Declarations:
163   INCLUDE "dimensions.h"
164   INCLUDE "iniprint.h"
165
166!------------------------------------------------------------------------------------------------------------------------------
167! Local variables
168   INTEGER, ALLOCATABLE :: hadv(:), vadv(:)                          !--- Horizontal/vertical transport scheme number
169#ifdef INCA
170   INTEGER, ALLOCATABLE :: had (:), hadv_inca(:), conv_flg_inca(:), &!--- Variables specific to INCA
171                           vad (:), vadv_inca(:),  pbl_flg_inca(:)
172   CHARACTER(LEN=8), ALLOCATABLE :: solsym_inca(:)                   !--- Tracers names for INCA
173   INTEGER :: nqINCA
174#endif
175   CHARACTER(LEN=2)      ::   suff(9)                                !--- Suffixes for schemes of order 3 or 4 (Prather)
176   CHARACTER(LEN=3)      :: descrq(30)                               !--- Advection scheme description tags
177   CHARACTER(LEN=maxlen) :: msg1                                     !--- String for messages
178   INTEGER :: fType                                                  !--- Tracers description file type ; 0: none
179                                                                     !--- 1/2/3: "traceur.def"/"tracer.def"/"tracer_*.def"
180   INTEGER :: nqtrue                                                 !--- Tracers nb from tracer.def (no higher order moments)
181   INTEGER :: iad                                                    !--- Advection scheme number
182   INTEGER :: ic, ip, np, iq, jq, it, nt, im, nm, ix, iz, nz, k      !--- Indexes and temporary variables
183   LOGICAL :: lerr, ll, lRepr
184   CHARACTER(LEN=1) :: p
185   TYPE(trac_type), ALLOCATABLE, TARGET :: ttr(:)
186   TYPE(trac_type), POINTER             :: t1, t(:)
187   TYPE(isot_type), POINTER             :: iso
188   INTEGER :: ierr
189
190   CHARACTER(LEN=*), PARAMETER :: modname="infotrac_init"
191!------------------------------------------------------------------------------------------------------------------------------
192! Initialization :
193!------------------------------------------------------------------------------------------------------------------------------
194   suff          = ['x ','y ','z ','xx','xy','xz','yy','yz','zz']
195   descrq( 1: 2) = ['LMV','BAK']
196   descrq(10:20) = ['VL1','VLP','FH1','FH2','VLH','   ','PPM','PPS','PPP','   ','SLP']
197   descrq(30)    =  'PRA'
198   
199   CALL msg('type_trac = "'//TRIM(type_trac)//'"', modname)
200   IF(strParse(type_trac, '|', types_trac, n=nt)) CALL abort_gcm(modname,'can''t parse "type_trac = '//TRIM(type_trac)//'"',1)
201
202   !---------------------------------------------------------------------------------------------------------------------------
203   DO it = 1, nt                                                          !--- nt>1=> "type_trac": coma-separated keywords list
204   !---------------------------------------------------------------------------------------------------------------------------
205      !--- MESSAGE ABOUT THE CHOSEN CONFIGURATION
206      msg1 = 'For type_trac = "'//TRIM(types_trac(it))//'":'
207      SELECT CASE(types_trac(it))
208         CASE('inca'); CALL msg(TRIM(msg1)//' coupling with INCA chemistry model, config_inca='//config_inca, modname)
209         CASE('inco'); CALL msg(TRIM(msg1)//' coupling jointly with INCA and CO2 cycle',  modname)
210         CASE('repr'); CALL msg(TRIM(msg1)//' coupling with REPROBUS chemistry model',    modname)
211         CASE('co2i'); CALL msg(TRIM(msg1)//' you have chosen to run with CO2 cycle',     modname)
212         CASE('coag'); CALL msg(TRIM(msg1)//' tracers are treated for COAGULATION tests', modname)
213         CASE('lmdz'); CALL msg(TRIM(msg1)//' tracers are treated in LMDZ only',          modname)
214         CASE DEFAULT; CALL abort_gcm(modname,'type_trac='//TRIM(types_trac(it))//' not possible yet.',1)
215      END SELECT
216
217      !--- COHERENCE TEST BETWEEN "type_trac" AND "config_inca"
218      IF(ANY(['inca', 'inco'] == types_trac(it)) .AND. ALL(['aero', 'aeNP', 'chem'] /= config_inca)) &
219         CALL abort_gcm(modname, 'Incoherence between type_trac and config_inca. Please modify "run.def"',1)
220
221      !--- COHERENCE TEST BETWEEN "type_trac" AND PREPROCESSING KEYS
222      SELECT CASE(types_trac(it))
223         CASE('inca', 'inco')
224#ifndef INCA
225            CALL abort_gcm(modname, 'You must add cpp key INCA and compile with INCA code', 1)
226#endif
227         CASE('repr')
228#ifndef REPROBUS
229            CALL abort_gcm(modname, 'You must add cpp key REPROBUS and compile with REPROBUS code', 1)
230#endif
231         CASE('coag')
232#ifndef CPP_StratAer
233            CALL abort_gcm(modname, 'You must add cpp key StratAer and compile with StratAer code', 1)
234#endif
235      END SELECT
236
237   !---------------------------------------------------------------------------------------------------------------------------
238   END DO
239   !---------------------------------------------------------------------------------------------------------------------------
240
241   !--- DISABLE "config_inca" OPTION FOR A RUN WITHOUT "INCA" IF IT DIFFERS FROM "none"
242   IF(fmsg('Setting config_inca="none" as you do not couple with INCA model', &
243         modname, ALL(types_trac /= 'inca') .AND. ALL(types_trac /= 'inco') .AND. config_inca /= 'none')) config_inca = 'none'
244
245   nqCO2 = COUNT( [ANY(types_trac == 'inco') .OR. (ANY(types_trac == 'co2i') .AND. ANY(types_trac == 'inca'))] )
246
247!==============================================================================================================================
248! 1) Get the numbers of: true (first order only) tracers "nqtrue", water tracers "nqo" (vapor/liquid/solid)
249!==============================================================================================================================
250   lRepr = ANY(types_trac(:) == 'repr')
251   IF(readTracersFiles(type_trac, fType, tracers, lRepr)) CALL abort_gcm(modname, 'problem with tracers file(s)',1)
252   !---------------------------------------------------------------------------------------------------------------------------
253   IF(fType == 0) CALL abort_gcm(modname, 'Missing tracers file: "traceur.def", "tracer.def" or "tracer_<keyword>.def file.',1)
254   !---------------------------------------------------------------------------------------------------------------------------
255   IF(fType == 1 .AND. ANY(['inca','inco'] == type_trac)) THEN       !=== FOUND OLD STYLE INCA "traceur.def" (single type_trac)
256   !---------------------------------------------------------------------------------------------------------------------------
257#ifdef INCA
258      nqo = SIZE(tracers) - nqCO2
259      CALL Init_chem_inca_trac(nqINCA)                               !--- Get nqINCA from INCA
260      nbtr = nqINCA + nqCO2                                          !--- Number of tracers passed to phytrac
261      nqtrue = nbtr + nqo                                            !--- Total number of "true" tracers
262      IF(ALL([2,3] /= nqo)) CALL abort_gcm(modname, 'Only 2 or 3 water phases allowed ; found nqo='//TRIM(int2str(nqo)), 1)
263      ALLOCATE(hadv(nqtrue), hadv_inca(nqINCA), conv_flg_inca(nqINCA), solsym_inca(nqINCA))
264      ALLOCATE(vadv(nqtrue), vadv_inca(nqINCA),  pbl_flg_inca(nqINCA))
265      CALL init_transport(solsym_inca, conv_flg_inca, pbl_flg_inca, hadv_inca, vadv_inca)
266      ALLOCATE(ttr(nqtrue))
267      ttr(1:nqo+nqCO2)                    = tracers
268      ttr(1    :      nqo   )%component   = 'lmdz'
269      ttr(1+nqo:nqCO2+nqo   )%component   = 'co2i'
270      ttr(1+nqo+nqCO2:nqtrue)%component   = 'inca'
271      ttr(1+nqo      :nqtrue)%name        = [('CO2     ', k=1, nqCO2), solsym_inca]
272      ttr(1+nqo+nqCO2:nqtrue)%parent      = tran0
273      ttr(1+nqo+nqCO2:nqtrue)%phase       = 'g'
274      lerr = getKey('hadv', had, ky=tracers(:)%keys)
275      lerr = getKey('vadv', vad, ky=tracers(:)%keys)
276      hadv(1:nqo) = had(:); hadv(nqo+1:nqtrue) = hadv_inca
277      vadv(1:nqo) = vad(:); vadv(nqo+1:nqtrue) = vadv_inca
278      CALL MOVE_ALLOC(FROM=ttr, TO=tracers)
279      CALL setGeneration(tracers)                                    !--- SET FIELDS %iGeneration, %gen0Name
280      DEALLOCATE(had, hadv_inca, vad, vadv_inca, conv_flg_inca, pbl_flg_inca, solsym_inca)
281#endif
282   !---------------------------------------------------------------------------------------------------------------------------
283   ELSE                                                              !=== OTHER CASES (OLD OR NEW FORMAT, NO INCA MODULE)
284   !---------------------------------------------------------------------------------------------------------------------------
285      nqo    =        COUNT(delPhase(tracers(:)%name)     == 'H2O' &
286                               .AND. tracers(:)%component == 'lmdz') !--- Number of water phases
287      nqtrue = SIZE(tracers)                                         !--- Total number of "true" tracers
288      nbtr   = nqtrue-COUNT(delPhase(tracers(:)%gen0Name) == 'H2O' &
289                               .AND. tracers(:)%component == 'lmdz') !--- Number of tracers passed to phytrac
290#ifdef INCA
291      nqINCA = COUNT(tracers(:)%component == 'inca')
292#endif
293      lerr = getKey('hadv', hadv, ky=tracers(:)%keys)
294      lerr = getKey('vadv', vadv, ky=tracers(:)%keys)
295   !---------------------------------------------------------------------------------------------------------------------------
296   END IF
297   !---------------------------------------------------------------------------------------------------------------------------
298
299   CALL getKey_init(tracers)
300
301   !--- Transfert the number of tracers to Reprobus
302#ifdef REPROBUS
303   CALL Init_chem_rep_trac(nbtr, nqo, tracers(:)%name)
304#endif
305
306!==============================================================================================================================
307! 2) Calculate nqtot, number of tracers needed (greater if advection schemes 20 or 30 have been chosen).
308!==============================================================================================================================
309   DO iq = 1, nqtrue
310      IF( hadv(iq)<20 .OR. (ANY(hadv(iq)==[20,30]) .AND. hadv(iq)==vadv(iq)) ) CYCLE
311      WRITE(msg1,'("The choice hadv=",i0,", vadv=",i0,a)')hadv(iq),vadv(iq),' for "'//TRIM(tracers(iq)%name)//'" is not available'
312      CALL abort_gcm(modname, TRIM(msg1), 1)
313   END DO
314   nqtot =    COUNT( hadv< 20 .AND. vadv< 20 ) &                     !--- No additional tracer
315         +  4*COUNT( hadv==20 .AND. vadv==20 ) &                     !--- 3  additional tracers
316         + 10*COUNT( hadv==30 .AND. vadv==30 )                       !--- 9  additional tracers
317
318   !--- More tracers due to the choice of advection scheme => assign total number of tracers
319   IF( nqtot /= nqtrue ) THEN
320      CALL msg('The choice of advection scheme for one or more tracers makes it necessary to add tracers')
321      CALL msg('The number of true tracers is '//TRIM(int2str(nqtrue)))
322      CALL msg('The total number of tracers needed is '//TRIM(int2str(nqtot)))
323   END IF
324
325!==============================================================================================================================
326! 3) Determine the advection scheme choice for water and tracers "iadv" and the fields long name, isAdvected.
327!     iadv = 1    "LMDZ-specific humidity transport" (for H2O vapour)          LMV
328!     iadv = 2    backward                           (for H2O liquid)          BAK
329!     iadv = 14   Van-Leer + specific humidity, modified by Francis Codron     VLH
330!     iadv = 10   Van-Leer (chosen for vapour and liquid water)                VL1
331!     iadv = 11   Van-Leer for hadv and PPM version (Monotonic) for vadv       VLP
332!     iadv = 12   Frederic Hourdin I                                           FH1
333!     iadv = 13   Frederic Hourdin II                                          FH2
334!     iadv = 16   Monotonic         PPM (Collela & Woodward 1984)              PPM
335!     iadv = 17   Semi-monotonic    PPM (overshoots allowed)                   PPS
336!     iadv = 18   Definite positive PPM (overshoots and undershoots allowed)   PPP
337!     iadv = 20   Slopes                                                       SLP
338!     iadv = 30   Prather                                                      PRA
339!
340!        In array q(ij,l,iq) : iq = 1/2[/3]    for vapour/liquid[/ice] water
341!        And optionaly:        iq = 3[4],nqtot for other tracers
342!==============================================================================================================================
343   ALLOCATE(ttr(nqtot))
344   jq = nqtrue+1; tracers(:)%iadv = -1
345   DO iq = 1, nqtrue
346      t1 => tracers(iq)
347
348      !--- VERIFY THE CHOICE OF ADVECTION SCHEME
349      iad = -1
350      IF(hadv(iq)     ==    vadv(iq)    ) iad = hadv(iq)
351      IF(hadv(iq)==10 .AND. vadv(iq)==16) iad = 11
352      WRITE(msg1,'("Bad choice of advection scheme for ",a,": hadv = ",i0,", vadv = ",i0)')TRIM(t1%name), hadv(iq), vadv(iq)
353      IF(iad == -1) CALL abort_gcm(modname, msg1, 1)
354
355      !--- SET FIELDS %longName, %iadv, %isAdvected, %isInPhysics
356      t1%longName   = t1%name; IF(iad > 0) t1%longName=TRIM(t1%name)//descrq(iad)
357      t1%iadv       = iad
358      t1%isAdvected = iad >= 0
359      t1%isInPhysics= delPhase(t1%gen0Name) /= 'H2O' &
360                          .OR. t1%component /= 'lmdz' !=== OTHER EXCEPTIONS TO BE ADDED: CO2i, SURSATURATED WATER CLOUD...
361      ttr(iq)       = t1
362
363      !--- DEFINE THE HIGHER ORDER TRACERS, IF ANY
364      nm = 0
365      IF(iad == 20) nm = 3                                           !--- 2nd order scheme
366      IF(iad == 30) nm = 9                                           !--- 3rd order scheme
367      IF(nm == 0) CYCLE                                              !--- No higher moments
368      ttr(jq+1:jq+nm)             = t1
369      ttr(jq+1:jq+nm)%name        = [(TRIM(t1%name)    //'-'//TRIM(suff(im)), im=1, nm) ]
370      ttr(jq+1:jq+nm)%parent      = [(TRIM(t1%parent)  //'-'//TRIM(suff(im)), im=1, nm) ]
371      ttr(jq+1:jq+nm)%longName    = [(TRIM(t1%longName)//'-'//TRIM(suff(im)), im=1, nm) ]
372      ttr(jq+1:jq+nm)%iadv        = [(-iad,    im=1, nm) ]
373      ttr(jq+1:jq+nm)%isAdvected  = [(.FALSE., im=1, nm) ]
374      jq = jq + nm
375   END DO
376   DEALLOCATE(hadv, vadv)
377   CALL MOVE_ALLOC(FROM=ttr, TO=tracers)
378
379   !--- SET FIELDS %iqParent, %nqChilds, %iGeneration, %iqDescen, %nqDescen
380   CALL indexUpdate(tracers)
381
382   !=== TEST ADVECTION SCHEME
383   DO iq=1,nqtot ; t1 => tracers(iq); iad = t1%iadv
384
385      !--- ONLY TESTED VALUES FOR TRACERS FOR NOW:               iadv = 14, 10 (and 0 for non-transported tracers)
386      IF(ALL([10,14,0] /= iad)) &
387         CALL abort_gcm(modname, 'Not tested for iadv='//TRIM(int2str(iad))//' ; 10 or 14 only are allowed !', 1)
388
389      !--- ONLY TESTED VALUES FOR PARENTS HAVING CHILDS FOR NOW: iadv = 14, 10 (PARENTS: TRACERS OF GENERATION 1)
390      IF(ALL([10,14] /= iad) .AND. t1%iGeneration == 1 .AND. ANY(tracers(:)%iGeneration > 1)) &
391         CALL abort_gcm(modname, 'iadv='//TRIM(int2str(iad))//' not implemented for parents ; 10 or 14 only are allowed !', 1)
392
393      !--- ONLY TESTED VALUES FOR CHILDS FOR NOW:                iadv = 10     (CHILDS:  TRACERS OF GENERATION GREATER THAN 1)
394      IF(fmsg('WARNING ! iadv='//TRIM(int2str(iad))//' not implemented for childs. Setting iadv=10 for "'//TRIM(t1%name)//'"',&
395         modname, iad /= 10 .AND. t1%iGeneration > 1)) t1%iadv = 10
396
397      !--- ONLY VALID SCHEME NUMBER FOR WATER VAPOUR:            iadv = 14
398      ll = t1%name /= addPhase('H2O','g')
399      IF(fmsg('WARNING ! iadv=14 is valid for water vapour only. Setting iadv=10 for "'//TRIM(t1%name)//'".', &
400         modname, iad == 14 .AND. ll))                 t1%iadv = 10
401   END DO
402
403   niso = 0; nzone=0; nphas=nqo; ntiso = 0; isoCheck=.FALSE.
404   IF(initIsotopes(tracers, isotopes)) CALL abort_gcm(modname, 'Problem when reading isotopes parameters', 1)
405   nbIso = SIZE(isotopes)
406   nqtottr = nqtot - COUNT(delPhase(tracers%gen0Name) == 'H2O' .AND. tracers%component == 'lmdz')
407   IF(nbIso/=0) THEN                        !--- ISOTOPES FOUND
408
409      !=== READ PHYSICAL PARAMETERS FROM "isotopes_params.def" FILE SPECIFIC TO WATER ISOTOPES
410      !    DONE HERE, AND NOT ONLY IN "infotrac_phy", BECAUSE SOME PHYSICAL PARAMS ARE NEEDED FOR RESTARTS (tnat, alpha_ideal)
411      CALL getKey_init(tracers, isotopes)
412      IF(isoSelect('H2O')) RETURN           !--- Select water isotopes ; finished if no water isotopes
413      iH2O = ixIso                          !--- Keep track of water family index
414      IF(getKey('tnat' , tnat,        isoName(1:niso))) CALL abort_gcm(modname, 'can''t read "tnat"', 1)
415      IF(getKey('alpha', alpha_ideal, isoName(1:niso))) CALL abort_gcm(modname, 'can''t read "alpha_ideal"', 1)
416
417      !=== MAKE SURE THE MEMBERS OF AN ISOTOPES FAMILY ARE PRESENT IN THE SAME PHASES
418      DO ix = 1, nbIso
419         iso => isotopes(ix)
420         !--- Check whether each isotope and tagging isotopic tracer is present in the same number of phases
421         DO it = 1, iso%ntiso
422            np = SUM([(COUNT(tracers(:)%name == addPhase(iso%trac(it), iso%phase(ip:ip))), ip=1, iso%nphas)])
423            IF(np == iso%nphas) CYCLE
424            WRITE(msg1,'("Found ",i0," phases for ",a," instead of ",i0)')np, TRIM(iso%trac(it)), iso%nphas
425            CALL abort_gcm(modname, msg1, 1)
426         END DO
427         DO it = 1, iso%niso
428            nz = SUM([(COUNT(iso%trac == TRIM(iso%trac(it))//'_'//iso%zone(iz)), iz=1, iso%nzone)])
429            IF(nz == iso%nzone) CYCLE
430            WRITE(msg1,'("Found ",i0," tagging zones for ",a," instead of ",i0)')nz, TRIM(iso%trac(it)), iso%nzone
431            CALL abort_gcm(modname, msg1, 1)
432         END DO
433      END DO
434   END IF
435
436   !--- Convection / boundary layer activation for all tracers
437   ALLOCATE(conv_flg(nbtr)); conv_flg(1:nbtr) = 1
438   ALLOCATE( pbl_flg(nbtr));  pbl_flg(1:nbtr) = 1
439
440   !--- Note: nqtottr can differ from nbtr when nmom/=0
441   IF(COUNT(tracers%iso_iName == 0) - COUNT(delPhase(tracers%name) == 'H2O' .AND. tracers%component == 'lmdz') /= nqtottr) &
442      CALL abort_gcm('infotrac_init', 'pb dans le calcul de nqtottr', 1)
443
444   !=== DISPLAY THE RESULTS
445   CALL msg('nqo    = '//TRIM(int2str(nqo)),    modname)
446   CALL msg('nbtr   = '//TRIM(int2str(nbtr)),   modname)
447   CALL msg('nqtrue = '//TRIM(int2str(nqtrue)), modname)
448   CALL msg('nqtot  = '//TRIM(int2str(nqtot)),  modname)
449   CALL msg('niso   = '//TRIM(int2str(niso)),   modname)
450   CALL msg('ntiso  = '//TRIM(int2str(ntiso)),  modname)
451#ifdef INCA
452   CALL msg('nqCO2  = '//TRIM(int2str(nqCO2)),  modname)
453   CALL msg('nqINCA = '//TRIM(int2str(nqINCA)), modname)
454#endif
455   t => tracers
456   CALL msg('Information stored in infotrac :', modname)
457   IF(dispTable('isssssssssiiiiiiiii', &
458      ['iq    ', 'name  ', 'lName ', 'gen0N ', 'parent', 'type  ', 'phase ', 'compon', 'isAdv ', 'isPhy ', &
459       'iadv  ', 'iGen  ', 'iqPar ', 'nqDes ', 'nqChld', 'iGroup', 'iName ', 'iZone ', 'iPhase'],          &
460      cat(t%name, t%longName, t%gen0Name, t%parent, t%type, t%phase, t%component, bool2str(t%isAdvected),  &
461                                                                                  bool2str(t%isInPhysics)),&
462      cat([(iq, iq=1, nqtot)], t%iadv, t%iGeneration, t%iqParent, t%nqDescen, t%nqChilds, t%iso_iGroup,    &
463                  t%iso_iName, t%iso_iZone, t%iso_iPhase), nColMax=maxTableWidth, nHead=2, sub=modname))   &
464      CALL abort_gcm(modname, "problem with the tracers table content", 1)
465   IF(niso > 0) THEN
466      CALL msg('Where, for isotopes family "'//TRIM(isotope%parent)//'":', modname)
467      CALL msg('  isoKeys = '//strStack(isoKeys%name), modname)
468      CALL msg('  isoName = '//strStack(isoName),      modname)
469      CALL msg('  isoZone = '//strStack(isoZone),      modname)
470      CALL msg('  isoPhas = '//TRIM(isoPhas),          modname)
471   ELSE
472      CALL msg('No isotopes identified.', modname)
473   END IF
474   CALL msg('end', modname)
475
476END SUBROUTINE infotrac_init
477
478
479!==============================================================================================================================
480!=== THE ROUTINE isoSelect IS USED TO SWITCH FROM AN ISOTOPE FAMILY TO ANOTHER: ISOTOPES DEPENDENT PARAMETERS ARE UPDATED
481!     Single generic "isoSelect" routine, using the predefined index of the parent (fast version) or its name (first call).
482!==============================================================================================================================
483LOGICAL FUNCTION isoSelectByName(iName, lVerbose) RESULT(lerr)
484   IMPLICIT NONE
485   CHARACTER(LEN=*),  INTENT(IN) :: iName
486   LOGICAL, OPTIONAL, INTENT(IN) :: lVerbose
487   INTEGER :: iIso
488   LOGICAL :: lV
489   lV = .FALSE.; IF(PRESENT(lVerbose)) lV = lVerbose
490   iIso = strIdx(isotopes(:)%parent, iName)
491   lerr = iIso == 0
492   IF(lerr) THEN
493      niso = 0; ntiso = 0; nzone=0; nphas=nqo; isoCheck=.FALSE.
494      CALL msg('no isotope family named "'//TRIM(iName)//'"', ll=lV)
495      RETURN
496   END IF
497   lerr = isoSelectByIndex(iIso, lV)
498END FUNCTION isoSelectByName
499!==============================================================================================================================
500LOGICAL FUNCTION isoSelectByIndex(iIso, lVerbose) RESULT(lerr)
501   IMPLICIT NONE
502   INTEGER,           INTENT(IN) :: iIso
503   LOGICAL, OPTIONAL, INTENT(IN) :: lVerbose
504   LOGICAL :: lv
505   lv = .FALSE.; IF(PRESENT(lVerbose)) lv = lVerbose
506   lerr = .FALSE.
507   IF(iIso == ixIso) RETURN                                          !--- Nothing to do if the index is already OK
508   lerr = iIso<=0 .OR. iIso>nbIso
509   CALL msg('Inconsistent isotopes family index '//TRIM(int2str(iIso))//': should be > 0 and <= '//TRIM(int2str(nbIso))//'"',&
510            ll=lerr .AND. lV)
511   IF(lerr) RETURN
512   ixIso = iIso                                                      !--- Update currently selected family index
513   isotope  => isotopes(ixIso)                                       !--- Select corresponding component
514   isoKeys  => isotope%keys;     niso     = isotope%niso
515   isoName  => isotope%trac;     ntiso    = isotope%ntiso
516   isoZone  => isotope%zone;     nzone    = isotope%nzone
517   isoPhas  => isotope%phase;    nphas    = isotope%nphas
518   itZonIso => isotope%itZonIso; isoCheck = isotope%check
519   iqIsoPha => isotope%iqIsoPha
520END FUNCTION isoSelectByIndex
521!==============================================================================================================================
522
523END MODULE infotrac
Note: See TracBrowser for help on using the repository browser.