My Project
 All Classes Files Functions Variables Macros
advect_p.F
Go to the documentation of this file.
1 !
2 ! $Header$
3 !
4  SUBROUTINE advect_p(ucov,vcov,teta,w,massebx,masseby,du,dv,dteta)
5  USE parallel
6  USE write_field_p
7  IMPLICIT NONE
8 c=======================================================================
9 c
10 c Auteurs: P. Le Van , Fr. Hourdin .
11 c -------
12 c
13 c Objet:
14 c ------
15 c
16 c *************************************************************
17 c .... calcul des termes d'advection vertic.pour u,v,teta,q ...
18 c *************************************************************
19 c ces termes sont ajoutes a du,dv,dteta et dq .
20 c Modif F.Forget 03/94 : on retire q de advect
21 c
22 c=======================================================================
23 c-----------------------------------------------------------------------
24 c Declarations:
25 c -------------
26 
27 #include "dimensions.h"
28 #include "paramet.h"
29 #include "comconst.h"
30 #include "comvert.h"
31 #include "comgeom.h"
32 #include "logic.h"
33 #include "ener.h"
34 
35 c Arguments:
36 c ----------
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 
42 c Local:
43 c ------
44 
45  REAL uav(ip1jmp1,llm),vav(ip1jm,llm),wsur2(ip1jmp1)
46  REAL unsaire2(ip1jmp1), ge(ip1jmp1)
47  REAL deuxjour, ww, gt, uu, vv
48 
49  INTEGER ij,l,ijb,ije
50 
51  EXTERNAL ssum
52  REAL ssum
53 
54 c-----------------------------------------------------------------------
55 c 2. Calculs preliminaires:
56 c -------------------------
57 
58  IF (conser) THEN
59  deuxjour = 2. * daysec
60 
61  DO 1 ij = 1, ip1jmp1
62  unsaire2(ij) = unsaire(ij) * unsaire(ij)
63  1 CONTINUE
64  END IF
65 
66 
67 c------------------ -yy ----------------------------------------------
68 c . Calcul de u
69 
70  DO l=1,llm
71 
72  ijb=ij_begin
73  ije=ij_end
74  if (pole_nord) ijb=ijb+iip1
75  if (pole_sud) ije=ije-iip1
76 
77 c DO ij = iip2, ip1jmp1
78 c uav(ij,l) = 0.25 * ( ucov(ij,l) + ucov(ij-iip1,l) )
79 c ENDDO
80 
81 c DO ij = iip2, ip1jm
82 c uav(ij,l) = uav(ij,l) + uav(ij+iip1,l)
83 c ENDDO
84 
85  DO ij = ijb, ije
86 
87  uav(ij,l)=0.25*(ucov(ij,l)+ucov(ij-iip1,l))
88  . +0.25*(ucov(ij+iip1,l)+ucov(ij,l))
89  ENDDO
90 
91  if (pole_nord) then
92  DO ij = 1, iip1
93  uav(ij ,l) = 0.
94  ENDDO
95  endif
96 
97  if (pole_sud) then
98  DO ij = 1, iip1
99  uav(ip1jm+ij,l) = 0.
100  ENDDO
101  endif
102 
103  ENDDO
104 
105 c call write_field3d_p('uav',reshape(uav,(/iip1,jjp1,llm/)))
106 
107 c------------------ -xx ----------------------------------------------
108 c . Calcul de v
109 
110  ijb=ij_begin
111  ije=ij_end
112  if (pole_sud) ije=ij_end-iip1
113 
114  DO l=1,llm
115 
116  DO ij = ijb+1, ije
117  vav(ij,l) = 0.25 * ( vcov(ij,l) + vcov(ij-1,l) )
118  ENDDO
119 
120  DO ij = ijb,ije,iip1
121  vav(ij,l) = vav(ij+iim,l)
122  ENDDO
123 
124 
125  DO ij = ijb, ije-1
126  vav(ij,l) = vav(ij,l) + vav(ij+1,l)
127  ENDDO
128 
129  DO ij = ijb, ije, iip1
130  vav(ij+iim,l) = vav(ij,l)
131  ENDDO
132 
133  ENDDO
134 c call write_field3d_p('vav',reshape(vav,(/iip1,jjm,llm/)))
135 c-----------------------------------------------------------------------
136 
137 
138 
139  DO 20 l = 1, llmm1
140 
141 
142 c ...... calcul de - w/2. au niveau l+1 .......
143  ijb=ij_begin
144  ije=ij_end+iip1
145  if (pole_sud) ije=ij_end
146 
147  DO 5 ij = ijb, ije
148  wsur2( ij ) = - 0.5 * w( ij,l+1 )
149  5 CONTINUE
150 
151 
152 c ..................... calcul pour du ..................
153 
154  ijb=ij_begin
155  ije=ij_end
156  if (pole_nord) ijb=ijb+iip1
157  if (pole_sud) ije=ije-iip1
158 
159  DO 6 ij = ijb ,ije-1
160  ww = wsur2( ij ) + wsur2( ij+1 )
161  uu = 0.5 * ( ucov(ij,l) + ucov(ij,l+1) )
162  du(ij,l) = du(ij,l) - ww * ( uu - uav(ij, l ) )/massebx(ij, l )
163  du(ij,l+1)= du(ij,l+1) + ww * ( uu - uav(ij,l+1) )/massebx(ij,l+1)
164  6 CONTINUE
165 
166 c ..... correction pour du(iip1,j,l) ........
167 c ..... du(iip1,j,l)= du(1,j,l) .....
168 
169 CDIR$ IVDEP
170  DO 7 ij = ijb+iip1-1, ije, iip1
171  du( ij, l ) = du( ij -iim, l )
172  du( ij,l+1 ) = du( ij -iim,l+1 )
173  7 CONTINUE
174 
175 c ................. calcul pour dv .....................
176  ijb=ij_begin
177  ije=ij_end
178  if (pole_sud) ije=ij_end-iip1
179 
180  DO 8 ij = ijb, ije
181  ww = wsur2( ij+iip1 ) + wsur2( ij )
182  vv = 0.5 * ( vcov(ij,l) + vcov(ij,l+1) )
183  dv(ij,l) = dv(ij, l ) - ww * (vv - vav(ij, l ) )/masseby(ij, l )
184  dv(ij,l+1)= dv(ij,l+1) + ww * (vv - vav(ij,l+1) )/masseby(ij,l+1)
185  8 CONTINUE
186 
187 c
188 
189 c ............................................................
190 c ............... calcul pour dh ...................
191 c ............................................................
192 
193 c ---z
194 c calcul de - d( teta * w ) qu'on ajoute a dh
195 c ...............
196  ijb=ij_begin
197  ije=ij_end
198 
199  DO 15 ij = ijb, ije
200  ww = wsur2(ij) * (teta(ij,l) + teta(ij,l+1) )
201  dteta(ij, l ) = dteta(ij, l ) - ww
202  dteta(ij,l+1) = dteta(ij,l+1) + ww
203  15 CONTINUE
204 
205 c ym ---> conser a voir plus tard
206 
207 c IF( conser) THEN
208 c
209 c DO 17 ij = 1,ip1jmp1
210 c ge(ij) = wsur2(ij) * wsur2(ij) * unsaire2(ij)
211 c 17 CONTINUE
212 c gt = SSUM( ip1jmp1,ge,1 )
213 c gtot(l) = deuxjour * SQRT( gt/ip1jmp1 )
214 c END IF
215 
216  20 CONTINUE
217 
218  RETURN
219  END