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

Last change on this file since 4158 was 4158, checked in by dcugnet, 2 years ago
  • include new type "tracer.def" tracers description files in DefLists/?
  • Remove the old code (lOldCode=.TRUE.)
  • Remove the hard-coded part in isotopes_mod that was ensuring the convergence => this version will not converge with previous isotopic physics, but the

difference, only due to round off errors on "tdifrel", is acceptable.

  • 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.3 KB
Line 
1!$Id: infotrac.F90 4158 2022-05-19 08:47:23Z 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, addPhase, indexUpdate,  nphases, ancestor,  &
8                                isot_type, old2newName,      delPhase,               getKey_init, tran0, &
9                                keys_type, initIsotopes,     getPhase, known_phases, getKey, setGeneration
10   IMPLICIT NONE
11
12   PRIVATE
13
14   !=== FOR TRACERS:
15   PUBLIC :: infotrac_init                                 !--- Initialization of the tracers
16   PUBLIC :: tracers, type_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
111
112   !=== DERIVED TYPES EMBEDDING MOST INFORMATIONS ABOUT TRACERS AND ISOTOPES FAMILIES
113   TYPE(trac_type), TARGET, SAVE, ALLOCATABLE ::  tracers(:)    !=== TRACERS DESCRIPTORS VECTOR
114   TYPE(isot_type), TARGET, SAVE, ALLOCATABLE :: isotopes(:)    !=== ISOTOPES PARAMETERS VECTOR
115
116   !=== ALIASES FOR CURRENTLY SELECTED ISOTOPES FAMILY OF VARIABLES EMBEDDED IN THE VECTOR "isotopes"
117   TYPE(isot_type),         SAVE, POINTER :: isotope            !--- CURRENTLY SELECTED ISOTOPES FAMILY DESCRIPTOR
118   INTEGER,                 SAVE          :: ixIso, iH2O        !--- Index of the selected isotopes family and H2O family
119   LOGICAL,                 SAVE          :: isoCheck           !--- Flag to trigger the checking routines
120   TYPE(keys_type),         SAVE, POINTER :: isoKeys(:)         !--- ONE SET OF KEYS FOR EACH ISOTOPE (LISTED IN isoName)
121   CHARACTER(LEN=maxlen),   SAVE, POINTER :: isoName(:),   &    !--- ISOTOPES NAMES FOR THE CURRENTLY SELECTED FAMILY
122                                             isoZone(:),   &    !--- TAGGING ZONES  FOR THE CURRENTLY SELECTED FAMILY
123                                             isoPhas            !--- USED PHASES    FOR THE CURRENTLY SELECTED FAMILY
124   INTEGER,                 SAVE          ::  niso, nzone, &    !--- NUMBER OF ISOTOPES, TAGGING ZONES AND PHASES
125                                             nphas, ntiso       !--- NUMBER OF PHASES AND ISOTOPES + ISOTOPIC TAGGING TRACERS
126   INTEGER,                 SAVE, POINTER ::itZonIso(:,:), &    !--- INDEX IN "isoTrac" AS f(tagging zone idx,  isotope idx)
127                                            iqIsoPha(:,:)       !--- INDEX IN "qx"      AS f(isotopic tracer idx, phase idx)
128
129   !=== VARIABLES FOR ISOTOPES INITIALIZATION AND FOR INCA
130   REAL,                SAVE, ALLOCATABLE ::     tnat(:), &     !--- Natural relative abundance of water isotope        (niso)
131                                          alpha_ideal(:)        !--- Ideal fractionning coefficient (for initial state) (niso)
132   INTEGER,             SAVE, ALLOCATABLE :: conv_flg(:), &     !--- Convection     activation ; needed for INCA        (nbtr)
133                                              pbl_flg(:)        !--- Boundary layer activation ; needed for INCA        (nbtr)
134
135CONTAINS
136
137SUBROUTINE infotrac_init
138   USE control_mod, ONLY: planet_type, config_inca
139#ifdef REPROBUS
140   USE CHEM_REP,    ONLY: Init_chem_rep_trac
141#endif
142   IMPLICIT NONE
143!==============================================================================================================================
144!
145!   Auteur:  P. Le Van /L. Fairhead/F.Hourdin
146!   -------
147!
148!   Modifications:
149!   --------------
150!   05/94: F.Forget      Modif special traceur
151!   02/02: M-A Filiberti Lecture de traceur.def
152!   01/22: D. Cugnet     Nouveaux tracer.def et tracer_*.def + encapsulation (types trac_type et isot_type)
153!
154!   Objet:
155!   ------
156!   GCM LMD nouvelle grille
157!
158!==============================================================================================================================
159!   ... modification de l'integration de q ( 26/04/94 ) ....
160!------------------------------------------------------------------------------------------------------------------------------
161! Declarations:
162   INCLUDE "dimensions.h"
163   INCLUDE "iniprint.h"
164
165!------------------------------------------------------------------------------------------------------------------------------
166! Local variables
167   INTEGER, ALLOCATABLE :: hadv(:), vadv(:)                          !--- Horizontal/vertical transport scheme number
168#ifdef INCA
169   INTEGER, ALLOCATABLE :: had (:), hadv_inca(:), conv_flg_inca(:), &!--- Variables specific to INCA
170                           vad (:), vadv_inca(:),  pbl_flg_inca(:)
171   CHARACTER(LEN=8), ALLOCATABLE :: solsym_inca(:)                   !--- Tracers names for INCA
172   INTEGER :: nqINCA
173#endif
174   CHARACTER(LEN=2)      ::   suff(9)                                !--- Suffixes for schemes of order 3 or 4 (Prather)
175   CHARACTER(LEN=3)      :: descrq(30)                               !--- Advection scheme description tags
176   CHARACTER(LEN=maxlen) :: msg1                                     !--- String for messages
177   CHARACTER(LEN=maxlen), ALLOCATABLE :: str(:)                      !--- Temporary storage
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
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
189   CHARACTER(LEN=maxlen), ALLOCATABLE :: tnom_0(:), tnom_transp(:)   !--- Tracer short name + transporting fluid name
190   CHARACTER(LEN=maxlen)              :: tchaine
191   INTEGER :: ierr
192
193   CHARACTER(LEN=*), PARAMETER :: modname="infotrac_init"
194!------------------------------------------------------------------------------------------------------------------------------
195! Initialization :
196!------------------------------------------------------------------------------------------------------------------------------
197   suff          = ['x ','y ','z ','xx','xy','xz','yy','yz','zz']
198   descrq( 1: 2) = ['LMV','BAK']
199   descrq(10:20) = ['VL1','VLP','FH1','FH2','VLH','   ','PPM','PPS','PPP','   ','SLP']
200   descrq(30)    =  'PRA'
201   
202   CALL msg('type_trac = "'//TRIM(type_trac)//'"', modname)
203   IF(strParse(type_trac, ',', str, n=nt)) CALL abort_gcm(modname,'can''t parse "type_trac = '//TRIM(type_trac)//'"',1)
204
205   !---------------------------------------------------------------------------------------------------------------------------
206   DO it = 1, nt                                                          !--- nt>1=> "type_trac": coma-separated keywords list
207   !---------------------------------------------------------------------------------------------------------------------------
208      !--- MESSAGE ABOUT THE CHOSEN CONFIGURATION
209      msg1 = 'For type_trac = "'//TRIM(str(it))//'":'
210      SELECT CASE(str(it))
211         CASE('inca'); CALL msg(TRIM(msg1)//' coupling with INCA chemistry model, config_inca='//config_inca, modname)
212         CASE('inco'); CALL msg(TRIM(msg1)//' coupling jointly with INCA and CO2 cycle',  modname)
213         CASE('repr'); CALL msg(TRIM(msg1)//' coupling with REPROBUS chemistry model',    modname)
214         CASE('co2i'); CALL msg(TRIM(msg1)//' you have chosen to run with CO2 cycle',     modname)
215         CASE('coag'); CALL msg(TRIM(msg1)//' tracers are treated for COAGULATION tests', modname)
216         CASE('lmdz'); CALL msg(TRIM(msg1)//' tracers are treated in LMDZ only',          modname)
217         CASE DEFAULT; CALL abort_gcm(modname,'type_trac='//TRIM(str(it))//' not possible yet.',1)
218      END SELECT
219
220      !--- COHERENCE TEST BETWEEN "type_trac" AND "config_inca"
221      IF(ANY(['inca', 'inco'] == str(it)) .AND. ALL(['aero', 'aeNP', 'chem'] /= config_inca)) &
222         CALL abort_gcm(modname, 'Incoherence between type_trac and config_inca. Please modify "run.def"',1)
223
224      !--- COHERENCE TEST BETWEEN "type_trac" AND PREPROCESSING KEYS
225      SELECT CASE(str(it))
226         CASE('inca', 'inco')
227#ifndef INCA
228            CALL abort_gcm(modname, 'You must add cpp key INCA and compile with INCA code', 1)
229#endif
230         CASE('repr')
231#ifndef REPROBUS
232            CALL abort_gcm(modname, 'You must add cpp key REPROBUS and compile with REPROBUS code', 1)
233#endif
234         CASE('coag')
235#ifndef CPP_StratAer
236            CALL abort_gcm(modname, 'You must add cpp key StratAer and compile with StratAer code', 1)
237#endif
238      END SELECT
239
240      !--- DISABLE "config_inca" OPTION FOR A RUN WITHOUT "INCA" IF IT DIFFERS FROM "none"
241      IF(fmsg('Setting config_inca="none" as you do not couple with INCA model', &
242              modname, ALL(['inca', 'inco'] /= str(it)) .AND. config_inca /= 'none')) config_inca = 'none'
243
244   !---------------------------------------------------------------------------------------------------------------------------
245   END DO
246   !---------------------------------------------------------------------------------------------------------------------------
247
248   nqCO2 = 0; IF(ANY(str == 'inco')) nqCO2 = 1
249
250!==============================================================================================================================
251! 1) Get the numbers of: true (first order only) tracers "nqtrue", water tracers "nqo" (vapor/liquid/solid)
252!==============================================================================================================================
253   IF(readTracersFiles(type_trac, fType, tracers)) CALL abort_gcm(modname, 'problem with tracers file(s)',1)
254   IF(fType == 0) CALL abort_gcm(modname, 'Missing tracers file: "traceur.def", "tracer.def" or "tracer_<keyword>.def file.',1)
255   !---------------------------------------------------------------------------------------------------------------------------
256   IF(fType == 1 .AND. ANY(['inca','inco'] == type_trac)) THEN       !=== FOUND OLD STYLE INCA "traceur.def" (single type_trac)
257   !---------------------------------------------------------------------------------------------------------------------------
258#ifdef INCA
259      nqo = SIZE(tracers)
260      IF(nqCO2==1 .AND. nqo==4) nqo = 3                              !--- Force nqo to 3 (ThL)
261      CALL Init_chem_inca_trac(nqINCA)                               !--- Get nqINCA from INCA
262      nbtr = nqINCA + nqCO2                                          !--- Number of tracers passed to phytrac
263      nqtrue = nbtr + nqo                                            !--- Total number of "true" tracers
264      IF(ALL([2,3] /= nqo)) CALL abort_gcm(modname, 'Only 2 or 3 water phases allowed ; found nqo='//TRIM(int2str(nqo)), 1)
265      ALLOCATE(hadv(nqtrue), hadv_inca(nqINCA), conv_flg_inca(nqINCA), solsym_inca(nqINCA))
266      ALLOCATE(vadv(nqtrue), vadv_inca(nqINCA),  pbl_flg_inca(nqINCA))
267      CALL init_transport(solsym_inca, conv_flg_inca, pbl_flg_inca, hadv_inca, vadv_inca)
268      ALLOCATE(ttr(nqtrue))
269      ttr(1:nqo+nqCO2)                    = tracers
270      ttr(1    :      nqo   )%component   = 'lmdz'
271      ttr(1+nqo:nqCO2+nqo   )%component   = 'co2i'
272      ttr(1+nqo+nqCO2:nqtrue)%component   = 'inca'
273      ttr(1+nqo      :nqtrue)%name        = [('CO2     ', k=1, nqCO2), solsym_inca]
274      ttr(1+nqo+nqCO2:nqtrue)%parent      = tran0
275      ttr(1+nqo+nqCO2:nqtrue)%phase       = 'g'
276      lerr = getKey('hadv', had, ky=ttr(:)%keys); hadv(:) = [had, hadv_inca]
277      lerr = getKey('vadv', vad, ky=ttr(:)%keys); vadv(:) = [vad, 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                                                              !=== FOUND NEW STYLE TRACERS CONFIGURATION FILE(S)
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      lerr = getKey('hadv', hadv, ky=tracers(:)%keys)
291      lerr = getKey('vadv', vadv, ky=tracers(:)%keys)
292   !---------------------------------------------------------------------------------------------------------------------------
293   END IF
294   !---------------------------------------------------------------------------------------------------------------------------
295
296   CALL getKey_init(tracers)
297
298   !--- Transfert the number of tracers to Reprobus
299#ifdef REPROBUS
300   CALL Init_chem_rep_trac(nbtr, nqo, tracers(:)%name)
301#endif
302
303!==============================================================================================================================
304! 2) Calculate nqtot, number of tracers needed (greater if advection schemes 20 or 30 have been chosen).
305!==============================================================================================================================
306   DO iq = 1, nqtrue
307      IF( hadv(iq)<20 .OR. (ANY(hadv(iq)==[20,30]) .AND. hadv(iq)==vadv(iq)) ) CYCLE
308      WRITE(msg1,'("The choice hadv=",i0,", vadv=",i0,a)')hadv(iq),vadv(iq),' for "'//TRIM(tracers(iq)%name)//'" is not available'
309      CALL abort_gcm(modname, TRIM(msg1), 1)
310   END DO
311   nqtot =    COUNT( hadv< 20 .AND. vadv< 20 ) &                     !--- No additional tracer
312         +  4*COUNT( hadv==20 .AND. vadv==20 ) &                     !--- 3  additional tracers
313         + 10*COUNT( hadv==30 .AND. vadv==30 )                       !--- 9  additional tracers
314
315   !--- More tracers due to the choice of advection scheme => assign total number of tracers
316   IF( nqtot /= nqtrue ) THEN
317      CALL msg('The choice of advection scheme for one or more tracers makes it necessary to add tracers')
318      CALL msg('The number of true tracers is '//TRIM(int2str(nqtrue)))
319      CALL msg('The total number of tracers needed is '//TRIM(int2str(nqtot)))
320   END IF
321
322!==============================================================================================================================
323! 3) Determine the advection scheme choice for water and tracers "iadv" and the fields long name, isAdvected.
324!     iadv = 1    "LMDZ-specific humidity transport" (for H2O vapour)          LMV
325!     iadv = 2    backward                           (for H2O liquid)          BAK
326!     iadv = 14   Van-Leer + specific humidity, modified by Francis Codron     VLH
327!     iadv = 10   Van-Leer (chosen for vapour and liquid water)                VL1
328!     iadv = 11   Van-Leer for hadv and PPM version (Monotonic) for vadv       VLP
329!     iadv = 12   Frederic Hourdin I                                           FH1
330!     iadv = 13   Frederic Hourdin II                                          FH2
331!     iadv = 16   Monotonic         PPM (Collela & Woodward 1984)              PPM
332!     iadv = 17   Semi-monotonic    PPM (overshoots allowed)                   PPS
333!     iadv = 18   Definite positive PPM (overshoots and undershoots allowed)   PPP
334!     iadv = 20   Slopes                                                       SLP
335!     iadv = 30   Prather                                                      PRA
336!
337!        In array q(ij,l,iq) : iq = 1/2[/3]    for vapour/liquid[/ice] water
338!        And optionaly:        iq = 3[4],nqtot for other tracers
339!==============================================================================================================================
340   ALLOCATE(ttr(nqtot))
341   jq = nqtrue+1; tracers(:)%iadv = -1
342   DO iq = 1, nqtrue
343      t1 => tracers(iq)
344
345      !--- VERIFY THE CHOICE OF ADVECTION SCHEME
346      iad = -1
347      IF(hadv(iq)     ==    vadv(iq)    ) iad = hadv(iq)
348      IF(hadv(iq)==10 .AND. vadv(iq)==16) iad = 11
349      WRITE(msg1,'("Bad choice of advection scheme for ",a,": hadv = ",i0,", vadv = ",i0)')TRIM(t1%name), hadv(iq), vadv(iq)
350      IF(iad == -1) CALL abort_gcm(modname, msg1, 1)
351
352      !--- SET FIELDS %longName, %iadv, %isAdvected, %isInPhysics
353      t1%longName   = t1%name; IF(iad > 0) t1%longName=TRIM(t1%name)//descrq(iad)
354      t1%iadv       = iad
355      t1%isAdvected = iad >= 0
356      t1%isInPhysics= delPhase(t1%gen0Name) /= 'H2O' &
357                          .OR. t1%component /= 'lmdz' !=== OTHER EXCEPTIONS TO BE ADDED: CO2i, SURSATURATED WATER CLOUD...
358      ttr(iq)       = t1
359
360      !--- DEFINE THE HIGHER ORDER TRACERS, IF ANY
361      nm = 0
362      IF(iad == 20) nm = 3                                           !--- 2nd order scheme
363      IF(iad == 30) nm = 9                                           !--- 3rd order scheme
364      IF(nm == 0) CYCLE                                              !--- No higher moments
365      ttr(jq+1:jq+nm)             = t1
366      ttr(jq+1:jq+nm)%name        = [(TRIM(t1%name)    //'-'//TRIM(suff(im)), im=1, nm) ]
367      ttr(jq+1:jq+nm)%parent      = [(TRIM(t1%parent)  //'-'//TRIM(suff(im)), im=1, nm) ]
368      ttr(jq+1:jq+nm)%longName    = [(TRIM(t1%longName)//'-'//TRIM(suff(im)), im=1, nm) ]
369      ttr(jq+1:jq+nm)%iadv        = [(-iad,    im=1, nm) ]
370      ttr(jq+1:jq+nm)%isAdvected  = [(.FALSE., im=1, nm) ]
371      jq = jq + nm
372   END DO
373   DEALLOCATE(hadv, vadv)
374   CALL MOVE_ALLOC(FROM=ttr, TO=tracers)
375
376   !--- SET FIELDS %iqParent, %nqChilds, %iGeneration, %iqDescen, %nqDescen
377   CALL indexUpdate(tracers)
378
379   !=== TEST ADVECTION SCHEME
380   DO iq=1,nqtot ; t1 => tracers(iq); iad = t1%iadv
381
382      !--- ONLY TESTED VALUES FOR TRACERS FOR NOW:               iadv = 14, 10 (and 0 for non-transported tracers)
383      IF(ALL([10,14,0] /= iad)) &
384         CALL abort_gcm(modname, 'Not tested for iadv='//TRIM(int2str(iad))//' ; 10 or 14 only are allowed !', 1)
385
386      !--- ONLY TESTED VALUES FOR PARENTS HAVING CHILDS FOR NOW: iadv = 14, 10 (PARENTS: TRACERS OF GENERATION 1)
387      IF(ALL([10,14] /= iad) .AND. t1%iGeneration == 1 .AND. ANY(tracers(:)%iGeneration > 1)) &
388         CALL abort_gcm(modname, 'iadv='//TRIM(int2str(iad))//' not implemented for parents ; 10 or 14 only are allowed !', 1)
389
390      !--- ONLY TESTED VALUES FOR CHILDS FOR NOW:                iadv = 10     (CHILDS:  TRACERS OF GENERATION GREATER THAN 1)
391      IF(fmsg('WARNING ! iadv='//TRIM(int2str(iad))//' not implemented for childs. Setting iadv=10 for "'//TRIM(t1%name)//'"',&
392         modname, iad /= 10 .AND. t1%iGeneration > 1)) t1%iadv = 10
393
394      !--- ONLY VALID SCHEME NUMBER FOR WATER VAPOUR:            iadv = 14
395      ll = t1%name /= addPhase('H2O','g')
396      IF(fmsg('WARNING ! iadv=14 is valid for water vapour only. Setting iadv=10 for "'//TRIM(t1%name)//'".', &
397         modname, iad == 14 .AND. ll))                 t1%iadv = 10
398   END DO
399
400   niso = 0; nzone=0; nphas=nqo; ntiso = 0; isoCheck=.FALSE.
401   CALL initIsotopes(tracers, isotopes)
402   nbIso = SIZE(isotopes)
403   nqtottr = nqtot - COUNT(tracers%gen0Name == 'H2O' .AND. tracers%component == 'lmdz')
404   IF(nbIso/=0) THEN                        !--- ISOTOPES FOUND
405
406      !=== READ PHYSICAL PARAMETERS FROM "isotopes_params.def" FILE SPECIFIC TO WATER ISOTOPES
407      !    DONE HERE, AND NOT ONLY IN "infotrac_phy", BECAUSE SOME PHYSICAL PARAMS ARE NEEDED FOR RESTARTS (tnat, alpha_ideal)
408      CALL getKey_init(tracers, isotopes)
409      IF(isoSelect('H2O')) RETURN           !--- Select water isotopes ; finished if no water isotopes
410      iH2O = ixIso                          !--- Keep track of water family index
411      IF(getKey('tnat' , tnat,        isoName(1:niso))) CALL abort_gcm(modname, 'can''t read "tnat"', 1)
412      IF(getKey('alpha', alpha_ideal, isoName(1:niso))) CALL abort_gcm(modname, 'can''t read "alpha_ideal"', 1)
413
414      !=== MAKE SURE THE MEMBERS OF AN ISOTOPES FAMILY ARE PRESENT IN THE SAME PHASES
415      DO ix = 1, nbIso
416         iso => isotopes(ix)
417         !--- Check whether each isotope and tagging isotopic tracer is present in the same number of phases
418         DO it = 1, iso%ntiso
419            np = SUM([(COUNT(tracers(:)%name == addPhase(iso%trac(it), iso%phase(ip:ip))), ip=1, iso%nphas)])
420            IF(np == iso%nphas) CYCLE
421            WRITE(msg1,'("Found ",i0," phases for ",a," instead of ",i0)')np, TRIM(iso%trac(it)), iso%nphas
422            CALL abort_gcm(modname, msg1, 1)
423         END DO
424         DO it = 1, iso%niso
425            nz = SUM([(COUNT(iso%trac == TRIM(iso%trac(it))//'_'//iso%zone(iz)), iz=1, iso%nzone)])
426            IF(nz == iso%nzone) CYCLE
427            WRITE(msg1,'("Found ",i0," tagging zones for ",a," instead of ",i0)')nz, TRIM(iso%trac(it)), iso%nzone
428            CALL abort_gcm(modname, msg1, 1)
429         END DO
430      END DO
431   END IF
432
433   !--- Convection / boundary layer activation for all tracers
434   ALLOCATE(conv_flg(nbtr)); conv_flg(1:nbtr) = 1
435   ALLOCATE( pbl_flg(nbtr));  pbl_flg(1:nbtr) = 1
436
437   !--- Note: nqtottr can differ from nbtr when nmom/=0
438!   IF(COUNT(tracers%iso_iName == 0) - COUNT(delPhase(tracers%name) == 'H2O' .AND. tracers%component == 'lmdz') /= nqtottr) &
439!      CALL abort_gcm('infotrac_init', 'pb dans le calcul de nqtottr', 1)
440
441   !=== DISPLAY THE RESULTS
442   CALL msg('nqo    = '//TRIM(int2str(nqo)),    modname)
443   CALL msg('nbtr   = '//TRIM(int2str(nbtr)),   modname)
444   CALL msg('nqtrue = '//TRIM(int2str(nqtrue)), modname)
445   CALL msg('nqtot  = '//TRIM(int2str(nqtot)),  modname)
446   CALL msg('niso   = '//TRIM(int2str(niso)),   modname)
447   CALL msg('ntiso  = '//TRIM(int2str(ntiso)),  modname)
448#ifdef INCA
449   CALL msg('nqCO2  = '//TRIM(int2str(nqCO2)),  modname)
450   CALL msg('nqINCA = '//TRIM(int2str(nqINCA)), modname)
451#endif
452   t => tracers
453   CALL msg('Information stored in infotrac :', modname)
454   IF(dispTable('isssssssssiiiiiiiii', &
455      ['iq    ', 'name  ', 'lName ', 'gen0N ', 'parent', 'type  ', 'phase ', 'compon', 'isAdv ', 'isPhy ', &
456       'iadv  ', 'iGen  ', 'iqPar ', 'nqDes ', 'nqChld', 'iGroup', 'iName ', 'iZone ', 'iPhase'],          &
457      cat(t%name, t%longName, t%gen0Name, t%parent, t%type, t%phase, t%component, bool2str(t%isAdvected),  &
458                                                                                  bool2str(t%isInPhysics)),&
459      cat([(iq, iq=1, nqtot)], t%iadv, t%iGeneration, t%iqParent, t%nqDescen, t%nqChilds, t%iso_iGroup,    &
460                                                    t%iso_iName, t%iso_iZone, t%iso_iPhase), sub=modname)) &
461      CALL abort_gcm(modname, "problem with the tracers table content", 1)
462   IF(niso > 0) THEN
463      CALL msg('Where, for isotopes family "'//TRIM(isotope%parent)//'":', modname)
464      CALL msg('  isoKeys = '//strStack(isoKeys%name), modname)
465      CALL msg('  isoName = '//strStack(isoName),      modname)
466      CALL msg('  isoZone = '//strStack(isoZone),      modname)
467      CALL msg('  isoPhas = '//TRIM(isoPhas),          modname)
468   ELSE
469      CALL msg('No isotopes identified.', modname)
470   END IF
471   CALL msg('end', modname)
472
473END SUBROUTINE infotrac_init
474
475
476!==============================================================================================================================
477!=== THE ROUTINE isoSelect IS USED TO SWITCH FROM AN ISOTOPE FAMILY TO ANOTHER: ISOTOPES DEPENDENT PARAMETERS ARE UPDATED
478!     Single generic "isoSelect" routine, using the predefined index of the parent (fast version) or its name (first call).
479!==============================================================================================================================
480LOGICAL FUNCTION isoSelectByName(iName, lVerbose) RESULT(lerr)
481   IMPLICIT NONE
482   CHARACTER(LEN=*),  INTENT(IN)  :: iName
483   LOGICAL, OPTIONAL, INTENT(IN) :: lVerbose
484   INTEGER :: iIso
485   LOGICAL :: lV
486   lV = .FALSE.; IF(PRESENT(lVerbose)) lV = lVerbose
487   iIso = strIdx(isotopes(:)%parent, iName)
488   lerr = iIso == 0
489   IF(lerr) THEN
490      niso = 0; ntiso = 0; nzone=0; nphas=nqo; isoCheck=.FALSE.
491      CALL msg('no isotope family named "'//TRIM(iName)//'"', ll=lV)
492      RETURN
493   END IF
494   lerr = isoSelectByIndex(iIso, lV)
495END FUNCTION isoSelectByName
496!==============================================================================================================================
497LOGICAL FUNCTION isoSelectByIndex(iIso, lVerbose) RESULT(lerr)
498   IMPLICIT NONE
499   INTEGER,           INTENT(IN) :: iIso
500   LOGICAL, OPTIONAL, INTENT(IN) :: lVerbose
501   LOGICAL :: lv
502   lv = .FALSE.; IF(PRESENT(lVerbose)) lv = lVerbose
503   lerr = .FALSE.
504   IF(iIso == ixIso) RETURN                                          !--- Nothing to do if the index is already OK
505   lerr = iIso<=0 .OR. iIso>nbIso
506   CALL msg('Inconsistent isotopes family index '//TRIM(int2str(iIso))//': should be > 0 and <= '//TRIM(int2str(nbIso))//'"',&
507            ll=lerr .AND. lV)
508   IF(lerr) RETURN
509   ixIso = iIso                                                      !--- Update currently selected family index
510   isotope  => isotopes(ixIso)                                       !--- Select corresponding component
511   isoKeys  => isotope%keys;     niso     = isotope%niso
512   isoName  => isotope%trac;     ntiso    = isotope%ntiso
513   isoZone  => isotope%zone;     nzone    = isotope%nzone
514   isoPhas  => isotope%phase;    nphas    = isotope%nphas
515   itZonIso => isotope%itZonIso; isoCheck = isotope%check
516   iqIsoPha => isotope%iqIsoPha
517END FUNCTION isoSelectByIndex
518!==============================================================================================================================
519
520END MODULE infotrac
Note: See TracBrowser for help on using the repository browser.