source: LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/vlspltgen_p.F @ 1114

Last change on this file since 1114 was 1114, checked in by jghattas, 15 years ago

Creation du module infotrac:

  • contient les variables de advtrac.h
  • contient la subroutine iniadvtrac renommer en infotrac_init
  • le nombre des traceurs est lu dans tracer.def en dynamique (ou par default ou recu par INCA)
  • ce module est utilise dans la dynamique et la physique
  • contient aussi la variable nbtr qui avant etait stockee dans dimphy

Le fichier advtrac.h n'existe plus.
La compilation ne prend plus en compte le nombre de traceur.

/JG
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 10.7 KB
Line 
1!
2! $Header$
3!
4       SUBROUTINE vlspltgen_p( q,iadv,pente_max,masse,w,pbaru,pbarv,pdt,
5     ,                                  p,pk,teta                 )
6c
7c     Auteurs:   P.Le Van, F.Hourdin, F.Forget, F.Codron
8c
9c    ********************************************************************
10c          Shema  d'advection " pseudo amont " .
11c      + test sur humidite specifique: Q advecte< Qsat aval
12c                   (F. Codron, 10/99)
13c    ********************************************************************
14c     q,pbaru,pbarv,w sont des arguments d'entree  pour le s-pg ....
15c
16c     pente_max facteur de limitation des pentes: 2 en general
17c                                                0 pour un schema amont
18c     pbaru,pbarv,w flux de masse en u ,v ,w
19c     pdt pas de temps
20c
21c     teta temperature potentielle, p pression aux interfaces,
22c     pk exner au milieu des couches necessaire pour calculer Qsat
23c   --------------------------------------------------------------------
24      USE parallel
25      USE mod_hallo
26      USE Write_Field_p
27      USE VAMPIR
28      USE infotrac, ONLY : nqtot
29      IMPLICIT NONE
30
31c
32#include "dimensions.h"
33#include "paramet.h"
34#include "logic.h"
35#include "comvert.h"
36#include "comconst.h"
37
38c
39c   Arguments:
40c   ----------
41      INTEGER iadv(nqtot)
42      REAL masse(ip1jmp1,llm),pente_max
43      REAL pbaru( ip1jmp1,llm ),pbarv( ip1jm,llm)
44      REAL q(ip1jmp1,llm,nqtot)
45      REAL w(ip1jmp1,llm),pdt
46      REAL p(ip1jmp1,llmp1),teta(ip1jmp1,llm),pk(ip1jmp1,llm)
47c
48c      Local
49c   ---------
50c
51      INTEGER ij,l
52c
53      REAL,SAVE :: qsat(ip1jmp1,llm)
54      REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: zm
55      REAL,SAVE :: mu(ip1jmp1,llm)
56      REAL,SAVE :: mv(ip1jm,llm)
57      REAL,SAVE :: mw(ip1jmp1,llm+1)
58      REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: zq
59      REAL zzpbar, zzw
60
61      REAL qmin,qmax
62      DATA qmin,qmax/0.,1.e33/
63
64c--pour rapport de melange saturant--
65
66      REAL rtt,retv,r2es,r3les,r3ies,r4les,r4ies,play
67      REAL ptarg,pdelarg,foeew,zdelta
68      REAL tempe(ip1jmp1)
69      INTEGER ijb,ije,iq
70      LOGICAL, SAVE :: firstcall=.TRUE.
71
72      type(request) :: MyRequest1
73      type(request) :: MyRequest2
74
75c    fonction psat(T)
76
77       FOEEW ( PTARG,PDELARG ) = EXP (
78     *          (R3LES*(1.-PDELARG)+R3IES*PDELARG) * (PTARG-RTT)
79     * / (PTARG-(R4LES*(1.-PDELARG)+R4IES*PDELARG)) )
80
81        r2es  = 380.11733
82        r3les = 17.269
83        r3ies = 21.875
84        r4les = 35.86
85        r4ies = 7.66
86        retv = 0.6077667
87        rtt  = 273.16
88
89c Allocate variables depending on dynamic variable nqtot
90         IF (firstcall) THEN
91            firstcall=.FALSE.
92            ALLOCATE(zm(ip1jmp1,llm,nqtot))
93            ALLOCATE(zq(ip1jmp1,llm,nqtot))
94         END IF
95
96c-- Calcul de Qsat en chaque point
97c-- approximation: au milieu des couches play(l)=(p(l)+p(l+1))/2
98c   pour eviter une exponentielle.
99
100      call SetTag(MyRequest1,100)
101      call SetTag(MyRequest2,101)
102
103       
104        ijb=ij_begin-iip1
105        ije=ij_end+iip1
106        if (pole_nord) ijb=ij_begin
107        if (pole_sud) ije=ij_end
108       
109c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
110        DO l = 1, llm
111         DO ij = ijb, ije
112          tempe(ij) = teta(ij,l) * pk(ij,l) /cpp
113         ENDDO
114         DO ij = ijb, ije
115          zdelta = MAX( 0., SIGN(1., rtt - tempe(ij)) )
116          play   = 0.5*(p(ij,l)+p(ij,l+1))
117          qsat(ij,l) = MIN(0.5, r2es* FOEEW(tempe(ij),zdelta) / play )
118          qsat(ij,l) = qsat(ij,l) / ( 1. - retv * qsat(ij,l) )
119         ENDDO
120        ENDDO
121c$OMP END DO NOWAIT
122c      PRINT*,'Debut vlsplt version debug sans vlyqs'
123
124        zzpbar = 0.5 * pdt
125        zzw    = pdt
126
127      ijb=ij_begin
128      ije=ij_end
129      if (pole_nord) ijb=ijb+iip1
130      if (pole_sud)  ije=ije-iip1
131
132c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)       
133      DO l=1,llm
134        DO ij = ijb,ije
135            mu(ij,l)=pbaru(ij,l) * zzpbar
136         ENDDO
137      ENDDO
138c$OMP END DO NOWAIT
139     
140      ijb=ij_begin-iip1
141      ije=ij_end
142      if (pole_nord) ijb=ij_begin
143      if (pole_sud)  ije=ij_end-iip1
144
145c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
146      DO l=1,llm
147         DO ij=ijb,ije
148            mv(ij,l)=pbarv(ij,l) * zzpbar
149         ENDDO
150      ENDDO
151c$OMP END DO NOWAIT
152
153      ijb=ij_begin
154      ije=ij_end
155
156c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
157      DO l=1,llm
158         DO ij=ijb,ije
159            mw(ij,l)=w(ij,l) * zzw
160         ENDDO
161      ENDDO
162c$OMP END DO NOWAIT
163
164c$OMP MASTER
165      DO ij=ijb,ije
166         mw(ij,llm+1)=0.
167      ENDDO
168c$OMP END MASTER
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
186c$OMP BARRIER           
187      DO iq=1,nqtot
188
189        if(iadv(iq) == 0) then
190       
191          cycle
192       
193        else if (iadv(iq)==10) then
194
195#ifdef _ADV_HALO       
196          call vlx_p(zq(1,1,iq),pente_max,zm(1,1,iq),mu,
197     &               ij_begin,ij_begin+2*iip1-1)
198          call vlx_p(zq(1,1,iq),pente_max,zm(1,1,iq),mu,
199     &               ij_end-2*iip1+1,ij_end)
200#else
201          call vlx_p(zq(1,1,iq),pente_max,zm(1,1,iq),mu,
202     &               ij_begin,ij_end)
203#endif
204
205c$OMP MASTER
206          call VTb(VTHallo)
207c$OMP END MASTER
208          call Register_Hallo(zq(1,1,iq),ip1jmp1,llm,2,2,2,2,MyRequest1)
209          call Register_Hallo(zm(1,1,iq),ip1jmp1,llm,1,1,1,1,MyRequest1)
210
211c$OMP MASTER
212          call VTe(VTHallo)
213c$OMP END MASTER
214        else if (iadv(iq)==14) then
215
216#ifdef _ADV_HALO           
217          call vlxqs_p(zq(1,1,iq),pente_max,zm(1,1,iq),mu,qsat,
218     &                 ij_begin,ij_begin+2*iip1-1)
219          call vlxqs_p(zq(1,1,iq),pente_max,zm(1,1,iq),mu,qsat,
220     &                 ij_end-2*iip1+1,ij_end)
221#else
222
223          call vlxqs_p(zq(1,1,iq),pente_max,zm(1,1,iq),mu,qsat,
224     &                 ij_begin,ij_end)
225#endif
226
227c$OMP MASTER
228          call VTb(VTHallo)
229c$OMP END MASTER
230
231          call Register_Hallo(zq(1,1,iq),ip1jmp1,llm,2,2,2,2,MyRequest1)
232          call Register_Hallo(zm(1,1,iq),ip1jmp1,llm,1,1,1,1,MyRequest1)
233
234c$OMP MASTER
235          call VTe(VTHallo)
236c$OMP END MASTER
237        else
238       
239          stop 'vlspltgen_p : schema non parallelise'
240     
241        endif
242     
243      enddo
244     
245     
246c$OMP BARRIER     
247c$OMP MASTER     
248      call VTb(VTHallo)
249c$OMP END MASTER
250
251      call SendRequest(MyRequest1)
252
253c$OMP MASTER
254      call VTe(VTHallo)
255c$OMP END MASTER       
256c$OMP BARRIER
257      do iq=1,nqtot
258
259        if(iadv(iq) == 0) then
260       
261          cycle
262       
263        else if (iadv(iq)==10) then
264
265#ifdef _ADV_HALLO
266          call vlx_p(zq(1,1,iq),pente_max,zm(1,1,iq),mu,
267     &               ij_begin+2*iip1,ij_end-2*iip1)
268#endif       
269        else if (iadv(iq)==14) then
270#ifdef _ADV_HALLO
271          call vlxqs_p(zq(1,1,iq),pente_max,zm(1,1,iq),mu,qsat,
272     &                 ij_begin+2*iip1,ij_end-2*iip1)
273#endif   
274        else
275       
276          stop 'vlspltgen_p : schema non parallelise'
277     
278        endif
279     
280      enddo
281c$OMP BARRIER     
282c$OMP MASTER
283      call VTb(VTHallo)
284c$OMP END MASTER
285
286!      call WaitRecvRequest(MyRequest1)
287!      call WaitSendRequest(MyRequest1)
288c$OMP BARRIER
289       call WaitRequest(MyRequest1)
290
291
292c$OMP MASTER
293      call VTe(VTHallo)
294c$OMP END MASTER
295c$OMP BARRIER
296 
297      do iq=1,nqtot
298
299        if(iadv(iq) == 0) then
300       
301          cycle
302       
303        else if (iadv(iq)==10) then
304       
305          call vly_p(zq(1,1,iq),pente_max,zm(1,1,iq),mv)
306 
307        else if (iadv(iq)==14) then
308     
309          call vlyqs_p(zq(1,1,iq),pente_max,zm(1,1,iq),mv,qsat)
310 
311        else
312       
313          stop 'vlspltgen_p : schema non parallelise'
314     
315        endif
316       
317       enddo
318
319
320      do iq=1,nqtot
321
322        if(iadv(iq) == 0) then
323         
324          cycle
325       
326        else if (iadv(iq)==10 .or. iadv(iq)==14 ) then
327
328c$OMP BARRIER       
329#ifdef _ADV_HALLO
330          call vlz_p(zq(1,1,iq),pente_max,zm(1,1,iq),mw,
331     &               ij_begin,ij_begin+2*iip1-1)
332          call vlz_p(zq(1,1,iq),pente_max,zm(1,1,iq),mw,
333     &               ij_end-2*iip1+1,ij_end)
334#else
335          call vlz_p(zq(1,1,iq),pente_max,zm(1,1,iq),mw,
336     &               ij_begin,ij_end)
337#endif
338c$OMP BARRIER
339
340c$OMP MASTER
341          call VTb(VTHallo)
342c$OMP END MASTER
343
344          call Register_Hallo(zq(1,1,iq),ip1jmp1,llm,2,2,2,2,MyRequest2)
345          call Register_Hallo(zm(1,1,iq),ip1jmp1,llm,1,1,1,1,MyRequest2)
346
347c$OMP MASTER
348          call VTe(VTHallo)
349c$OMP END MASTER       
350c$OMP BARRIER
351        else
352       
353          stop 'vlspltgen_p : schema non parallelise'
354     
355        endif
356     
357      enddo
358c$OMP BARRIER     
359
360c$OMP MASTER       
361      call VTb(VTHallo)
362c$OMP END MASTER
363
364      call SendRequest(MyRequest2)
365
366c$OMP MASTER
367      call VTe(VTHallo)
368c$OMP END MASTER       
369
370c$OMP BARRIER
371      do iq=1,nqtot
372
373        if(iadv(iq) == 0) then
374         
375          cycle
376       
377        else if (iadv(iq)==10 .or. iadv(iq)==14 ) then
378c$OMP BARRIER       
379
380#ifdef _ADV_HALLO
381          call vlz_p(zq(1,1,iq),pente_max,zm(1,1,iq),mw,
382     &               ij_begin+2*iip1,ij_end-2*iip1)
383#endif
384
385c$OMP BARRIER       
386        else
387       
388          stop 'vlspltgen_p : schema non parallelise'
389     
390        endif
391     
392      enddo
393
394c$OMP BARRIER
395c$OMP MASTER
396      call VTb(VTHallo)
397c$OMP END MASTER
398
399!      call WaitRecvRequest(MyRequest2)
400!      call WaitSendRequest(MyRequest2)
401c$OMP BARRIER
402       CALL WaitRequest(MyRequest2)
403
404c$OMP MASTER
405      call VTe(VTHallo)
406c$OMP END MASTER
407c$OMP BARRIER
408
409
410      do iq=1,nqtot
411
412        if(iadv(iq) == 0) then
413       
414          cycle
415       
416        else if (iadv(iq)==10) then
417       
418          call vly_p(zq(1,1,iq),pente_max,zm(1,1,iq),mv)
419 
420        else if (iadv(iq)==14) then
421     
422          call vlyqs_p(zq(1,1,iq),pente_max,zm(1,1,iq),mv,qsat)
423 
424        else
425       
426          stop 'vlspltgen_p : schema non parallelise'
427     
428        endif
429       
430       enddo
431
432      do iq=1,nqtot
433
434        if(iadv(iq) == 0) then
435         
436          cycle
437       
438        else if (iadv(iq)==10) then
439       
440          call vlx_p(zq(1,1,iq),pente_max,zm(1,1,iq),mu,
441     &               ij_begin,ij_end)
442 
443        else if (iadv(iq)==14) then
444     
445          call vlxqs_p(zq(1,1,iq),pente_max,zm(1,1,iq),mu,qsat,
446     &                 ij_begin,ij_end)
447 
448        else
449       
450          stop 'vlspltgen_p : schema non parallelise'
451     
452        endif
453       
454       enddo
455
456     
457      ijb=ij_begin
458      ije=ij_end
459c$OMP BARRIER     
460
461
462      DO iq=1,nqtot
463
464c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)         
465        DO l=1,llm
466           DO ij=ijb,ije
467c             print *,'zq-->',ij,l,iq,zq(ij,l,iq)
468c            print *,'q-->',ij,l,iq,q(ij,l,iq)
469             q(ij,l,iq)=zq(ij,l,iq)
470           ENDDO
471        ENDDO
472c$OMP END DO NOWAIT         
473
474c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
475        DO l=1,llm
476           DO ij=ijb,ije-iip1+1,iip1
477              q(ij+iim,l,iq)=q(ij,l,iq)
478           ENDDO
479        ENDDO
480c$OMP END DO NOWAIT 
481
482      ENDDO
483       
484     
485c$OMP BARRIER
486
487cc$OMP MASTER     
488c      call WaitSendRequest(MyRequest1)
489c      call WaitSendRequest(MyRequest2)
490cc$OMP END MASTER
491cc$OMP BARRIER
492
493      RETURN
494      END
Note: See TracBrowser for help on using the repository browser.