source: LMDZ6/trunk/libf/dyn3dmem/qminimum_loc.F @ 5084

Last change on this file since 5084 was 5001, checked in by dcugnet, 5 months ago
  • strings_mod:
    • remove "test()" function (was not very clear)
    • modifications of the "str2bool" function: result is O/1 for .FALSE./.TRUE. and -1 if the string was not a boolean.
    • more general "find()" function (for several numerical types)
    • more general "cat()" function (can append a 2D array with vectors, 1D arrays with scalars)
    • few simplifications (in "strParse") and minor changes
  • readTracFiles_mod:
    • remove internal usage of direct keys ("%" symbol) in favor of the "getKey" function. => moving toward a totally generic tracers derived type.
    • improve the internal management of the error return value "lerr".
    • remove "fGetKey", "fGetKeys", "setDirectKeys" functions
    • new functions to add/remove a phase: "addPhase", "delPhase"
    • more general "addKey(key[(:)], val[(:)], ky(:), [lOverWrite])" function: . input argument "val" can be string/integer/real/logical . (key, val, ky ): add the <key> =<val> pair to ky . (key, val(:), ky(:)): add the <key> =<val(i)> pair to ky(i) for 1<=i<=SIZE(ky) . (key(:), val(:), ky(:)): add the <key(i)>=<val(i)> pair to ky(i) for 1<=i<=SIZE(ky)
    • more general "getKey(key[(:)], val[(:)], itr [, ky(:)][, nam(:)][, def][, lDisp])" (tracer index version)

and "getKey(key[(:)], val[(:)], tname[, ky(:)]. [, def][, lDisp])" (tracer name version) functions:

. output argument "val" can be string/integer/real/logical
. if present, the default value <def> is retained if the corresponding key was not found.
. get values from "ky(:)" if present, otherwise from internal database "tracers(:)" or "isotope ».
. if "keyn" is a vector, try with each element in indices order until a value is found
. (key[(:)], val, itr/tname[,ky(:)][, ...]): get the value <val> of tracer nr. itr or named "tname"
. (key[(:)], val(:), itr/tname[,ky(:)][, ...]): same + parsing of the value with « , », then storage in <val(:)>
. (key[(:)], val(:)[, ky(:)][, nam(:)][, ...]): same for all tracers (optional names list <nam(:)>) of database.
. (key[(:)], val(:), tname(:)[, ky(:)][, ...]): same for the tracers named « tnames(:)"

  • more general "dispTraSection" function
  • much simplified "indexUpdate" function ; "ancestor*" and "idxAncestor" functions are removed.
  • "readIsotopesFile" is renamed to "processIsotopes" for more clarity
  • cosmetic changes
  • fix for isotopes: iq_val and iq_liq are usable for "q" only, not for "q_follow" and "zx_defau_diag" => use hardcoded indices (1 for vapor and 2 for liquid) for these variables
  • 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: 10.2 KB
