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

Last change on this file since 4124 was 4124, checked in by dcugnet, 2 years ago

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