source: LMDZ6/trunk/libf/dyn3dmem/qminimum_loc.f90 @ 5407

Last change on this file since 5407 was 5285, checked in by abarral, 7 weeks ago

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