source: LMDZ6/trunk/libf/dyn3dmem/advect_new_loc.f90 @ 5410

Last change on this file since 5410 was 5324, checked in by abarral, 5 weeks ago

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