Ignore:
Timestamp:
May 9, 2022, 12:35:40 PM (3 years ago)
Author:
dcugnet
Message:
  • Some variables are renamed or replaced by direct equivalents:
    • iso_indnum -> tracers(:)%iso_iName
    • niso_possibles -> niso
    • iqiso -> iqIsoPha ; index_trac -> itZonIso
    • ok_iso_verif -> isoCheck
    • ntraceurs_zone -> nzone ; ntraciso -> ntiso
    • qperemin -> min_qparent ; masseqmin -> min_qmass ; ratiomin -> min_ratio
  • Some renamed variables are only aliased with the older name (using USE <module>, ONLY: <oldName> => <newName>) in routines where they are repeated many times.
  • Few hard-coded indexes are now computed (examples: ilic, iso, ivap, irneb, iq_vap, iq_liq, iso_H2O, iso_HDO, iso_HTO, iso_O17, iso_O18).
  • The IF(isoCheck) test is now embedded in the check_isotopes_seq and check_isotopes_loc routines (lighter calling).
Location:
LMDZ6/trunk/libf/dyn3d
Files:
6 edited
1 moved

Legend:

Unmodified
Added
Removed
  • LMDZ6/trunk/libf/dyn3d/advtrac.F90

    r4064 r4143  
    1111   !            M.A Filiberti (04/2002)
    1212   !
    13    USE infotrac,     ONLY: nqtot, tracers,ok_iso_verif
     13   USE infotrac,     ONLY: nqtot, tracers, isoCheck
    1414   USE control_mod,  ONLY: iapp_tracvl, day_step
    1515   USE comconst_mod, ONLY: dtvr
     
    215215#endif
    216216
    217    IF(ok_iso_verif) THEN
    218       WRITE(*,*) 'advtrac 227'
    219       CALL check_isotopes_seq(q,ip1jmp1,'advtrac 162')
    220    END IF
     217   IF(isoCheck) WRITE(*,*) 'advtrac 227'
     218   CALL check_isotopes_seq(q,ip1jmp1,'advtrac 162')
    221219
    222220   !-------------------------------------------------------------------------
     
    346344   END DO
    347345
    348    IF(ok_iso_verif) then
    349       WRITE(*,*) 'advtrac 402'
    350       CALL check_isotopes_seq(q,ip1jmp1,'advtrac 397')
    351    END IF
     346   IF(isoCheck) WRITE(*,*) 'advtrac 402'
     347   CALL check_isotopes_seq(q,ip1jmp1,'advtrac 397')
    352348
    353349   !-------------------------------------------------------------------------
  • LMDZ6/trunk/libf/dyn3d/check_isotopes.F90

    r4142 r4143  
    1         subroutine check_isotopes_seq(q,ip1jmp1,err_msg)
    2         USE infotrac, ONLY: nqtot, nqo, niso, ntraciso, nzone,
    3      &                     use_iso,
    4      &                     iqiso, index_trac,indnum_fn_num, tnat
    5         implicit none
     1SUBROUTINE check_isotopes_seq(q, ip1jmp1, err_msg)
     2   USE strings_mod, ONLY: maxlen, msg, strIdx, strStack, int2str, real2str
     3   USE infotrac,    ONLY: nqtot, niso, nphas, isotope, isoCheck, iqIsoPha, isoSelect, &
     4                          ntiso, iH2O, nzone, tracers, isoName,  itZonIso, tnat
     5   IMPLICIT NONE
     6   include "dimensions.h"
     7   REAL,             INTENT(INOUT) :: q(ip1jmp1,llm,nqtot)
     8   INTEGER,          INTENT(IN)    :: ip1jmp1
     9   CHARACTER(LEN=*), INTENT(IN)    :: err_msg    !--- Error message to display
     10   CHARACTER(LEN=maxlen) :: modname, msg1, nm(2)
     11   INTEGER :: ixt, ipha, k, i, iq, iiso, izon, ieau, iqeau, iqpar
     12   INTEGER, ALLOCATABLE :: ix(:)
     13   REAL    :: xtractot, xiiso, deltaD, q1, q2
     14   REAL, PARAMETER :: borne     = 1e19,  &
     15                      errmax    = 1e-8,  &       !--- Max. absolute error
     16                      errmaxrel = 1e-3,  &       !--- Max. relative error
     17                      qmin      = 1e-11, &
     18                      deltaDmax =1000.0, &
     19                      deltaDmin =-999.0, &
     20                      ridicule  = 1e-12
     21   INTEGER, SAVE :: iso_eau, iso_HDO, iso_O18, &
     22                             iso_O17, iso_HTO
     23   LOGICAL, SAVE :: first=.TRUE.
    624
    7 #include "dimensions.h"
     25   modname='check_isotopes'
     26   IF(.NOT.isoCheck)    RETURN                   !--- No need to check => finished
     27   IF(isoSelect('H2O')) RETURN                   !--- No H2O isotopes group found
     28   IF(niso == 0)        RETURN                   !--- No isotopes => finished
     29   IF(first) THEN
     30      iso_eau = strIdx(isoName,'H2[16]O')
     31      iso_HDO = strIdx(isoName,'H[2]HO')
     32      iso_O18 = strIdx(isoName,'H2[18]O')
     33      iso_O17 = strIdx(isoName,'H2[17]O')
     34      iso_HTO = strIdx(isoName,'H[3]HO')
     35      first = .FALSE.
     36   END IF
     37   CALL msg('31: err_msg='//TRIM(err_msg), modname)
    838
    9         ! inputs
    10         integer ip1jmp1
    11         real q(ip1jmp1,llm,nqtot)
    12         character*(*) err_msg ! message d''erreur à afficher
     39   !--- CHECK FOR NaNs (FOR ALL THE ISOTOPES, INCLUDING GEOGRAPHIC TAGGING TRACERS)
     40   modname = 'check_isotopes:iso_verif_noNaN'
     41   DO ixt = 1, ntiso
     42      DO ipha = 1, nphas
     43         iq = iqIsoPha(ixt,ipha)
     44         DO k = 1, llm
     45            DO i = 1, ip1jmp1
     46               IF(ABS(q(i,k,iq)) < borne) CYCLE
     47               WRITE(msg1,'(s,"(",i0,",",i0,",",i0,") = ",ES12.4)')TRIM(isoName(ixt)),i,k,iq,q(i,k,iq)
     48               CALL msg(msg1, modname)
     49               CALL abort_gcm(modname, 'Error with isotopes: '//TRIM(err_msg), 1)
     50            END DO
     51         END DO
     52      END DO
     53   END DO
    1354
    14         ! locals
    15         integer ixt,phase,k,i,iq,iiso,izone,ieau,iqeau
    16         real xtractot,xiiso
    17         real borne
    18         real qmin
    19         real errmax ! erreur maximale en absolu.
    20         real errmaxrel ! erreur maximale en relatif autorisée
    21         real deltaDmax,deltaDmin
    22         real ridicule
    23         parameter (borne=1e19)
    24         parameter (errmax=1e-8)
    25         parameter (errmaxrel=1e-3)
    26         parameter (qmin=1e-11)
    27         parameter (deltaDmax=1000.0,deltaDmin=-999.0)
    28         parameter (ridicule=1e-12)
    29         real deltaD
     55   !--- CHECK CONSERVATION (MAIN ISOTOPE AND PARENT CONCENTRATIONS MUST BE EQUAL)
     56   modname = 'check_isotopes:iso_verif_egalite'
     57   ixt = iso_eau
     58   IF(ixt /= 0) THEN
     59      DO ipha = 1, nphas
     60         iq = iqIsoPha(ixt,ipha)
     61         iqpar = tracers(iq)%iqParent
     62         DO k = 1, llm
     63            DO i = 1, ip1jmp1
     64               q1 = q(i,k,iqpar)
     65               q2 = q(i,k,iq)
     66!--- IMPROVEMENT in case at least one isotope is not negligible compared to the main isotopic form.
     67!    This would be probably required to sum from smallest to highest concentrations ; the corresponding
     68!    indices vector can be computed once only (in the initializations stage), using mean concentrations.
     69!              q2 = SUM(q(i,k,tracers(iqPar)%iqDesc), DIM=3)
     70               IF(ABS(q1-q2) <= errmax .OR. ABS(q1-q2)/MAX(MAX(ABS(q1),ABS(q2)),1e-18) <= errmaxrel) THEN
     71                  q(i,k,iq) = q1                 !--- Bidouille pour convergence
     72!                 q(i,k,tracers(iqPar)%iqDesc) = q(i,k,tracers(iqPar)%iqDesc) * q1 / q2
     73                  CYCLE
     74               END IF
     75               CALL msg('ixt, iq = '//TRIM(strStack(int2str([ixt,iq]))), modname)
     76               msg1 = '('//TRIM(strStack(int2str([i,k])))//')'
     77               CALL msg(TRIM(tracers(iqpar)%name)//TRIM(msg1)//' = '//TRIM(real2str(q1)), modname)
     78               CALL msg(TRIM(tracers(iq   )%name)//TRIM(msg1)//' = '//TRIM(real2str(q2)), modname)
     79               CALL abort_gcm(modname, 'Error with isotopes: '//TRIM(err_msg), 1)
     80            END DO
     81         END DO
     82      END DO
     83   END IF
    3084
    31         if (niso > 0) then
     85   !--- CHECK DELTA ANOMALIES
     86   modname = 'check_isotopes:iso_verif_aberrant'
     87   ix = [ iso_HDO  ,   iso_O18 ]
     88   nm = ['deltaD  ', 'deltaO18']
     89   DO iiso = 1, SIZE(ix)
     90      ixt = ix(iiso)
     91      IF(ixt  == 0) CYCLE
     92      DO ipha = 1, nphas
     93         iq = iqIsoPha(ixt,ipha)
     94         iqpar = tracers(iq)%iqParent
     95         DO k = 1, llm
     96            DO i = 1, ip1jmp1
     97               q1 = q(i,k,iqpar)
     98               q2 = q(i,k,iq)
     99!--- IMPROVEMENT in case at least one isotope is not negligible compared to the main isotopic form.
     100!    This would be probably required to sum from smallest to highest concentrations ; the corresponding
     101!    indices vector can be computed once only (in the initializations stage), using mean concentrations.
     102!              q2 = SUM(q(i,k,tracers(iqPar)%iqDesc), DIM=3)
     103               IF(q2 <= qmin) CYCLE
     104               deltaD = (q2/q1/tnat(ixt)-1.)*1000.
     105               IF(deltaD <= deltaDmax .AND. deltaD >= deltaDmin) CYCLE
     106               CALL msg('ixt, iq = '//TRIM(strStack(int2str([ixt,iq]))), modname)
     107               msg1 = '('//TRIM(strStack(int2str([i,k])))//')'
     108               CALL msg(TRIM(tracers(iqpar)%name)//TRIM(msg1)//' = '//TRIM(real2str(q1)), modname)
     109               CALL msg(TRIM(tracers(iq   )%name)//TRIM(msg1)//' = '//TRIM(real2str(q2)), modname)
     110               CALL msg(TRIM(nm(iiso))//TRIM(real2str(deltaD)), modname)
     111               CALL abort_gcm(modname, 'Error with isotopes: '//TRIM(err_msg), 1)
     112            END DO
     113         END DO
     114      END DO
     115   END DO
    32116
    33         write(*,*) 'check_isotopes 31: err_msg=',err_msg
    34         ! verifier que rien n'est NaN
    35         do ixt=1,ntraciso
    36           do phase=1,nqo
    37             iq=iqiso(ixt,phase)
    38             do k=1,llm
    39               DO i = 1,ip1jmp1
    40                 if ((q(i,k,iq).gt.-borne).and.
    41      :            (q(i,k,iq).lt.borne)) then
    42                 else !if ((x(ixt,i,j).gt.-borne).and.
    43                   write(*,*) 'erreur detectee par iso_verif_noNaN:'
    44                   write(*,*) err_msg
    45                   write(*,*) 'q,i,k,iq=',q(i,k,iq),i,k,iq
    46                   write(*,*) 'borne=',borne
    47                   call abort_gcm('check_isotopes_loc','plantage iso',0)
    48                 endif  !if ((x(ixt,i,j).gt.-borne).and.
    49               enddo !DO i = 1,ip1jmp1
    50             enddo !do k=1,llm
    51           enddo !do phase=1,nqo
    52         enddo !do ixt=1,ntraciso
     117   IF(nzone == 0) RETURN
    53118
    54         !write(*,*) 'check_isotopes 52'
    55         ! verifier que l'eau normale est OK
    56         if (use_iso(1)) then
    57           ixt=indnum_fn_num(1)
    58           do phase=1,nqo
    59             iq=iqiso(ixt,phase)
    60             do k=1,llm
    61             DO i = 1,ip1jmp1 
    62               if ((abs(q(i,k,phase)-q(i,k,iq)).gt.errmax).and.
    63      :          (abs((q(i,k,phase)-q(i,k,iq))/
    64      :           max(max(abs(q(i,k,phase)),abs(q(i,k,iq))),1e-18))
    65      :           .gt.errmaxrel)) then
    66                   write(*,*) 'erreur detectee par iso_verif_egalite:'
    67                   write(*,*) err_msg
    68                   write(*,*) 'ixt,phase=',ixt,phase
    69                   write(*,*) 'q,iq,i,k=',q(i,k,iq),iq,i,k
    70                   write(*,*) 'q(i,k,phase)=',q(i,k,phase)
    71                   call abort_gcm('check_isotopes_loc','plantage iso',0)
    72               endif !if ((abs(q(i,k,phase)-q(i,k,iq)).gt.errmax).and.
    73               ! bidouille pour éviter divergence:
    74               q(i,k,iq)= q(i,k,phase)
    75             enddo ! DO i = 1,ip1jmp1
    76             enddo !do k=1,llm
    77           enddo ! do phase=1,nqo
    78         endif !if (use_iso(1)) then
    79        
    80         !write(*,*) 'check_isotopes 78'
    81         ! verifier que HDO est raisonable
    82         if (use_iso(2)) then
    83           ixt=indnum_fn_num(2)
    84           do phase=1,nqo
    85             iq=iqiso(ixt,phase)
    86             do k=1,llm
    87             DO i = 1,ip1jmp1
    88             if (q(i,k,iq).gt.qmin) then
    89              deltaD=(q(i,k,iq)/q(i,k,phase)/tnat(2)-1)*1000
    90              if ((deltaD.gt.deltaDmax).or.(deltaD.lt.deltaDmin)) then
    91                   write(*,*) 'erreur detectee par iso_verif_aberrant:'
    92                   write(*,*) err_msg
    93                   write(*,*) 'ixt,phase=',ixt,phase
    94                   write(*,*) 'q,iq,i,k,=',q(i,k,iq),iq,i,k
    95                   write(*,*) 'q=',q(i,k,:)
    96                   write(*,*) 'deltaD=',deltaD
    97                   call abort_gcm('check_isotopes_loc','plantage iso',0)
    98              endif !if ((deltaD.gt.deltaDmax).or.(deltaD.lt.deltaDmin)) then
    99             endif !if (q(i,k,iq).gt.qmin) then
    100             enddo !DO i = 1,ip1jmp1
    101             enddo !do k=1,llm
    102           enddo ! do phase=1,nqo
    103         endif !if (use_iso(2)) then
     119   !--- CHECK FOR TAGGING TRACERS DELTAD ANOMALIES
     120   modname = 'check_isotopes:iso_verif_aberrant'
     121   IF(iso_eau /= 0 .AND. iso_HDO /= 0) THEN
     122      DO izon = 1, nzone
     123         ixt  = itZonIso(izon, iso_HDO)
     124         ieau = itZonIso(izon, iso_eau)
     125         DO ipha = 1, nphas
     126            iq    = iqIsoPha(ixt,  ipha)
     127            iqeau = iqIsoPha(ieau, ipha)
     128            DO k = 1, llm
     129               DO i = 1, ip1jmp1
     130                  q1 = q(i,k,iqeau)
     131                  q2 = q(i,k,iq)
     132                  IF(q2<=qmin) CYCLE
     133                  deltaD = (q2/q1/tnat(iso_HDO)-1.)*1000.
     134                  IF(deltaD <= deltaDmax .AND. deltaD >= deltaDmin) CYCLE
     135                  CALL msg('izon, ipha = '//TRIM(strStack(int2str([izon, ipha]))), modname)
     136                  CALL msg( 'ixt, ieau = '//TRIM(strStack(int2str([ ixt, ieau]))), modname)
     137                  msg1 = '('//TRIM(strStack(int2str([i,k])))//')'
     138                  CALL msg(TRIM(tracers(iqeau)%name)//TRIM(msg1)//' = '//TRIM(real2str(q1)), modname)
     139                  CALL msg(TRIM(tracers(iq   )%name)//TRIM(msg1)//' = '//TRIM(real2str(q2)), modname)
     140                  CALL msg('deltaD = '//TRIM(real2str(deltaD)), modname)
     141                  CALL abort_gcm(modname, 'Error with isotopes: '//TRIM(err_msg), 1)
     142               END DO
     143            END DO
     144         END DO
     145      END DO
     146   END IF
    104147
    105         !write(*,*) 'check_isotopes 103'
    106         ! verifier que O18 est raisonable
    107         if (use_iso(3)) then
    108           ixt=indnum_fn_num(3)
    109           do phase=1,nqo
    110             iq=iqiso(ixt,phase)
    111             do k=1,llm
    112             DO i = 1,ip1jmp1
    113             if (q(i,k,iq).gt.qmin) then
    114              deltaD=(q(i,k,iq)/q(i,k,phase)/tnat(3)-1)*1000
    115              if ((deltaD.gt.deltaDmax).or.(deltaD.lt.deltaDmin)) then
    116                   write(*,*) 'erreur detectee iso_verif_aberrant O18:'
    117                   write(*,*) err_msg
    118                   write(*,*) 'ixt,phase=',ixt,phase
    119                   write(*,*) 'q,iq,i,k,=',q(i,k,phase),iq,i,k
    120                   write(*,*) 'xt=',q(i,k,:)
    121                   write(*,*) 'deltaO18=',deltaD
    122                   call abort_gcm('check_isotopes_loc','plantage iso',0)
    123              endif !if ((deltaD.gt.deltaDmax).or.(deltaD.lt.deltaDmin)) then
    124             endif !if (q(i,k,iq).gt.qmin) then
    125             enddo !DO i = 1,ip1jmp1
    126             enddo !do k=1,llm
    127           enddo ! do phase=1,nqo
    128         endif !if (use_iso(2)) then
     148   !--- CHECK FOR TAGGING TRACERS CONSERVATION (PARENT AND TAGGING TRACERS SUM OVER ALL REGIONS MUST BE EQUAL)
     149   DO iiso = 1, niso
     150      DO ipha = 1, nphas
     151         iq = iqIsoPha(iiso, ipha)
     152         DO k = 1, llm
     153            DO i = 1, ip1jmp1
     154               xiiso = q(i,k,iq)
     155               xtractot = SUM(q(i, k, iqIsoPha(itZonIso(1:nzone,iiso), ipha)))
     156               IF(ABS(xtractot-xiiso) > errmax .AND. ABS(xtractot-xiiso)/MAX(MAX(ABS(xtractot),ABS(xiiso)),1e-18) > errmaxrel) THEN
     157                  CALL msg('Error in iso_verif_aberrant trac: '//TRIM(err_msg))
     158                  CALL msg('iiso, ipha = '//TRIM(strStack(int2str([iiso, ipha]))), modname)
     159                  CALL msg('q('//TRIM(strStack(int2str([i,k])))//',:) = '//TRIM(strStack(real2str(q(i,k,:)))), modname)
     160                  CALL abort_gcm(modname, 'Error with isotopes: '//TRIM(err_msg), 1)
     161               END IF
     162               IF(ABS(xtractot) <= ridicule) CYCLE
     163               DO izon = 1, nzone
     164                  q(i,k,iq) = q(i,k,iq) / xtractot * xiiso !--- Bidouille pour convergence
     165               END DO
     166            END DO
     167         END DO
     168      END DO
     169   END DO
    129170
     171END SUBROUTINE check_isotopes_seq
    130172
    131         !write(*,*) 'check_isotopes 129'
    132         if (nzone > 0) then
    133 
    134           if (use_iso(2).and.use_iso(1)) then
    135             do izone=1,nzone
    136              ixt=index_trac(izone,indnum_fn_num(2))
    137              ieau=index_trac(izone,indnum_fn_num(1))
    138              do phase=1,nqo
    139                iq=iqiso(ixt,phase)
    140                iqeau=iqiso(ieau,phase)
    141                do k=1,llm
    142                 DO i = 1,ip1jmp1
    143                 if (q(i,k,iq).gt.qmin) then
    144                  deltaD=(q(i,k,iq)/q(i,k,iqeau)/tnat(2)-1)*1000
    145                  if ((deltaD.gt.deltaDmax).or.
    146      &                   (deltaD.lt.deltaDmin)) then
    147                   write(*,*) 'erreur dans iso_verif_aberrant trac:'
    148                   write(*,*) err_msg
    149                   write(*,*) 'izone,phase=',izone,phase
    150                   write(*,*) 'ixt,ieau=',ixt,ieau
    151                   write(*,*) 'q,iq,i,k,=',q(i,k,iq),iq,i,k
    152                   write(*,*) 'deltaD=',deltaD
    153                   call abort_gcm('check_isotopes_loc','plantage iso',0)
    154                  endif !if ((deltaD.gt.deltaDmax).or.
    155                 endif !if (q(i,k,iq).gt.qmin) then
    156                 enddo !DO i = 1,ip1jmp1
    157                 enddo  ! do k=1,llm
    158               enddo ! do phase=1,nqo   
    159             enddo !do izone=1,nzone
    160           endif !if (use_iso(2).and.use_iso(1)) then
    161 
    162           do iiso=1,niso
    163            do phase=1,nqo
    164               iq=iqiso(iiso,phase)
    165               do k=1,llm
    166                 DO i = 1,ip1jmp1
    167                    xtractot=0.0
    168                    xiiso=q(i,k,iq)
    169                    do izone=1,nzone
    170                       iq=iqiso(index_trac(izone,iiso),phase)
    171                       xtractot=xtractot+ q(i,k,iq)
    172                    enddo !do izone=1,ntraceurs_zone
    173                    if ((abs(xtractot-xiiso).gt.errmax).and.
    174      :                  (abs(xtractot-xiiso)/
    175      :                  max(max(abs(xtractot),abs(xiiso)),1e-18)
    176      :                  .gt.errmaxrel)) then
    177                   write(*,*) 'erreur detectee par iso_verif_traceurs:'
    178                   write(*,*) err_msg
    179                   write(*,*) 'iiso,phase=',iiso,phase
    180                   write(*,*) 'i,k,=',i,k
    181                   write(*,*) 'q(i,k,:)=',q(i,k,:)
    182                   call abort_gcm('check_isotopes_loc','plantage iso',0)
    183                  endif !if ((abs(q(i,k,phase)-q(i,k,iq)).gt.errmax).and.
    184                  
    185                  ! bidouille pour éviter divergence:
    186                  if (abs(xtractot).gt.ridicule) then
    187                    do izone=1,nzone
    188                      ixt=index_trac(izone,iiso)
    189                      q(i,k,iq)=q(i,k,iq)/xtractot*xiiso
    190                    enddo !do izone=1,nzone               
    191                   endif !if ((abs(xtractot).gt.ridicule) then
    192                 enddo !DO i = 1,ip1jmp1
    193               enddo !do k=1,llm
    194            enddo !do phase=1,nqo
    195           enddo !do iiso=1,niso
    196 
    197         endif !if (nzone > 0)
    198 
    199         endif ! if (niso > 0)
    200         !write(*,*) 'check_isotopes 198'
    201        
    202         end
    203 
  • LMDZ6/trunk/libf/dyn3d/dynetat0.F90

    r4124 r4143  
    66! Purpose: Initial state reading.
    77!-------------------------------------------------------------------------------
    8   USE infotrac,    ONLY: nqtot, tracers, niso, iqiso, iso_indnum, iso_num, tnat, alpha_ideal, iH2O
     8  USE infotrac,    ONLY: nqtot, tracers, niso, iqIsoPha, tnat, alpha_ideal, iH2O
    99  USE strings_mod, ONLY: maxlen, msg, strStack, real2str
    1010  USE netcdf,      ONLY: NF90_OPEN,  NF90_NOWRITE, NF90_INQ_VARID, &
     
    145145#endif
    146146    ELSE IF(tracers(iq)%iso_iGroup == iH2O .AND. niso > 0) THEN                          !=== WATER ISOTOPES
    147 !     iName    = tracers(iq)%iso_iName  ! (next commit)
    148       iName    = iso_num(iq)
     147      iName    = tracers(iq)%iso_iName
    149148      iPhase   = tracers(iq)%iso_iPhase
    150149      iqParent = tracers(iq)%iqParent
     
    154153      ELSE
    155154         CALL msg('Tracer <'//TRIM(var)//'> is missing => initialized to its parent isotope concentration.', modname)
    156          q(:,:,:,iq) = q(:,:,:,iqiso(iso_indnum(iq),iPhase))
     155         q(:,:,:,iq) = q(:,:,:,iqIsoPha(iName,iPhase))
    157156      END IF
    158157    !--------------------------------------------------------------------------------------------------------------------------
  • LMDZ6/trunk/libf/dyn3d/iniacademic.F90

    r4124 r4143  
    55
    66  USE filtreg_mod, ONLY: inifilr
    7   USE infotrac,    ONLY: nqtot, niso, niso_possibles, ok_iso_verif, tnat, alpha_ideal, &
    8                          iqiso, tracers, iso_indnum, iso_num
     7  USE infotrac,    ONLY: nqtot, niso, tnat, alpha_ideal, iqIsoPha, tracers
    98  USE control_mod, ONLY: day_step,planet_type
    109  use exner_hyb_m, only: exner_hyb
     
    282281              ! CRisi: init des isotopes
    283282              ! distill de Rayleigh très simplifiée
    284 !             iName    = tracers(iq)%iso_iName  ! (next commit)
    285               iName    = iso_num(iq)
     283              iName    = tracers(iq)%iso_iName
    286284              if (niso <= 0 .OR. iName <= 0) CYCLE
    287285              iPhase   = tracers(iq)%iso_iPhase
     
    290288                 q(:,:,iq) = q(:,:,iqParent)*tnat(iName)*(q(:,:,iqParent)/30.e-3)**(alpha_ideal(iName)-1.)
    291289              ELSE
    292                  q(:,:,iq) = q(:,:,iqiso(iso_indnum(iq),iPhase))
     290                 q(:,:,iq) = q(:,:,iqIsoPha(iName,iPhase))
    293291              END IF
    294292           enddo
     
    297295        endif ! of if (planet_type=="earth")
    298296
    299         if (ok_iso_verif) call check_isotopes_seq(q,1,ip1jmp1,'iniacademic_loc')
     297        call check_isotopes_seq(q,1,ip1jmp1,'iniacademic_loc')
    300298
    301299        ! add random perturbation to temperature
  • LMDZ6/trunk/libf/dyn3d/leapfrog.F

    r4120 r4143  
    1111      use IOIPSL
    1212#endif
    13       USE infotrac, ONLY: nqtot,ok_iso_verif
     13      USE infotrac, ONLY: nqtot, isoCheck
    1414      USE guide_mod, ONLY : guide_main
    1515      USE write_field, ONLY: writefield
     
    2626      USE temps_mod, ONLY: jD_ref,jH_ref,itaufin,day_ini,day_ref,
    2727     &                        start_time,dt
     28      USE strings_mod, ONLY: msg
    2829
    2930      IMPLICIT NONE
     
    237238      jH_cur = jH_cur - int(jH_cur)
    238239
    239         if (ok_iso_verif) then
    240            call check_isotopes_seq(q,ip1jmp1,'leapfrog 321')
    241         endif !if (ok_iso_verif) then
     240      call check_isotopes_seq(q,ip1jmp1,'leapfrog 321')
    242241
    243242#ifdef CPP_IOIPSL
     
    271270!      CALL filtreg ( finvmaold ,jjp1, llm, -2,2, .TRUE., 1 )
    272271
    273         if (ok_iso_verif) then
    274            call check_isotopes_seq(q,ip1jmp1,'leapfrog 400')
    275         endif !if (ok_iso_verif) then
     272      call check_isotopes_seq(q,ip1jmp1,'leapfrog 400')
    276273
    277274   2  CONTINUE ! Matsuno backward or leapfrog step begins here
     
    324321
    325322
    326         if (ok_iso_verif) then
    327            call check_isotopes_seq(q,ip1jmp1,'leapfrog 589')
    328         endif !if (ok_iso_verif) then
     323      call check_isotopes_seq(q,ip1jmp1,'leapfrog 589')
    329324
    330325c-----------------------------------------------------------------------
     
    345340c   -------------------------------------------------------------
    346341
    347         if (ok_iso_verif) then
    348            call check_isotopes_seq(q,ip1jmp1,
     342      call check_isotopes_seq(q,ip1jmp1,
    349343     &           'leapfrog 686: avant caladvtrac')
    350         endif !if (ok_iso_verif) then
    351344
    352345      IF( forward. OR . leapf )  THEN
     
    376369c   ----------------------------------
    377370
    378         if (ok_iso_verif) then
    379            write(*,*) 'leapfrog 720'
    380            call check_isotopes_seq(q,ip1jmp1,'leapfrog 756')
    381         endif !if (ok_iso_verif) then
     371       CALL msg('720', modname, isoCheck)
     372       call check_isotopes_seq(q,ip1jmp1,'leapfrog 756')
    382373       
    383374       CALL integrd ( nqtot,vcovm1,ucovm1,tetam1,psm1,massem1 ,
     
    385376!     $              finvmaold                                    )
    386377
    387        if (ok_iso_verif) then
    388           write(*,*) 'leapfrog 724'
    389            call check_isotopes_seq(q,ip1jmp1,'leapfrog 762')
    390         endif !if (ok_iso_verif) then
     378       CALL msg('724', modname, isoCheck)
     379       call check_isotopes_seq(q,ip1jmp1,'leapfrog 762')
    391380
    392381c .P.Le Van (26/04/94  ajout de  finvpold dans l'appel d'integrd)
     
    552541        CALL massdair(p,masse)
    553542
    554         if (ok_iso_verif) then
    555            call check_isotopes_seq(q,ip1jmp1,'leapfrog 1196')
    556         endif !if (ok_iso_verif) then
     543        call check_isotopes_seq(q,ip1jmp1,'leapfrog 1196')
    557544
    558545c-----------------------------------------------------------------------
     
    639626c   preparation du pas d'integration suivant  ......
    640627
    641         if (ok_iso_verif) then
    642            call check_isotopes_seq(q,ip1jmp1,'leapfrog 1509')
    643         endif !if (ok_iso_verif) then
     628      call check_isotopes_seq(q,ip1jmp1,'leapfrog 1509')
    644629
    645630      IF ( .NOT.purmats ) THEN
     
    703688            ENDIF ! of IF((MOD(itau,iperiod).EQ.0).OR.(itau.EQ.itaufin))
    704689
    705         if (ok_iso_verif) then
    706            call check_isotopes_seq(q,ip1jmp1,'leapfrog 1584')
    707         endif !if (ok_iso_verif) then
     690            call check_isotopes_seq(q,ip1jmp1,'leapfrog 1584')
    708691
    709692c-----------------------------------------------------------------------
     
    790773      ELSE ! of IF (.not.purmats)
    791774
    792         if (ok_iso_verif) then
    793            call check_isotopes_seq(q,ip1jmp1,'leapfrog 1664')
    794         endif !if (ok_iso_verif) then
     775            call check_isotopes_seq(q,ip1jmp1,'leapfrog 1664')
    795776
    796777c       ........................................................
     
    817798            ELSE ! of IF(forward) i.e. backward step
    818799 
    819         if (ok_iso_verif) then
    820            call check_isotopes_seq(q,ip1jmp1,'leapfrog 1698')
    821         endif !if (ok_iso_verif) then 
     800              call check_isotopes_seq(q,ip1jmp1,'leapfrog 1698')
    822801
    823802              IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin) THEN
  • LMDZ6/trunk/libf/dyn3d/qminimum.F

    r4124 r4143  
    44      SUBROUTINE qminimum( q,nqtot,deltap )
    55
    6       USE infotrac, ONLY: niso, ntraciso,iqiso,ok_iso_verif
     6      USE infotrac, ONLY: niso, ntiso,iqIsoPha, tracers
     7      USE strings_mod, ONLY: strIdx
     8      USE readTracFiles_mod, ONLY: addPhase
    79      IMPLICIT none
    810c
     
    1618      REAL q(ip1jmp1,llm,nqtot), deltap(ip1jmp1,llm)
    1719c
    18       INTEGER iq_vap, iq_liq
    19       PARAMETER ( iq_vap = 1 ) ! indice pour l'eau vapeur
    20       PARAMETER ( iq_liq = 2 ) ! indice pour l'eau liquide
    21       REAL seuil_vap, seuil_liq
    22       PARAMETER ( seuil_vap = 1.0e-10 ) ! seuil pour l'eau vapeur
    23       PARAMETER ( seuil_liq = 1.0e-11 ) ! seuil pour l'eau liquide
     20      LOGICAL, SAVE :: first=.TRUE.
     21      INTEGER, SAVE :: iq_vap, iq_liq        ! indices pour l'eau vapeur/liquide
     22      REAL, PARAMETER :: seuil_vap = 1.0e-10 ! seuil pour l'eau vapeur
     23      REAL, PARAMETER :: seuil_liq = 1.0e-11 ! seuil pour l'eau liquide
    2424c
    2525c  NB. ....( Il est souhaitable mais non obligatoire que les valeurs des
     
    4343      !INTEGER nb_pump
    4444      INTEGER ixt
     45
     46      IF(first) THEN
     47         iq_vap = strIdx(tracers(:)%name, addPhase('H2O', 'g'))
     48         iq_liq = strIdx(tracers(:)%name, addPhase('H2O', 'l'))
     49         first = .FALSE.
     50      END IF
    4551c
    4652c Quand l'eau liquide est trop petite (ou negative), on prend
     
    4955c
    5056
    51         if (ok_iso_verif) then
    52            call check_isotopes_seq(q,ip1jmp1,'qminimum 52')   
    53         endif !if (ok_iso_verif) then     
     57      call check_isotopes_seq(q,ip1jmp1,'qminimum 52')   
    5458
    5559      zx_defau_diag(:,:,:)=0.0
     
    127131          if (zx_defau_diag(i,k,iq_vap).gt.0.0) then             
    128132              ! on ajoute la vapeur en k             
    129               do ixt=1,ntraciso
    130                q(i,k,iqiso(ixt,iq_vap))=q(i,k,iqiso(ixt,iq_vap))
    131      :              +zx_defau_diag(i,k,iq_vap)
    132      :              *q(i,k-1,iqiso(ixt,iq_vap))/q_follow(i,k-1,iq_vap)
     133              do ixt=1,ntiso
     134               q(i,k,iqIsoPha(ixt,iq_vap))=q(i,k,iqIsoPha(ixt,iq_vap))
     135     :           +zx_defau_diag(i,k,iq_vap)
     136     :           *q(i,k-1,iqIsoPha(ixt,iq_vap))/q_follow(i,k-1,iq_vap)
    133137               
    134138              ! et on la retranche en k-1
    135                q(i,k-1,iqiso(ixt,iq_vap))=q(i,k-1,iqiso(ixt,iq_vap))
     139               q(i,k-1,iqIsoPha(ixt,iq_vap))=
     140     :            q(i,k-1,iqIsoPha(ixt,iq_vap))
    136141     :              -zx_defau_diag(i,k,iq_vap)
    137142     :              *deltap(i,k)/deltap(i,k-1)
    138      :              *q(i,k-1,iqiso(ixt,iq_vap))/q_follow(i,k-1,iq_vap)
     143     :              *q(i,k-1,iqIsoPha(ixt,iq_vap))
     144     :              /q_follow(i,k-1,iq_vap)
    139145
    140146              enddo !do ixt=1,niso
     
    148154       enddo !do k=2,llm
    149155
    150         if (ok_iso_verif) then     
    151            call check_isotopes_seq(q,ip1jmp1,'qminimum 168')
    152         endif !if (ok_iso_verif) then
     156       call check_isotopes_seq(q,ip1jmp1,'qminimum 168')
    153157       
    154158     
     
    160164
    161165              ! on ajoute eau liquide en k en k             
    162               do ixt=1,ntraciso
    163                q(i,k,iqiso(ixt,iq_liq))=q(i,k,iqiso(ixt,iq_liq))
     166              do ixt=1,ntiso
     167               q(i,k,iqIsoPha(ixt,iq_liq))=q(i,k,iqIsoPha(ixt,iq_liq))
    164168     :              +zx_defau_diag(i,k,iq_liq)
    165      :              *q(i,k,iqiso(ixt,iq_vap))/q_follow(i,k,iq_vap)
     169     :              *q(i,k,iqIsoPha(ixt,iq_vap))/q_follow(i,k,iq_vap)
    166170              ! et on la retranche à la vapeur en k
    167                q(i,k,iqiso(ixt,iq_vap))=q(i,k,iqiso(ixt,iq_vap))
     171               q(i,k,iqIsoPha(ixt,iq_vap))=q(i,k,iqIsoPha(ixt,iq_vap))
    168172     :              -zx_defau_diag(i,k,iq_liq)
    169      :              *q(i,k,iqiso(ixt,iq_vap))/q_follow(i,k,iq_vap)   
     173     :              *q(i,k,iqIsoPha(ixt,iq_vap))/q_follow(i,k,iq_vap)   
    170174              enddo !do ixt=1,niso
    171175              q_follow(i,k,iq_liq)=   q_follow(i,k,iq_liq)
     
    177181       enddo !do k=2,llm 
    178182
    179         if (ok_iso_verif) then
    180            call check_isotopes_seq(q,ip1jmp1,'qminimum 197')
    181         endif !if (ok_iso_verif) then
     183       call check_isotopes_seq(q,ip1jmp1,'qminimum 197')
    182184
    183185      endif !if (niso > 0) then
  • LMDZ6/trunk/libf/dyn3d/vlsplt.F

    r4064 r4143  
    125125      RECURSIVE SUBROUTINE vlx(q,pente_max,masse,u_m,iq)
    126126      USE infotrac, ONLY : nqtot,tracers, ! CRisi
    127      &                     qperemin,masseqmin,ratiomin ! MVals et CRisi
     127     &                     min_qParent,min_qMass,min_ratio ! MVals et CRisi
    128128
    129129c     Auteurs:   P.Le Van, F.Hourdin, F.Forget
     
    428428            !Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)
    429429            !Mvals: veiller a ce qu'on n'ait pas de denominateur nul
    430             masseq(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),masseqmin)
    431             if (q(ij,l,iq).gt.qperemin) then
     430            masseq(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),min_qMass)
     431            if (q(ij,l,iq).gt.min_qParent) then
    432432              Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)
    433433            else
    434               Ratio(ij,l,iq2)=ratiomin
     434              Ratio(ij,l,iq2)=min_ratio
    435435            endif
    436436          enddo   
     
    449449         DO ij=iip2+1,ip1jm
    450450            !MVals: veiller a ce qu'on ait pas de denominateur nul
    451             new_m=max(masse(ij,l,iq)+u_m(ij-1,l)-u_m(ij,l),masseqmin)
     451            new_m=max(masse(ij,l,iq)+u_m(ij-1,l)-u_m(ij,l),min_qMass)
    452452            q(ij,l,iq)=(q(ij,l,iq)*masse(ij,l,iq)+
    453453     &      u_mq(ij-1,l)-u_mq(ij,l))
     
    485485      RECURSIVE SUBROUTINE vly(q,pente_max,masse,masse_adv_v,iq)
    486486      USE infotrac, ONLY : nqtot,tracers, ! CRisi
    487      &                     qperemin,masseqmin,ratiomin ! MVals et CRisi
     487     &                     min_qParent,min_qMass,min_ratio ! MVals et CRisi
    488488c
    489489c     Auteurs:   P.Le Van, F.Hourdin, F.Forget
     
    752752            !Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)     
    753753            !MVals: veiller a ce qu'on n'ait pas de denominateur nul
    754             masseq(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),masseqmin)
    755             if (q(ij,l,iq).gt.qperemin) then
     754            masseq(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),min_qMass)
     755            if (q(ij,l,iq).gt.min_qParent) then
    756756              Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)
    757757            else
    758               Ratio(ij,l,iq2)=ratiomin
     758              Ratio(ij,l,iq2)=min_ratio
    759759            endif
    760760          enddo   
     
    848848      RECURSIVE SUBROUTINE vlz(q,pente_max,masse,w,iq)
    849849      USE infotrac, ONLY : nqtot,tracers, ! CRisi
    850      &                     qperemin,masseqmin,ratiomin ! MVals et CRisi
     850     &                     min_qParent,min_qMass,min_ratio ! MVals et CRisi
    851851c
    852852c     Auteurs:   P.Le Van, F.Hourdin, F.Forget
     
    977977            !Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)       
    978978            !MVals: veiller a ce qu'on n'ait pas de denominateur nul
    979             masseq(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),masseqmin)
    980             if (q(ij,l,iq).gt.qperemin) then
     979            masseq(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),min_qMass)
     980            if (q(ij,l,iq).gt.min_qParent) then
    981981              Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)
    982982            else
    983               Ratio(ij,l,iq2)=ratiomin
     983              Ratio(ij,l,iq2)=min_ratio
    984984            endif     
    985985          enddo   
Note: See TracChangeset for help on using the changeset viewer.