source: LMDZ6/branches/Amaury_dev/libf/dyn3dmem/vlspltgen_loc.F90 @ 5133

Last change on this file since 5133 was 5128, checked in by abarral, 5 months ago

Correct bug in vlspltqs_loc.f90 from r2270 where we call SSUM with incorrect arguments.
Merge the three different versions of abort_gcm into one
Fix seq, para 3D compilation broken from r5107 onwards
(lint) usual + Remove uneeded fixed-form continuations
  • 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
File size: 11.4 KB
Line 
1! $Header$
2
3SUBROUTINE vlspltgen_loc(q, pente_max, masse, w, pbaru, pbarv, &
4        pdt, p, pk, teta)
5
6  !
7  ! Auteurs:   P.Le Van, F.Hourdin, F.Forget, F.Codron
8  !
9  !    ********************************************************************
10  !      Schema  d'advection " pseudo amont " .
11  !  + test sur humidite specifique: Q advecte< Qsat aval
12  !               (F. Codron, 10/99)
13  !    ********************************************************************
14  ! q,pbaru,pbarv,w sont des arguments d'entree  pour le s-pg ....
15  !
16  ! pente_max facteur de limitation des pentes: 2 en general
17  !                                            0 pour un schema amont
18  ! pbaru,pbarv,w flux de masse en u ,v ,w
19  ! pdt pas de temps
20  !
21  ! teta temperature potentielle, p pression aux interfaces,
22  ! pk exner au milieu des couches necessaire pour calculer Qsat
23  !   --------------------------------------------------------------------
24  USE parallel_lmdz
25  USE mod_hallo
26  USE write_field_loc, ONLY: WriteField_u, WriteField_v
27  USE lmdz_vampir
28  ! CRisi: on rajoute variables utiles d'infotrac
29  USE infotrac, ONLY: nqtot, tracers, isoCheck
30  USE vlspltgen_mod
31  USE comconst_mod, ONLY: cpp
32  USE logic_mod, ONLY: adv_qsat_liq
33
34
35  IMPLICIT NONE
36
37  !
38  include "dimensions.h"
39  include "paramet.h"
40
41  !
42  !   Arguments:
43  !   ----------
44  REAL :: masse(ijb_u:ije_u, llm), pente_max
45  REAL :: pbaru(ijb_u:ije_u, llm), pbarv(ijb_v:ije_v, llm)
46  REAL :: q(ijb_u:ije_u, llm, nqtot)
47  REAL :: w(ijb_u:ije_u, llm), pdt
48  REAL :: p(ijb_u:ije_u, llmp1), teta(ijb_u:ije_u, llm)
49  REAL :: pk(ijb_u:ije_u, llm)
50  !
51  !  Local
52  !   ---------
53  !
54  INTEGER :: ij, l
55  !
56  REAL :: zzpbar, zzw
57
58  REAL :: qmin, qmax
59  DATA qmin, qmax/0., 1.e33/
60
61  !--pour rapport de melange saturant--
62
63  REAL :: rtt, retv, r2es, r3les, r3ies, r4les, r4ies, play
64  REAL :: ptarg, pdelarg, foeew, zdelta
65  REAL :: tempe(ijb_u:ije_u)
66  INTEGER :: ijb, ije, iq, iq2, ifils
67  LOGICAL, SAVE :: firstcall = .TRUE.
68  !$OMP THREADPRIVATE(firstcall)
69  type(request), SAVE :: MyRequest1
70  !$OMP THREADPRIVATE(MyRequest1)
71  type(request), SAVE :: MyRequest2
72  !$OMP THREADPRIVATE(MyRequest2)
73  !    fonction psat(T)
74
75  FOEEW (PTARG, PDELARG) = EXP (&
76          (R3LES * (1. - PDELARG) + R3IES * PDELARG) * (PTARG - RTT) &
77                  / (PTARG - (R4LES * (1. - PDELARG) + R4IES * PDELARG)))
78
79  r2es = 380.11733
80  r3les = 17.269
81  r3ies = 21.875
82  r4les = 35.86
83  r4ies = 7.66
84  retv = 0.6077667
85  rtt = 273.16
86
87  ! Allocate variables depending on dynamic variable nqtot
88
89  IF (firstcall) THEN
90    firstcall = .FALSE.
91  END IF
92  !-- Calcul de Qsat en chaque point
93  !-- approximation: au milieu des couches play(l)=(p(l)+p(l+1))/2
94  !   pour eviter une exponentielle.
95
96  CALL SetTag(MyRequest1, 100)
97  CALL SetTag(MyRequest2, 101)
98
99  ijb = ij_begin - iip1
100  ije = ij_end + iip1
101  IF (pole_nord) ijb = ij_begin
102  IF (pole_sud) ije = ij_end
103
104  !$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
105  DO l = 1, llm
106    DO ij = ijb, ije
107      tempe(ij) = teta(ij, l) * pk(ij, l) / cpp
108    ENDDO
109    DO ij = ijb, ije
110      IF (adv_qsat_liq) THEN
111        zdelta = 0.
112      ELSE
113        zdelta = MAX(0., SIGN(1., rtt - tempe(ij)))
114      ENDIF
115      play = 0.5 * (p(ij, l) + p(ij, l + 1))
116      qsat(ij, l) = MIN(0.5, r2es * FOEEW(tempe(ij), zdelta) / play)
117      qsat(ij, l) = qsat(ij, l) / (1. - retv * qsat(ij, l))
118    ENDDO
119  ENDDO
120  !$OMP END DO NOWAIT
121  ! PRINT*,'Debut vlsplt version debug sans vlyqs'
122
123  zzpbar = 0.5 * pdt
124  zzw = pdt
125
126  ijb = ij_begin
127  ije = ij_end
128  IF (pole_nord) ijb = ijb + iip1
129  IF (pole_sud)  ije = ije - iip1
130
131  !$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
132  DO l = 1, llm
133    DO ij = ijb, ije
134      mu(ij, l) = pbaru(ij, l) * zzpbar
135    ENDDO
136  ENDDO
137  !$OMP END DO NOWAIT
138
139  ijb = ij_begin - iip1
140  ije = ij_end
141  IF (pole_nord) ijb = ij_begin
142  IF (pole_sud)  ije = ij_end - iip1
143
144  !$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
145  DO l = 1, llm
146    DO ij = ijb, ije
147      mv(ij, l) = pbarv(ij, l) * zzpbar
148    ENDDO
149  ENDDO
150  !$OMP END DO NOWAIT
151
152  ijb = ij_begin
153  ije = ij_end
154
155  DO iq = 1, nqtot
156    !$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
157    DO l = 1, llm
158      DO ij = ijb, ije
159        mw(ij, l, iq) = w(ij, l) * zzw
160      ENDDO
161    ENDDO
162    !$OMP END DO NOWAIT
163  ENDDO
164
165  DO iq = 1, nqtot
166    !$OMP MASTER
167    DO ij = ijb, ije
168      mw(ij, llm + 1, iq) = 0.
169    ENDDO
170    !$OMP END MASTER
171  ENDDO
172
173  ijb = ij_begin
174  ije = ij_end
175
176  DO iq = 1, nqtot
177    !$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
178    DO l = 1, llm
179      zq(ijb:ije, l, iq) = q(ijb:ije, l, iq)
180      zm(ijb:ije, l, iq) = masse(ijb:ije, l)
181    ENDDO
182    !$OMP END DO NOWAIT
183  ENDDO
184
185  ! verif temporaire
186  ijb = ij_begin
187  ije = ij_end
188  CALL check_isotopes(zq, ijb, ije, 'vlspltgen_loc 191')
189
190  !$OMP BARRIER
191  DO iq = 1, nqtot
192    ! CRisi: on ne boucle que sur les pères= ceux qui sont transportés directement par l'air
193    IF(tracers(iq)%parent /= 'air') CYCLE
194    !WRITE(*,*) 'vlspltgen 192: iq,iadv=',iq,tracers(iq)%iadv
195    SELECT CASE(tracers(iq)%iadv)
196    CASE(0); CYCLE
197    CASE(10)
198      CALL vlx_loc(zq, pente_max, zm, mu, &
199              ij_begin, ij_end, iq)
200
201      !$OMP MASTER
202      CALL VTb(VTHallo)
203      !$OMP END MASTER
204      CALL Register_Hallo_u(zq(:, :, iq), llm, 2, 2, 2, 2, MyRequest1)
205      CALL Register_Hallo_u(zm(:, :, iq), llm, 1, 1, 1, 1, MyRequest1)
206      ! CRisi
207      do ifils = 1, tracers(iq)%nqDescen
208        iq2 = tracers(iq)%iqDescen(ifils)
209        CALL Register_Hallo_u(zq(:, :, iq2), llm, 2, 2, 2, 2, MyRequest1)
210        CALL Register_Hallo_u(zm(:, :, iq2), llm, 1, 1, 1, 1, MyRequest1)
211      enddo
212
213      !$OMP MASTER
214      CALL VTe(VTHallo)
215      !$OMP END MASTER
216    CASE(14)
217      CALL vlxqs_loc(zq, pente_max, zm, mu, &
218              qsat, ij_begin, ij_end, iq)
219
220      !$OMP MASTER
221      CALL VTb(VTHallo)
222      !$OMP END MASTER
223
224      CALL Register_Hallo_u(zq(:, :, iq), llm, 2, 2, 2, 2, MyRequest1)
225      CALL Register_Hallo_u(zm(:, :, iq), llm, 1, 1, 1, 1, MyRequest1)
226      do ifils = 1, tracers(iq)%nqDescen
227        iq2 = tracers(iq)%iqDescen(ifils)
228        CALL Register_Hallo_u(zq(:, :, iq2), llm, 2, 2, 2, 2, MyRequest1)
229        CALL Register_Hallo_u(zm(:, :, iq2), llm, 1, 1, 1, 1, MyRequest1)
230      enddo
231
232      !$OMP MASTER
233      CALL VTe(VTHallo)
234      !$OMP END MASTER
235    CASE DEFAULT
236      CALL abort_gcm("vlspltgen_loc", "schema non parallelise", 1)
237    END SELECT
238
239  enddo !DO iq=1,nqtot
240
241
242  !$OMP BARRIER
243  !$OMP MASTER
244  CALL VTb(VTHallo)
245  !$OMP END MASTER
246
247  CALL SendRequest(MyRequest1)
248
249  !$OMP MASTER
250  CALL VTe(VTHallo)
251  !$OMP END MASTER
252  !$OMP BARRIER
253
254  ! verif temporaire
255  ijb = ij_begin - 2 * iip1
256  ije = ij_end + 2 * iip1
257  IF (pole_nord) ijb = ij_begin
258  IF (pole_sud)  ije = ij_end
259  CALL check_isotopes(zq, ij_begin, ij_end, 'vlspltgen_loc 280')
260
261  do iq = 1, nqtot
262    IF(tracers(iq)%parent /= 'air') CYCLE
263    !WRITE(*,*) 'vlspltgen 279: iq=',iq
264
265    SELECT CASE(tracers(iq)%iadv)
266    CASE(0); CYCLE
267    CASE(10)
268    CASE(14)
269    CASE DEFAULT
270      CALL abort_gcm("vlspltgen_p", "schema non parallelise", 1)
271    END SELECT
272
273  enddo
274  !$OMP BARRIER
275  !$OMP MASTER
276  CALL VTb(VTHallo)
277  !$OMP END MASTER
278
279  ! CALL WaitRecvRequest(MyRequest1)
280  ! CALL WaitSendRequest(MyRequest1)
281  !$OMP BARRIER
282  CALL WaitRequest(MyRequest1)
283
284
285  !$OMP MASTER
286  CALL VTe(VTHallo)
287  !$OMP END MASTER
288  !$OMP BARRIER
289
290  IF(isoCheck) THEN
291    CALL check_isotopes(zq, ij_begin, ij_end, 'vlspltgen_loc 326')
292    ijb = ij_begin - 2 * iip1
293    ije = ij_end + 2 * iip1
294    IF (pole_nord) ijb = ij_begin
295    IF (pole_sud)  ije = ij_end
296    CALL check_isotopes(zq, ijb, ije, 'vlspltgen_loc 336')
297  END IF
298
299  do iq = 1, nqtot
300    IF(tracers(iq)%parent /= 'air') CYCLE
301
302    SELECT CASE(tracers(iq)%iadv)
303    CASE(0); CYCLE
304    CASE(10); CALL   vly_loc(zq, pente_max, zm, mv, iq)
305    CASE(14); CALL vlyqs_loc(zq, pente_max, zm, mv, qsat, iq)
306    CASE DEFAULT
307      CALL abort_gcm("vlspltgen_p", "schema non parallelise", 1)
308    END SELECT
309
310  enddo
311
312  CALL check_isotopes(zq, ij_begin, ij_end, 'vlspltgen_loc 357')
313
314  do iq = 1, nqtot
315    IF(tracers(iq)%parent /= 'air') CYCLE
316    SELECT CASE(tracers(iq)%iadv)
317    CASE(0); CYCLE
318    CASE(10, 14)
319      !$OMP BARRIER
320      CALL vlz_loc(zq, pente_max, zm, mw, &
321              ij_begin, ij_end, iq)
322      !$OMP BARRIER
323
324      !$OMP MASTER
325      CALL VTb(VTHallo)
326      !$OMP END MASTER
327
328      CALL Register_Hallo_u(zq(:, :, iq), llm, 2, 2, 2, 2, MyRequest2)
329      CALL Register_Hallo_u(zm(:, :, iq), llm, 1, 1, 1, 1, MyRequest2)
330      ! CRisi
331      do ifils = 1, tracers(iq)%nqDescen
332        iq2 = tracers(iq)%iqDescen(ifils)
333        CALL Register_Hallo_u(zq(:, :, iq2), llm, 2, 2, 2, 2, MyRequest2)
334        CALL Register_Hallo_u(zm(:, :, iq2), llm, 1, 1, 1, 1, MyRequest2)
335      enddo
336      !$OMP MASTER
337      CALL VTe(VTHallo)
338      !$OMP END MASTER
339      !$OMP BARRIER
340    CASE DEFAULT
341
342      CALL abort_gcm("vlspltgen_p", "schema non parallelise", 1)
343    END SELECT
344
345  enddo
346  !$OMP BARRIER
347
348  !$OMP MASTER
349  CALL VTb(VTHallo)
350  !$OMP END MASTER
351
352  CALL SendRequest(MyRequest2)
353
354  !$OMP MASTER
355  CALL VTe(VTHallo)
356  !$OMP END MASTER
357
358  CALL check_isotopes(zq, ij_begin, ij_end, 'vlspltgen_loc 429')
359
360  !$OMP BARRIER
361  do iq = 1, nqtot
362    IF(tracers(iq)%parent /= 'air') CYCLE
363    !WRITE(*,*) 'vlspltgen 409: iq=',iq
364
365    SELECT CASE(tracers(iq)%iadv)
366    CASE(0); CYCLE
367    CASE(10, 14)
368      !$OMP BARRIER
369    CASE DEFAULT
370      CALL abort_gcm("vlspltgen_p", "schema non parallelise", 1)
371    END SELECT
372
373  enddo
374  !WRITE(*,*) 'vlspltgen_loc 476'
375
376  !$OMP BARRIER
377  !WRITE(*,*) 'vlspltgen_loc 477'
378  !$OMP MASTER
379  CALL VTb(VTHallo)
380  !$OMP END MASTER
381
382  ! CALL WaitRecvRequest(MyRequest2)
383  ! CALL WaitSendRequest(MyRequest2)
384  !$OMP BARRIER
385  CALL WaitRequest(MyRequest2)
386
387  !$OMP MASTER
388  CALL VTe(VTHallo)
389  !$OMP END MASTER
390  !$OMP BARRIER
391
392
393  !WRITE(*,*) 'vlspltgen_loc 494'
394  CALL check_isotopes(zq, ij_begin, ij_end, 'vlspltgen_loc 461')
395
396  do iq = 1, nqtot
397    IF(tracers(iq)%parent /= 'air') CYCLE
398    SELECT CASE(tracers(iq)%iadv)
399    CASE(0); CYCLE
400    CASE(10); CALL   vly_loc(zq, pente_max, zm, mv, iq)
401    CASE(14); CALL vlyqs_loc(zq, pente_max, zm, mv, qsat, iq)
402    CASE DEFAULT
403      CALL abort_gcm("vlspltgen_p", "schema non parallelise", 1)
404    END SELECT
405
406  enddo !do iq=1,nqtot
407
408  CALL check_isotopes(zq, ij_begin, ij_end, 'vlspltgen_loc 493')
409
410  do iq = 1, nqtot
411    IF(tracers(iq)%parent /= 'air') CYCLE
412    SELECT CASE(tracers(iq)%iadv)
413    CASE(0); CYCLE
414    CASE(10); CALL   vlx_loc(zq, pente_max, zm, mu, &
415            ij_begin, ij_end, iq)
416    CASE(14); CALL vlxqs_loc(zq, pente_max, zm, mu, &
417            qsat, ij_begin, ij_end, iq)
418    CASE DEFAULT
419      CALL abort_gcm("vlspltgen_p", "schema non parallelise", 1)
420    END SELECT
421
422  enddo !do iq=1,nqtot
423
424  !WRITE(*,*) 'vlspltgen 550: apres derniere serie de CALL vlx'
425  CALL check_isotopes(zq, ij_begin, ij_end, 'vlspltgen_loc 521')
426
427  ijb = ij_begin
428  ije = ij_end
429  !$OMP BARRIER
430
431  DO iq = 1, nqtot
432    !$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
433    DO l = 1, llm
434      DO ij = ijb, ije
435        ! print *,'zq-->',ij,l,iq,zq(ij,l,iq)
436        ! print *,'q-->',ij,l,iq,q(ij,l,iq)
437        q(ij, l, iq) = zq(ij, l, iq)
438      ENDDO
439    ENDDO
440    !$OMP END DO NOWAIT
441
442    !$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
443    DO l = 1, llm
444      DO ij = ijb, ije - iip1 + 1, iip1
445        q(ij + iim, l, iq) = q(ij, l, iq)
446      ENDDO
447    ENDDO
448    !$OMP END DO NOWAIT
449  ENDDO !DO iq=1,nqtot
450
451  CALL check_isotopes(q, ij_begin, ij_end, 'vlspltgen_loc 557')
452
453  !$OMP BARRIER
454
455  !c$OMP MASTER
456  ! CALL WaitSendRequest(MyRequest1)
457  ! CALL WaitSendRequest(MyRequest2)
458  !c$OMP END MASTER
459  !c$OMP BARRIER
460
461  !WRITE(*,*) 'vlspltgen 597: sortie'
462
463END SUBROUTINE vlspltgen_loc
Note: See TracBrowser for help on using the repository browser.