Line 
1!
2!     $Id: qminimum_loc.F 5001 2024-07-01 09:25:05Z fairhead $
3!
4      SUBROUTINE qminimum_loc( q,nqtot,deltap )
5      USE parallel_lmdz
6      USE infotrac, ONLY: niso, ntiso, iqIsoPha, tracers,
7     &                    isoCheck, min_qParent
8      USE strings_mod, ONLY: strIdx
9      USE readTracFiles_mod, ONLY: addPhase
10      IMPLICIT none
11c
12c  -- Objet : Traiter les valeurs trop petites (meme negatives)
13c             pour l'eau vapeur et l'eau liquide
14c
15      include "dimensions.h"
16      include "paramet.h"
17      include "iniprint.h"
18c
19      INTEGER nqtot ! CRisi: on remplace nq par nqtot
20      REAL q(ijb_u:ije_u,llm,nqtot), deltap(ijb_u:ije_u,llm)
21c
22      LOGICAL, SAVE :: first=.TRUE.
23      INTEGER, SAVE :: iq_vap, iq_liq        ! indices pour l'eau vapeur/liquide
24c$OMP THREADPRIVATE(iq_vap, iq_liq, first)
25      REAL, PARAMETER :: seuil_vap = 1.0e-10 ! seuil pour l'eau vapeur
26      REAL, PARAMETER :: seuil_liq = 1.0e-11 ! seuil pour l'eau liquide
27c
28c  NB. ....( Il est souhaitable mais non obligatoire que les valeurs des
29c            parametres seuil_vap, seuil_liq soient pareilles a celles
30c            qui  sont utilisees dans la routine    ADDFI       )
31c     .................................................................
32c
33cDC iq_val and iq_liq are usable for q only, NOT for q_follow
34c   and zx_defau_diag (crash if iq_val/liq==3) => vapor/liquid
35c   water at hardcoded indices 1/2 in these variables
36      INTEGER i, k, iq
37      REAL zx_defau, zx_abc, zx_pump(ijb_u:ije_u), pompe
38
39      real zx_defau_diag(ijb_u:ije_u,llm,2)
40      real q_follow(ijb_u:ije_u,llm,2)
41c
42      REAL SSUM
43      EXTERNAL SSUM
44c
45      INTEGER imprim
46      SAVE imprim
47      DATA imprim /0/
48c$OMP THREADPRIVATE(imprim)
49      INTEGER ijb,ije
50      INTEGER Index_pump(ij_end-ij_begin+1)
51      INTEGER nb_pump
52      INTEGER ixt
53      INTEGER iso_verif_noNaN_nostop
54
55c$OMP BARRIER
56
57      !write(lunout,*) 'qminimum 52: entree'
58      IF(first) THEN
59         iq_vap = strIdx(tracers(:)%name, addPhase('H2O', 'g'))
60         iq_liq = strIdx(tracers(:)%name, addPhase('H2O', 'l'))
61         first = .FALSE.
62      END IF
63c
64c Quand l'eau liquide est trop petite (ou negative), on prend
65c l'eau vapeur de la meme couche et la convertit en eau liquide
66c (sans changer la temperature !)
67c
68
69      call check_isotopes(q,ij_begin,ij_end,'qminimum 52')   
70
71      ijb=ij_begin
72      ije=ij_end
73
74      DO k = 1, llm
75c$OMP DO SCHEDULE(STATIC)
76        DO i = ijb, ije
77          zx_defau_diag(i,k,1)=0.0
78          zx_defau_diag(i,k,2)=0.0
79          q_follow(i,k,1)=q(i,k,iq_vap)
80          q_follow(i,k,2)=q(i,k,iq_liq)
81        ENDDO
82c$OMP END DO NOWAIT
83      ENDDO
84
85      !write(lunout,*) 'qminimum 57'
86      DO k = 1, llm
87c$OMP DO SCHEDULE(STATIC)       
88        DO i = ijb, ije
89          if (seuil_liq - q(i,k,iq_liq) .gt. 0.d0 ) then
90
91            if (niso > 0) zx_defau_diag(i,k,2)=AMAX1
92     :               ( seuil_liq - q(i,k,iq_liq), 0.0 )
93
94            q(i,k,iq_vap) = q(i,k,iq_vap) + q(i,k,iq_liq) - seuil_liq
95            q(i,k,iq_liq) = seuil_liq
96          endif
97        END DO
98c$OMP END DO NOWAIT
99      END DO
100
101c
102c Quand l'eau vapeur est trop faible (ou negative), on complete
103c le defaut en prennant de l'eau vapeur de la couche au-dessous.
104c
105      !write(lunout,*) 'qminimum 81'
106      DO k = llm, 2, -1
107ccc      zx_abc = dpres(k) / dpres(k-1)
108c$OMP DO SCHEDULE(STATIC)
109        DO i = ijb, ije
110
111          if ( seuil_vap - q(i,k,iq_vap) .gt. 0.d0 ) then
112
113            if (niso > 0) zx_defau_diag(i,k,1)
114     &           = AMAX1( seuil_vap - q(i,k,iq_vap), 0.0 )
115
116            q(i,k-1,iq_vap) = q(i,k-1,iq_vap) - (seuil_vap
117     &           -q(i,k,iq_vap)) * deltap(i,k)/deltap(i,k-1)
118            q(i,k,iq_vap)   =  seuil_vap 
119
120          endif
121        ENDDO
122c$OMP END DO NOWAIT
123      ENDDO
124
125c
126c Quand il s'agit de la premiere couche au-dessus du sol, on
127c doit imprimer un message d'avertissement (saturation possible).
128c
129      !write(lunout,*) 'qminimum 106'
130      nb_pump=0
131c$OMP DO SCHEDULE(STATIC)
132      DO i = ijb, ije
133         zx_pump(i) = AMAX1( 0.0, seuil_vap - q(i,1,iq_vap) )
134         q(i,1,iq_vap)  = AMAX1( q(i,1,iq_vap), seuil_vap )
135         IF (zx_pump(i) > 0.0) THEN
136            nb_pump = nb_pump+1
137            Index_pump(nb_pump)=i
138         ENDIF
139      ENDDO
140c$OMP END DO NOWAIT
141!      pompe = SSUM(ije-ijb+1,zx_pump(ijb),1)
142
143      IF (imprim.LE.100 .AND. nb_pump .GT. 0 ) THEN
144         PRINT *, 'ATT!:on pompe de l eau au sol'
145         DO i = 1, nb_pump
146               imprim = imprim + 1
147               PRINT*,'  en ',index_pump(i),zx_pump(index_pump(i))
148         ENDDO
149      ENDIF
150
151      !write(lunout,*) 'qminimum 128'
152      if (niso > 0) then
153              !write(lunout,*) 'qminimum 140'
154      ! CRisi: traiter de même les traceurs d'eau
155      ! Mais il faut les prendre à l'envers pour essayer de conserver la
156      ! masse.
157      ! 1) pompage dans le sol 
158      ! On suppose que ce pompage se fait sans isotopes -> on ne modifie
159      ! rien ici et on croise les doigts pour que ça ne soit pas trop
160      ! génant
161      ! en fait, si, c'est genant quand les isotopes doivent eux même transporter des
162      ! traceurs -> apporter aussi un peu d'isotopes... Combien?
163      ! Essayer tnat/2 = -500 permil? C'est déjà mieux que -1000
164      ! permil...
165      ! pb: que faire pour les traceurs?
166c$OMP DO SCHEDULE(STATIC)     
167      DO i = ijb, ije
168        if (zx_pump(i).gt.0.0) then
169          q_follow(i,1,1)=q_follow(i,1,1)+zx_pump(i)
170        endif !if (zx_pump(i).gt.0.0) then
171      enddo !DO i = ijb, ije 
172c$OMP END DO NOWAIT
173
174      ! 2) transfert de vap vers les couches plus hautes
175      !write(lunout,*) 'qminimum 158'
176      do k=2,llm
177c$OMP DO SCHEDULE(STATIC)     
178        DO i = ijb, ije
179          if (zx_defau_diag(i,k,1).gt.0.0) then             
180              ! on ajoute la vapeur en k     
181!              write(lunout,*) 'i,k,q_follow(i,k-1,ivap)=',
182!     :                 i,k,q_follow(i,k-1,1)         
183              if (q_follow(i,k-1,1).lt.min_qParent) then
184                write(lunout,*) 'tmp qmin: on stoppe'
185                write(lunout,*) 'zx_pump(i)=',zx_pump(i)
186                write(lunout,*) 'q_follow(i,:,ivap)=',
187     :                   q_follow(i,:,1)
188                write(lunout,*) 'k=',k
189                call abort_gcm("qminimum","not enough vapor",1)
190              endif 
191            do ixt=1,ntiso
192!                write(lunout,*) 'qmin 168: ixt=',ixt
193!                write(lunout,*) 'q(i,k,iqIsoPha(ixt,iq_vap))=',
194!     :             q(i,k,iqIsoPha(ixt,iq_vap))
195!                write(lunout,*) 'zx_defau_diag(i,k,ivap)=',
196!     :                  zx_defau_diag(i,k,1)
197!                write(lunout,*) 'q(i,k-1,iqIsoPha(ixt,iq_vap))=',
198!     :                   q(i,k-1,iqIsoPha(ixt,iq_vap))     
199
200               q(i,k,iqIsoPha(ixt,iq_vap))=q(i,k,iqIsoPha(ixt,iq_vap))
201     :           +zx_defau_diag(i,k,1)
202     :           *q(i,k-1,iqIsoPha(ixt,iq_vap))/q_follow(i,k-1,1)
203               
204              if (isoCheck) then
205                if(iso_verif_noNaN_nostop(q(i,k,iqIsoPha(ixt,iq_vap)),
206     :                   'qminimum 155').eq.1) then
207                   write(*,*) 'i,k,ixt=',i,k,ixt
208                   write(*,*) 'q_follow(i,k-1,ivap)=',
209     :                   q_follow(i,k-1,1)
210                   write(*,*) 'q(i,k,iqIsoPha(ixt,iq_vap))=',
211     :                   q(i,k,iqIsoPha(ixt,iq_vap))
212                   write(*,*) 'zx_defau_diag(i,k,ivap)=',
213     :                   zx_defau_diag(i,k,1)
214                   write(*,*) 'q(i,k-1,iqIsoPha(ixt,iq_vap))=',
215     :                   q(i,k-1,iqIsoPha(ixt,iq_vap))
216                CALL abort_gcm("qminimum_loc","stopped",1)
217                endif
218              endif
219
220              ! et on la retranche en k-1
221               q(i,k-1,iqIsoPha(ixt,iq_vap)) =
222     :            q(i,k-1,iqIsoPha(ixt,iq_vap))
223     :              -zx_defau_diag(i,k,1)
224     :              *deltap(i,k)/deltap(i,k-1)
225     :              *q(i,k-1,iqIsoPha(ixt,iq_vap))
226     :              /q_follow(i,k-1,1)
227
228               if (isoCheck) then
229                if (iso_verif_noNaN_nostop(
230     :              q(i,k-1,iqIsoPha(ixt,iq_vap)),
231     :                   'qminimum 175').eq.1) then
232                   write(*,*) 'k,i,ixt=',k,i,ixt
233                   write(*,*) 'q_follow(i,k-1,ivap)=',
234     :                   q_follow(i,k-1,1)
235                   write(*,*) 'q(i,k,iqIsoPha(ixt,iq_vap))=',
236     :                   q(i,k,iqIsoPha(ixt,iq_vap))
237                   write(*,*) 'zx_defau_diag(i,k,ivap)=',
238     :                   zx_defau_diag(i,k,1)
239                   write(*,*) 'q(i,k-1,iqIsoPha(ixt,iq_vap))=',
240     :                   q(i,k-1,iqIsoPha(ixt,iq_vap))
241                   CALL abort_gcm("qminimum_loc","stopped",1)
242                endif
243              endif
244
245              enddo !do ixt=1,niso
246              q_follow(i,k,1)=   q_follow(i,k,1)
247     :               +zx_defau_diag(i,k,1)
248              q_follow(i,k-1,1)=   q_follow(i,k-1,1)
249     :               -zx_defau_diag(i,k,1)
250     :              *deltap(i,k)/deltap(i,k-1)
251          endif !if (zx_defau_diag(i,k,1).gt.0.0) then
252        enddo !DO i = 1, ip1jmp1       
253c$OMP END DO NOWAIT
254        enddo !do k=2,llm
255
256        call check_isotopes(q,ijb,ije,'qminimum 168')
257       
258     
259        ! 3) transfert d'eau de la vapeur au liquide
260        !write(*,*) 'qminimum 164'
261        do k=1,llm
262c$OMP DO SCHEDULE(STATIC)
263        DO i = ijb, ije
264          if (zx_defau_diag(i,k,2).gt.0.0) then
265
266              ! on ajoute eau liquide en k en k             
267              do ixt=1,ntiso
268               q(i,k,iqIsoPha(ixt,iq_liq))=q(i,k,iqIsoPha(ixt,iq_liq))
269     :              +zx_defau_diag(i,k,2)
270     :              *q(i,k,iqIsoPha(ixt,iq_vap))/q_follow(i,k,1)
271              ! et on la retranche à la vapeur en k
272               q(i,k,iqIsoPha(ixt,iq_vap))=q(i,k,iqIsoPha(ixt,iq_vap))
273     :              -zx_defau_diag(i,k,2)
274     :              *q(i,k,iqIsoPha(ixt,iq_vap))/q_follow(i,k,1)   
275              enddo !do ixt=1,niso
276              q_follow(i,k,2)=   q_follow(i,k,2)
277     :               +zx_defau_diag(i,k,2)
278              q_follow(i,k,1)=   q_follow(i,k,1)
279     :               -zx_defau_diag(i,k,2)
280          endif !if (zx_defau_diag(i,k,1).gt.0.0) then
281        enddo !DO i = ijb, ije
282c$OMP END DO NOWAIT       
283       enddo !do k=2,llm 
284
285       call check_isotopes(q,ijb,ije,'qminimum 197')
286
287      endif !if (niso > 0) then
288      !write(*,*) 'qminimum 188'
289c$OMP BARRIER
290
291c
292      RETURN
293      END
Note: See TracBrowser for help on using the repository browser.