source: LMDZ5/branches/testing/libf/dyn3dpar/advect_new_p.F @ 5312

Last change on this file since 5312 was 2641, checked in by Laurent Fairhead, 8 years ago

Merged trunk changes r2593:2640 into testing branch
  • 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:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 7.1 KB
RevLine 
[774]1!
2! $Header$
3!
[763]4      SUBROUTINE advect_new_p(ucov,vcov,teta,w,massebx,masseby,
5     &                        du,dv,dteta)
[1864]6      USE parallel_lmdz
[763]7      USE write_field_p
[2641]8      USE comconst_mod, ONLY: daysec
9      USE logic_mod, ONLY: conser
10     
[763]11      IMPLICIT NONE
12c=======================================================================
13c
14c   Auteurs:  P. Le Van , Fr. Hourdin  .
15c   -------
16c
17c   Objet:
18c   ------
19c
20c   *************************************************************
21c   .... calcul des termes d'advection vertic.pour u,v,teta,q ...
22c   *************************************************************
23c        ces termes sont ajoutes a du,dv,dteta et dq .
24c  Modif F.Forget 03/94 : on retire q de advect
25c
26c=======================================================================
27c-----------------------------------------------------------------------
28c   Declarations:
29c   -------------
30
31#include "dimensions.h"
32#include "paramet.h"
33#include "comgeom.h"
34
35c   Arguments:
36c   ----------
37
38      REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm),teta(ip1jmp1,llm)
39      REAL massebx(ip1jmp1,llm),masseby(ip1jm,llm),w(ip1jmp1,llm)
40      REAL dv(ip1jm,llm),du(ip1jmp1,llm),dteta(ip1jmp1,llm)
41      REAL,SAVE :: dv1(ip1jm,llm),du1(ip1jmp1,llm),dteta1(ip1jmp1,llm)
42      REAL,SAVE :: dv2(ip1jm,llm),du2(ip1jmp1,llm),dteta2(ip1jmp1,llm)
43c   Local:
44c   ------
45
46      REAL,SAVE :: uav(ip1jmp1,llm),vav(ip1jm,llm)
47      REAL wsur2(ip1jmp1)
48      REAL unsaire2(ip1jmp1), ge(ip1jmp1)
49      REAL deuxjour, ww, gt, uu, vv
50
51      INTEGER  ij,l,ijb,ije
52
53      EXTERNAL  SSUM
54      REAL      SSUM
55
56c-----------------------------------------------------------------------
57c   2. Calculs preliminaires:
58c   -------------------------
59
60      IF (conser)  THEN
61         deuxjour = 2. * daysec
62
63         DO   1  ij   = 1, ip1jmp1
64         unsaire2(ij) = unsaire(ij) * unsaire(ij)
65   1     CONTINUE
66      END IF
67
68
69c------------------  -yy ----------------------------------------------
70c   .  Calcul de     u
71
72c$OMP MASTER
73      ijb=ij_begin
74      ije=ij_end
75      if (pole_nord) ijb=ijb+iip1
76      if (pole_sud)  ije=ije-iip1
77
78      DO ij=ijb,ije
79        du2(ij,1)=0.
80      ENDDO
81     
82      ijb=ij_begin
83      ije=ij_end
84      if (pole_sud)  ije=ij_end-iip1
85     
86      DO ij=ijb,ije
87        dv2(ij,1)=0.
88      ENDDO
89     
90      ijb=ij_begin
91      ije=ij_end
92
93      DO ij=ijb,ije
94        dteta2(ij,1)=0.
95      ENDDO
96c$OMP END MASTER
97
98 
99c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)   
100      DO  l=1,llm
101         
102         ijb=ij_begin
103         ije=ij_end
104         if (pole_nord) ijb=ijb+iip1
105         if (pole_sud)  ije=ije-iip1
106         
107c         DO    ij     = iip2, ip1jmp1
108c            uav(ij,l) = 0.25 * ( ucov(ij,l) + ucov(ij-iip1,l) )
109c         ENDDO
110
111c         DO    ij     = iip2, ip1jm
112c            uav(ij,l) = uav(ij,l) + uav(ij+iip1,l)
113c         ENDDO
114         
115         DO    ij     = ijb, ije
116                 
117           uav(ij,l)=0.25*(ucov(ij,l)+ucov(ij-iip1,l))
[2641]118     .                     +0.25*(ucov(ij+iip1,l)+ucov(ij,l))
[763]119         ENDDO
120         
121         if (pole_nord) then
122           DO      ij         = 1, iip1
123              uav(ij      ,l) = 0.
124           ENDDO
125         endif
126         
127         if (pole_sud) then
128           DO      ij         = 1, iip1
129              uav(ip1jm+ij,l) = 0.
130           ENDDO
131         endif
132         
133      ENDDO
134c$OMP END DO     
135c      call write_field3d_p('uav',reshape(uav,(/iip1,jjp1,llm/)))
136     
137c------------------  -xx ----------------------------------------------
138c   .  Calcul de     v
139     
140      ijb=ij_begin
141      ije=ij_end
142      if (pole_sud)  ije=ij_end-iip1
143
144c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
145      DO  l=1,llm
146         
147         DO    ij   = ijb+1, ije
148           vav(ij,l) = 0.25 * ( vcov(ij,l) + vcov(ij-1,l) )
149         ENDDO
150         
151         DO    ij   = ijb,ije,iip1
152          vav(ij,l) = vav(ij+iim,l)
153         ENDDO
154         
155         
156         DO    ij   = ijb, ije-1
157          vav(ij,l) = vav(ij,l) + vav(ij+1,l)
158         ENDDO
159         
160         DO    ij       = ijb, ije, iip1
161          vav(ij+iim,l) = vav(ij,l)
162         ENDDO
163         
164      ENDDO
165c$OMP END DO
166c       call write_field3d_p('vav',reshape(vav,(/iip1,jjm,llm/)))
167
168c-----------------------------------------------------------------------
169c$OMP BARRIER
170
171c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
172      DO 20 l = 1, llmm1
173
174
175c       ......   calcul de  - w/2.    au niveau  l+1   .......
176      ijb=ij_begin
177      ije=ij_end+iip1
178      if (pole_sud)  ije=ij_end
179     
180      DO 5   ij   = ijb, ije
181      wsur2( ij ) = - 0.5 * w( ij,l+1 )
182   5  CONTINUE
183
184
185c     .....................     calcul pour  du     ..................
186     
187      ijb=ij_begin
188      ije=ij_end
189      if (pole_nord) ijb=ijb+iip1
190      if (pole_sud)  ije=ije-iip1
191         
192      DO 6 ij = ijb ,ije-1
193      ww        = wsur2 (  ij  )     + wsur2( ij+1 )
194      uu        = 0.5 * ( ucov(ij,l) + ucov(ij,l+1) )
195      du1(ij,l)  =  ww * ( uu - uav(ij, l ) )/massebx(ij, l )
196      du2(ij,l+1)=  ww * ( uu - uav(ij,l+1) )/massebx(ij,l+1)
197   6  CONTINUE
198
199c     .................    calcul pour   dv      .....................
200      ijb=ij_begin
201      ije=ij_end
202      if (pole_sud)  ije=ij_end-iip1
203     
204      DO 8 ij = ijb, ije
205      ww        = wsur2( ij+iip1 )   + wsur2( ij )
206      vv        = 0.5 * ( vcov(ij,l) + vcov(ij,l+1) )
207      dv1(ij,l)  =  ww * (vv - vav(ij, l ) )/masseby(ij, l )
208      dv2(ij,l+1)=  ww * (vv - vav(ij,l+1) )/masseby(ij,l+1)
209   8  CONTINUE
210
211c
212
213c     ............................................................
214c     ...............    calcul pour   dh      ...................
215c     ............................................................
216
217c                       ---z
218c       calcul de  - d( teta  * w )      qu'on ajoute a   dh
219c                   ...............
220        ijb=ij_begin
221        ije=ij_end
222       
223        DO 15 ij = ijb, ije
224         ww            = wsur2(ij) * (teta(ij,l) + teta(ij,l+1) )
225         dteta1(ij, l ) =   ww
226         dteta2(ij,l+1) =   ww
227  15    CONTINUE
228
229c ym ---> conser a voir plus tard
230
231c      IF( conser)  THEN
232c       
233c        DO 17 ij = 1,ip1jmp1
234c        ge(ij)   = wsur2(ij) * wsur2(ij) * unsaire2(ij)
235c  17    CONTINUE
236c        gt       = SSUM( ip1jmp1,ge,1 )
237c        gtot(l)  = deuxjour * SQRT( gt/ip1jmp1 )
238c      END IF
239
240  20  CONTINUE
241c$OMP END DO
242
243      ijb=ij_begin
244      ije=ij_end
245      if (pole_nord) ijb=ijb+iip1
246      if (pole_sud)  ije=ije-iip1
247     
248c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)       
249      DO l=1,llm
250        DO ij=ijb,ije-1
[2641]251          du(ij,l)=du(ij,l)+du2(ij,l)-du1(ij,l)
252        ENDDO
[763]253
254        DO   ij   = ijb+iip1-1, ije, iip1
255         du( ij, l  ) = du( ij -iim, l  )
256        ENDDO
257      ENDDO
258c$OMP END DO NOWAIT
259     
260      ijb=ij_begin
261      ije=ij_end
262      if (pole_sud)  ije=ij_end-iip1
263
264c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
265      DO l=1,llm
266        DO ij=ijb,ije
[2641]267          dv(ij,l)=dv(ij,l)+dv2(ij,l)-dv1(ij,l)
268        ENDDO
[763]269      ENDDO
270c$OMP END DO NOWAIT     
271      ijb=ij_begin
272      ije=ij_end
273
274c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)   
275      DO l=1,llm
276        DO ij=ijb,ije
[2641]277          dteta(ij,l)=dteta(ij,l)+dteta2(ij,l)-dteta1(ij,l)
278        ENDDO
[763]279      ENDDO
280c$OMP END DO NOWAIT     
281
282      RETURN
283      END
Note: See TracBrowser for help on using the repository browser.