source: LMDZ6/branches/SETHET_DECOUPLE/libf/dyn3dmem/vlspltgen_loc.F @ 5437

Last change on this file since 5437 was 2603, checked in by Ehouarn Millour, 8 years ago

Cleanup in the dynamics: turn logic.h into module logic_mod.F90
EM
  • 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: 15.2 KB
Line 
1!
2! $Header$
3!
4       SUBROUTINE vlspltgen_loc( q,iadv,pente_max,masse,w,pbaru,pbarv,
5     &                           pdt, p,pk,teta                 )
6     
7c
8c     Auteurs:   P.Le Van, F.Hourdin, F.Forget, F.Codron
9c
10c    ********************************************************************
11c          Shema  d'advection " pseudo amont " .
12c      + test sur humidite specifique: Q advecte< Qsat aval
13c                   (F. Codron, 10/99)
14c    ********************************************************************
15c     q,pbaru,pbarv,w sont des arguments d'entree  pour le s-pg ....
16c
17c     pente_max facteur de limitation des pentes: 2 en general
18c                                                0 pour un schema amont
19c     pbaru,pbarv,w flux de masse en u ,v ,w
20c     pdt pas de temps
21c
22c     teta temperature potentielle, p pression aux interfaces,
23c     pk exner au milieu des couches necessaire pour calculer Qsat
24c   --------------------------------------------------------------------
25      USE parallel_lmdz
26      USE mod_hallo
27      USE Write_Field_loc
28      USE VAMPIR
29      ! CRisi: on rajoute variables utiles d'infotrac 
30      USE infotrac, ONLY : nqtot,nqperes,nqdesc,nqfils,iqfils,
31     &    ok_iso_verif
32      USE vlspltgen_mod
33      USE comconst_mod, ONLY: cpp
34      IMPLICIT NONE
35
36c
37      include "dimensions.h"
38      include "paramet.h"
39
40c
41c   Arguments:
42c   ----------
43      INTEGER iadv(nqtot)
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)
50c
51c      Local
52c   ---------
53c
54      INTEGER ij,l
55c
56      REAL zzpbar, zzw
57
58      REAL qmin,qmax
59      DATA qmin,qmax/0.,1.e33/
60
61c--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)
73c    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
87c Allocate variables depending on dynamic variable nqtot
88
89         IF (firstcall) THEN
90            firstcall=.FALSE.
91         END IF
92c-- Calcul de Qsat en chaque point
93c-- approximation: au milieu des couches play(l)=(p(l)+p(l+1))/2
94c   pour eviter une exponentielle.
95
96      call SetTag(MyRequest1,100)
97      call SetTag(MyRequest2,101)
98
99       
100        ijb=ij_begin-iip1
101        ije=ij_end+iip1
102        if (pole_nord) ijb=ij_begin
103        if (pole_sud) ije=ij_end
104       
105c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
106        DO l = 1, llm
107         DO ij = ijb, ije
108          tempe(ij) = teta(ij,l) * pk(ij,l) /cpp
109         ENDDO
110         DO ij = ijb, ije
111          zdelta = MAX( 0., SIGN(1., rtt - tempe(ij)) )
112          play   = 0.5*(p(ij,l)+p(ij,l+1))
113          qsat(ij,l) = MIN(0.5, r2es* FOEEW(tempe(ij),zdelta) / play )
114          qsat(ij,l) = qsat(ij,l) / ( 1. - retv * qsat(ij,l) )
115         ENDDO
116        ENDDO
117c$OMP END DO NOWAIT
118c      PRINT*,'Debut vlsplt version debug sans vlyqs'
119
120        zzpbar = 0.5 * pdt
121        zzw    = pdt
122
123      ijb=ij_begin
124      ije=ij_end
125      if (pole_nord) ijb=ijb+iip1
126      if (pole_sud)  ije=ije-iip1
127
128c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)       
129      DO l=1,llm
130        DO ij = ijb,ije
131            mu(ij,l)=pbaru(ij,l) * zzpbar
132         ENDDO
133      ENDDO
134c$OMP END DO NOWAIT
135     
136      ijb=ij_begin-iip1
137      ije=ij_end
138      if (pole_nord) ijb=ij_begin
139      if (pole_sud)  ije=ij_end-iip1
140
141c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
142      DO l=1,llm
143         DO ij=ijb,ije
144            mv(ij,l)=pbarv(ij,l) * zzpbar
145         ENDDO
146      ENDDO
147c$OMP END DO NOWAIT
148
149      ijb=ij_begin
150      ije=ij_end
151
152      DO iq=1,nqtot
153c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
154      DO l=1,llm
155         DO ij=ijb,ije
156            mw(ij,l,iq)=w(ij,l) * zzw
157         ENDDO
158      ENDDO
159c$OMP END DO NOWAIT
160      ENDDO
161
162      DO iq=1,nqtot 
163c$OMP MASTER
164      DO ij=ijb,ije
165         mw(ij,llm+1,iq)=0.
166      ENDDO
167c$OMP END MASTER
168      ENDDO
169
170c      CALL SCOPY(ijp1llm,q,1,zq,1)
171c      CALL SCOPY(ijp1llm,masse,1,zm,1)
172
173       ijb=ij_begin
174       ije=ij_end
175
176      DO iq=1,nqtot       
177c$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
182c$OMP END DO NOWAIT
183      ENDDO
184
185#ifdef DEBUG_IO     
186       CALL WriteField_u('mu',mu)
187       CALL WriteField_v('mv',mv)
188       CALL WriteField_u('mw',mw)
189       CALL WriteField_u('qsat',qsat)
190#endif
191
192      ! verif temporaire
193      ijb=ij_begin
194      ije=ij_end 
195      if (ok_iso_verif) then
196        call check_isotopes(zq,ijb,ije,'vlspltgen_loc 191')
197      endif !if (ok_iso_verif) then   
198
199c$OMP BARRIER           
200!      DO iq=1,nqtot
201      DO iq=1,nqperes ! CRisi: on ne boucle que sur les pères= ceux qui sont transportés directement par l'air
202       !write(*,*) 'vlspltgen 192: iq,iadv=',iq,iadv(iq)
203#ifdef DEBUG_IO   
204       CALL WriteField_u('zq',zq(:,:,iq))
205       CALL WriteField_u('zm',zm(:,:,iq))
206#endif
207        if(iadv(iq) == 0) then
208       
209          cycle
210       
211        else if (iadv(iq)==10) then
212
213#ifdef _ADV_HALO       
214! CRisi: on ajoute les nombres de fils et tableaux des fils
215! On suppose qu'on ne peut advecter les fils que par le schéma 10. 
216          call vlx_loc(zq,pente_max,zm,mu,
217     &                     ij_begin,ij_begin+2*iip1-1,iq)
218          call vlx_loc(zq,pente_max,zm,mu,
219     &               ij_end-2*iip1+1,ij_end,iq)
220#else
221          call vlx_loc(zq,pente_max,zm,mu,
222     &                     ij_begin,ij_end,iq)
223#endif
224
225c$OMP MASTER
226          call VTb(VTHallo)
227c$OMP END MASTER
228          call Register_Hallo_u(zq(:,:,iq),llm,2,2,2,2,MyRequest1)
229          call Register_Hallo_u(zm(:,:,iq),llm,1,1,1,1,MyRequest1)
230! CRisi
231          do ifils=1,nqdesc(iq)
232            iq2=iqfils(ifils,iq)
233            call Register_Hallo_u(zq(:,:,iq2),llm,2,2,2,2,MyRequest1)
234            call Register_Hallo_u(zm(:,:,iq2),llm,1,1,1,1,MyRequest1)
235          enddo
236
237c$OMP MASTER
238          call VTe(VTHallo)
239c$OMP END MASTER
240        else if (iadv(iq)==14) then
241
242#ifdef _ADV_HALO           
243          call vlxqs_loc(zq,pente_max,zm,mu,
244     &                   qsat,ij_begin,ij_begin+2*iip1-1,iq)
245          call vlxqs_loc(zq,pente_max,zm,mu,
246     &                   qsat,ij_end-2*iip1+1,ij_end,iq)
247#else
248          call vlxqs_loc(zq,pente_max,zm,mu,
249     &                   qsat,ij_begin,ij_end,iq)
250#endif
251
252c$OMP MASTER
253          call VTb(VTHallo)
254c$OMP END MASTER
255
256          call Register_Hallo_u(zq(:,:,iq),llm,2,2,2,2,MyRequest1)
257          call Register_Hallo_u(zm(:,:,iq),llm,1,1,1,1,MyRequest1)
258          do ifils=1,nqdesc(iq)
259            iq2=iqfils(ifils,iq)
260            call Register_Hallo_u(zq(:,:,iq2),llm,2,2,2,2,MyRequest1)
261            call Register_Hallo_u(zm(:,:,iq2),llm,1,1,1,1,MyRequest1)
262          enddo
263
264c$OMP MASTER
265          call VTe(VTHallo)
266c$OMP END MASTER
267        else
268       
269          stop 'vlspltgen_p : schema non parallelise'
270     
271        endif
272     
273      enddo !DO iq=1,nqperes
274     
275     
276c$OMP BARRIER     
277c$OMP MASTER     
278      call VTb(VTHallo)
279c$OMP END MASTER
280
281      call SendRequest(MyRequest1)
282
283c$OMP MASTER
284      call VTe(VTHallo)
285c$OMP END MASTER       
286c$OMP BARRIER
287
288      ! verif temporaire
289      ijb=ij_begin-2*iip1
290      ije=ij_end+2*iip1 
291      if (pole_nord) ijb=ij_begin
292      if (pole_sud)  ije=ij_end 
293      if (ok_iso_verif) then
294           call check_isotopes(zq,ij_begin,ij_end,'vlspltgen_loc 280')
295      endif !if (ok_iso_verif) then
296
297      do iq=1,nqperes
298        !write(*,*) 'vlspltgen 279: iq=',iq
299
300        if(iadv(iq) == 0) then
301       
302          cycle
303       
304        else if (iadv(iq)==10) then
305
306#ifdef _ADV_HALLO
307          call vlx_loc(zq,pente_max,zm,mu,
308     &                 ij_begin+2*iip1,ij_end-2*iip1,iq)
309#endif       
310        else if (iadv(iq)==14) then
311#ifdef _ADV_HALLO
312          call vlxqs_loc(zq,pente_max,zm,mu,
313     &                    qsat,ij_begin+2*iip1,ij_end-2*iip1,iq)
314#endif   
315        else
316       
317          stop 'vlspltgen_p : schema non parallelise'
318     
319        endif
320     
321      enddo
322c$OMP BARRIER     
323c$OMP MASTER
324      call VTb(VTHallo)
325c$OMP END MASTER
326
327!      call WaitRecvRequest(MyRequest1)
328!      call WaitSendRequest(MyRequest1)
329c$OMP BARRIER
330       call WaitRequest(MyRequest1)
331
332
333c$OMP MASTER
334      call VTe(VTHallo)
335c$OMP END MASTER
336c$OMP BARRIER
337
338     
339      if (ok_iso_verif) then
340           call check_isotopes(zq,ij_begin,ij_end,'vlspltgen_loc 326')
341      endif !if (ok_iso_verif) then       
342      if (ok_iso_verif) then
343           ijb=ij_begin-2*iip1
344           ije=ij_end+2*iip1
345           if (pole_nord) ijb=ij_begin
346           if (pole_sud)  ije=ij_end
347           call check_isotopes(zq,ijb,ije,'vlspltgen_loc 336')
348      endif !if (ok_iso_verif) then 
349
350      do iq=1,nqperes
351       !write(*,*) 'vlspltgen 321: iq=',iq
352#ifdef DEBUG_IO   
353       CALL WriteField_u('zq',zq(:,:,iq))
354       CALL WriteField_u('zm',zm(:,:,iq))
355#endif
356
357        if(iadv(iq) == 0) then
358       
359          cycle
360       
361        else if (iadv(iq)==10) then
362       
363          call vly_loc(zq,pente_max,zm,mv,iq)
364 
365        else if (iadv(iq)==14) then
366     
367          call vlyqs_loc(zq,pente_max,zm,mv,
368     &                   qsat,iq)
369 
370        else
371       
372          stop 'vlspltgen_p : schema non parallelise'
373     
374        endif
375       
376       enddo
377
378      if (ok_iso_verif) then
379           call check_isotopes(zq,ij_begin,ij_end,'vlspltgen_loc 357')
380      endif !if (ok_iso_verif) then
381
382      do iq=1,nqperes
383      !write(*,*) 'vlspltgen 349: iq=',iq
384#ifdef DEBUG_IO   
385       CALL WriteField_u('zq',zq(:,:,iq))
386       CALL WriteField_u('zm',zm(:,:,iq))
387#endif
388        if(iadv(iq) == 0) then
389         
390          cycle
391       
392        else if (iadv(iq)==10 .or. iadv(iq)==14 ) then
393
394c$OMP BARRIER       
395#ifdef _ADV_HALLO
396          call vlz_loc(zq,pente_max,zm,mw,
397     &               ij_begin,ij_begin+2*iip1-1,iq)
398          call vlz_loc(zq,pente_max,zm,mw,
399     &               ij_end-2*iip1+1,ij_end,iq)
400#else
401          call vlz_loc(zq,pente_max,zm,mw,
402     &               ij_begin,ij_end,iq)
403#endif
404c$OMP BARRIER
405
406c$OMP MASTER
407          call VTb(VTHallo)
408c$OMP END MASTER
409
410          call Register_Hallo_u(zq(:,:,iq),llm,2,2,2,2,MyRequest2)
411          call Register_Hallo_u(zm(:,:,iq),llm,1,1,1,1,MyRequest2)
412          ! CRisi
413          do ifils=1,nqdesc(iq)
414            iq2=iqfils(ifils,iq)
415            call Register_Hallo_u(zq(:,:,iq2),llm,2,2,2,2,MyRequest2)
416            call Register_Hallo_u(zm(:,:,iq2),llm,1,1,1,1,MyRequest2)
417          enddo     
418c$OMP MASTER
419          call VTe(VTHallo)
420c$OMP END MASTER       
421c$OMP BARRIER
422        else
423       
424          stop 'vlspltgen_p : schema non parallelise'
425     
426        endif
427     
428      enddo
429c$OMP BARRIER     
430
431c$OMP MASTER       
432      call VTb(VTHallo)
433c$OMP END MASTER
434
435      call SendRequest(MyRequest2)
436
437c$OMP MASTER
438      call VTe(VTHallo)
439c$OMP END MASTER       
440
441
442      if (ok_iso_verif) then
443           call check_isotopes(zq,ij_begin,ij_end,'vlspltgen_loc 429')
444      endif !if (ok_iso_verif) then
445
446c$OMP BARRIER
447      do iq=1,nqperes
448      !write(*,*) 'vlspltgen 409: iq=',iq
449
450        if(iadv(iq) == 0) then
451         
452          cycle
453       
454        else if (iadv(iq)==10 .or. iadv(iq)==14 ) then
455c$OMP BARRIER       
456
457#ifdef _ADV_HALLO
458          call vlz_loc(zq,pente_max,zm,mw,
459     &               ij_begin+2*iip1,ij_end-2*iip1,iq)
460#endif
461
462c$OMP BARRIER       
463        else
464       
465          stop 'vlspltgen_p : schema non parallelise'
466     
467        endif
468     
469      enddo
470      !write(*,*) 'vlspltgen_loc 476'
471
472c$OMP BARRIER
473      !write(*,*) 'vlspltgen_loc 477'
474c$OMP MASTER
475      call VTb(VTHallo)
476c$OMP END MASTER
477
478!      call WaitRecvRequest(MyRequest2)
479!      call WaitSendRequest(MyRequest2)
480c$OMP BARRIER
481       CALL WaitRequest(MyRequest2)
482
483c$OMP MASTER
484      call VTe(VTHallo)
485c$OMP END MASTER
486c$OMP BARRIER
487
488
489      !write(*,*) 'vlspltgen_loc 494'
490      if (ok_iso_verif) then
491           call check_isotopes(zq,ij_begin,ij_end,'vlspltgen_loc 461')
492      endif !if (ok_iso_verif) then
493
494      do iq=1,nqperes
495      !write(*,*) 'vlspltgen 449: iq=',iq
496#ifdef DEBUG_IO   
497       CALL WriteField_u('zq',zq(:,:,iq))
498       CALL WriteField_u('zm',zm(:,:,iq))
499#endif
500        if(iadv(iq) == 0) then
501       
502          cycle
503       
504        else if (iadv(iq)==10) then
505       
506          call vly_loc(zq,pente_max,zm,mv,iq)
507 
508        else if (iadv(iq)==14) then
509     
510          call vlyqs_loc(zq,pente_max,zm,mv,
511     &                   qsat,iq)
512 
513        else
514       
515          stop 'vlspltgen_p : schema non parallelise'
516     
517        endif
518       
519       enddo !do iq=1,nqperes
520
521      if (ok_iso_verif) then
522           call check_isotopes(zq,ij_begin,ij_end,'vlspltgen_loc 493')
523      endif !if (ok_iso_verif) then
524
525      do iq=1,nqperes
526      !write(*,*) 'vlspltgen 477: iq=',iq
527#ifdef DEBUG_IO   
528       CALL WriteField_u('zq',zq(:,:,iq))
529       CALL WriteField_u('zm',zm(:,:,iq))
530#endif
531        if(iadv(iq) == 0) then
532         
533          cycle
534       
535        else if (iadv(iq)==10) then
536       
537          call vlx_loc(zq,pente_max,zm,mu,
538     &               ij_begin,ij_end,iq)
539 
540        else if (iadv(iq)==14) then
541     
542          call vlxqs_loc(zq,pente_max,zm,mu,
543     &                 qsat, ij_begin,ij_end,iq)
544 
545        else
546       
547          stop 'vlspltgen_p : schema non parallelise'
548     
549        endif
550       
551       enddo !do iq=1,nqperes
552
553      !write(*,*) 'vlspltgen 550: apres derniere serie de call vlx'
554      if (ok_iso_verif) then
555           call check_isotopes(zq,ij_begin,ij_end,'vlspltgen_loc 521')
556      endif !if (ok_iso_verif) then
557     
558      ijb=ij_begin
559      ije=ij_end
560      !write(*,*) 'vlspltgen_loc 557'
561c$OMP BARRIER     
562
563      !write(*,*) 'vlspltgen_loc 559' 
564      DO iq=1,nqtot
565       !write(*,*) 'vlspltgen_loc 561, iq=',iq 
566#ifdef DEBUG_IO   
567       CALL WriteField_u('zq',zq(:,:,iq))
568       CALL WriteField_u('zm',zm(:,:,iq))
569#endif
570c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)         
571        DO l=1,llm
572           DO ij=ijb,ije
573c             print *,'zq-->',ij,l,iq,zq(ij,l,iq)
574c             print *,'q-->',ij,l,iq,q(ij,l,iq)
575             q(ij,l,iq)=zq(ij,l,iq)
576           ENDDO
577        ENDDO
578c$OMP END DO NOWAIT   
579      !write(*,*) 'vlspltgen_loc 575'     
580
581c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
582        DO l=1,llm
583           DO ij=ijb,ije-iip1+1,iip1
584              q(ij+iim,l,iq)=q(ij,l,iq)
585           ENDDO
586        ENDDO
587c$OMP END DO NOWAIT 
588      !write(*,*) 'vlspltgen_loc 583' 
589      ENDDO !DO iq=1,nqtot
590       
591      if (ok_iso_verif) then
592           call check_isotopes(q,ij_begin,ij_end,'vlspltgen_loc 557')
593      endif !if (ok_iso_verif) then
594
595c$OMP BARRIER
596
597cc$OMP MASTER     
598c      call WaitSendRequest(MyRequest1)
599c      call WaitSendRequest(MyRequest2)
600cc$OMP END MASTER
601cc$OMP BARRIER
602
603      !write(*,*) 'vlspltgen 597: sortie' 
604      RETURN
605      END
Note: See TracBrowser for help on using the repository browser.