/[lmdze]/trunk/libf/phylmd/aaam_bud.f
ViewVC logotype

Contents of /trunk/libf/phylmd/aaam_bud.f

Parent Directory Parent Directory | Revision Log Revision Log

Revision 7 - (show annotations)
Mon Mar 31 12:24:17 2008 UTC (16 years, 1 month ago) by guez
File size: 10729 byte(s) 
This revision is not in working order. Pending some moving of files.

Important changes. In the program "etat0_lim": ozone coefficients from
Mobidic are regridded in time instead of pressure ; consequences in
"etat0". In the program "gcm", ozone coefficients from Mobidic are
read once per day only for the current day and regridded in pressure ;
consequences in "o3_chem_m", "regr_pr_coefoz", "phytrac" and
"regr_pr_comb_coefoz_m".

NetCDF95 is a library and does not export NetCDF.

New variables "nag_gl_options", "nag_fcalls_options" and
"nag_cross_options" in "nag_tools.mk".

"check_coefoz.jnl" rewritten entirely for new version of
"coefoz_LMDZ.nc".

Target "obj_etat0_lim" moved from "GNUmakefile" to "nag_rules.mk".

Added some "intent" attributes in "calfis", "clmain", "clqh",
"cltrac", "cltracrn", "cvltr", "ini_undefSTD", "moy_undefSTD",
"nflxtr", "phystokenc", "phytrac", "readsulfate", "readsulfate_preind"
and "undefSTD".

In "dynetat0", "dynredem0" and "gcm", "phis" has rank 2 instead of
1. "phis" has assumed shape in "dynredem0".

Added module containing "dynredem0". Changed some calls with NetCDF
Fortran 77 interface to calls with NetCDF95 interface.

Replaced calls to "ssum" by calls to "sum" in "inigeom".

In "make.sh", new option "-c" to change compiler.

In "aaam_bud", argument "rjour" deleted.

In "physiq": renamed some variables; deleted variable "xjour".

In "phytrac": renamed some variables; new argument "lmt_pas".
1 subroutine aaam_bud (iam,nlon,nlev,rsec,
2 i rea,rg,ome,
3 i plat,plon,phis,
4 i dragu,liftu,phyu,
5 i dragv,liftv,phyv,
6 i p, u, v,
7 o aam, torsfc)
8 c
9 use dimens_m
10 use dimphy
11 implicit none
12 c======================================================================
13 c Auteur(s): F.Lott (LMD/CNRS) date: 20031020
14 c Object: Compute different terms of the axial AAAM Budget.
15 C No outputs, every AAM quantities are written on the IAM
16 C File.
17 c
18 c Modif : I.Musat (LMD/CNRS) date : 20041020
19 c Outputs : axial components of wind AAM "aam" and total surface torque "torsfc",
20 c but no write in the iam file.
21 c
22 C WARNING: Only valid for regular rectangular grids.
23 C REMARK: CALL DANS PHYSIQ AFTER lift_noro:
24 C CALL aaam_bud (27,klon,klev,rjourvrai,gmtime,
25 C C ra,rg,romega,
26 C C rlat,rlon,pphis,
27 C C zustrdr,zustrli,zustrph,
28 C C zvstrdr,zvstrli,zvstrph,
29 C C paprs,u,v)
30 C
31 C======================================================================
32 c Explicit Arguments:
33 c ==================
34 c iam-----input-I-File number where AAMs and torques are written
35 c It is a formatted file that has been opened
36 c in physiq.F
37 c nlon----input-I-Total number of horizontal points that get into physics
38 c nlev----input-I-Number of vertical levels
39 c rsec----input-R-Seconde de la journee
40 c rea-----input-R-Earth radius
41 c rg------input-R-gravity constant
42 c ome-----input-R-Earth rotation rate
43 c plat ---input-R-Latitude en degres
44 c plon ---input-R-Longitude en degres
45 c phis ---input-R-Geopotential at the ground
46 c dragu---input-R-orodrag stress (zonal)
47 c liftu---input-R-orolift stress (zonal)
48 c phyu----input-R-Stress total de la physique (zonal)
49 c dragv---input-R-orodrag stress (Meridional)
50 c liftv---input-R-orolift stress (Meridional)
51 c phyv----input-R-Stress total de la physique (Meridional)
52 c p-------input-R-Pressure (Pa) at model half levels
53 c u-------input-R-Horizontal wind (m/s)
54 c v-------input-R-Meridional wind (m/s)
55 c aam-----output-R-Axial Wind AAM (=raam(3))
56 c torsfc--output-R-Total surface torque (=tmou(3)+tsso(3)+tbls(3))
57 c
58 c Implicit Arguments:
59 c ===================
60 c
61 c iim--common-I: Number of longitude intervals
62 c jjm--common-I: Number of latitude intervals
63 c klon-common-I: Number of points seen by the physics
64 c iim*(jjm-1)+2 for instance
65 c klev-common-I: Number of vertical layers
66 c======================================================================
67 c Local Variables:
68 c ================
69 c dlat-----R: Latitude increment (Radians)
70 c dlon-----R: Longitude increment (Radians)
71 c raam ---R: Wind AAM (3 Components, 1 & 2 Equatoriales; 3 Axiale)
72 c oaam ---R: Mass AAM (3 Components, 1 & 2 Equatoriales; 3 Axiale)
73 c tmou-----R: Resolved Mountain torque (3 components)
74 c tsso-----R: Parameterised Moutain drag torque (3 components)
75 c tbls-----R: Parameterised Boundary layer torque (3 components)
76 c
77 c LOCAL ARRAY:
78 c ===========
79 c zs ---R: Topographic height
80 c ps ---R: Surface Pressure
81 c ub ---R: Barotropic wind zonal
82 c vb ---R: Barotropic wind meridional
83 c zlat ---R: Latitude in radians
84 c zlon ---R: Longitude in radians
85 c======================================================================
86
87 c
88 c ARGUMENTS
89 c
90 INTEGER iam,nlon,nlev
91 real, intent(in):: rsec
92 real rea
93 real, intent(in):: rg
94 real ome
95 REAL, intent(in):: plat(nlon),plon(nlon)
96 real phis(nlon)
97 REAL dragu(nlon),liftu(nlon),phyu(nlon)
98 REAL dragv(nlon),liftv(nlon),phyv(nlon)
99 REAL, intent(in):: p(nlon,nlev+1)
100 real u(nlon,nlev), v(nlon,nlev)
101 c
102 c Variables locales:
103 c
104 INTEGER i,j,k,l
105 REAL xpi,hadley,hadday
106 REAL dlat,dlon
107 REAL raam(3),oaam(3),tmou(3),tsso(3),tbls(3)
108 integer iax
109 cIM ajout aam, torsfc
110 c aam = composante axiale du Wind AAM raam
111 c torsfc = composante axiale de (tmou+tsso+tbls)
112 REAL aam, torsfc
113
114 REAL ZS(801,401),PS(801,401)
115 REAL UB(801,401),VB(801,401)
116 REAL SSOU(801,401),SSOV(801,401)
117 REAL BLSU(801,401),BLSV(801,401)
118 REAL ZLON(801),ZLAT(401)
119 C
120 C PUT AAM QUANTITIES AT ZERO:
121 C
122 if(iim+1.gt.801.or.jjm+1.gt.401)then
123 print *,' Pb de dimension dans aaam_bud'
124 stop
125 endif
126
127 xpi=acos(-1.)
128 hadley=1.e18
129 hadday=1.e18*24.*3600.
130 dlat=xpi/float(jjm)
131 dlon=2.*xpi/float(iim)
132
133 do iax=1,3
134 oaam(iax)=0.
135 raam(iax)=0.
136 tmou(iax)=0.
137 tsso(iax)=0.
138 tbls(iax)=0.
139 enddo
140
141 C MOUNTAIN HEIGHT, PRESSURE AND BAROTROPIC WIND:
142
143 C North pole values (j=1):
144
145 l=1
146
147 ub(1,1)=0.
148 vb(1,1)=0.
149 do k=1,nlev
150 ub(1,1)=ub(1,1)+u(l,k)*(p(l,k)-p(l,k+1))/rg
151 vb(1,1)=vb(1,1)+v(l,k)*(p(l,k)-p(l,k+1))/rg
152 enddo
153
154 zlat(1)=plat(l)*xpi/180.
155
156 do i=1,iim+1
157
158 zs(i,1)=phis(l)/rg
159 ps(i,1)=p(l,1)
160 ub(i,1)=ub(1,1)
161 vb(i,1)=vb(1,1)
162 ssou(i,1)=dragu(l)+liftu(l)
163 ssov(i,1)=dragv(l)+liftv(l)
164 blsu(i,1)=phyu(l)-dragu(l)-liftu(l)
165 blsv(i,1)=phyv(l)-dragv(l)-liftv(l)
166
167 enddo
168
169
170 do j = 2,jjm
171
172 C Values at Greenwich (Periodicity)
173
174 zs(iim+1,j)=phis(l+1)/rg
175 ps(iim+1,j)=p(l+1,1)
176 ssou(iim+1,j)=dragu(l+1)+liftu(l+1)
177 ssov(iim+1,j)=dragv(l+1)+liftv(l+1)
178 blsu(iim+1,j)=phyu(l+1)-dragu(l+1)-liftu(l+1)
179 blsv(iim+1,j)=phyv(l+1)-dragv(l+1)-liftv(l+1)
180 zlon(iim+1)=-plon(l+1)*xpi/180.
181 zlat(j)=plat(l+1)*xpi/180.
182
183 ub(iim+1,j)=0.
184 vb(iim+1,j)=0.
185 do k=1,nlev
186 ub(iim+1,j)=ub(iim+1,j)+u(l+1,k)*(p(l+1,k)-p(l+1,k+1))/rg
187 vb(iim+1,j)=vb(iim+1,j)+v(l+1,k)*(p(l+1,k)-p(l+1,k+1))/rg
188 enddo
189
190
191 do i=1,iim
192
193 l=l+1
194 zs(i,j)=phis(l)/rg
195 ps(i,j)=p(l,1)
196 ssou(i,j)=dragu(l)+liftu(l)
197 ssov(i,j)=dragv(l)+liftv(l)
198 blsu(i,j)=phyu(l)-dragu(l)-liftu(l)
199 blsv(i,j)=phyv(l)-dragv(l)-liftv(l)
200 zlon(i)=plon(l)*xpi/180.
201
202 ub(i,j)=0.
203 vb(i,j)=0.
204 do k=1,nlev
205 ub(i,j)=ub(i,j)+u(l,k)*(p(l,k)-p(l,k+1))/rg
206 vb(i,j)=vb(i,j)+v(l,k)*(p(l,k)-p(l,k+1))/rg
207 enddo
208
209 enddo
210
211 enddo
212
213
214 C South Pole
215
216 l=l+1
217 ub(1,jjm+1)=0.
218 vb(1,jjm+1)=0.
219 do k=1,nlev
220 ub(1,jjm+1)=ub(1,jjm+1)+u(l,k)*(p(l,k)-p(l,k+1))/rg
221 vb(1,jjm+1)=vb(1,jjm+1)+v(l,k)*(p(l,k)-p(l,k+1))/rg
222 enddo
223 zlat(jjm+1)=plat(l)*xpi/180.
224
225 do i=1,iim+1
226 zs(i,jjm+1)=phis(l)/rg
227 ps(i,jjm+1)=p(l,1)
228 ssou(i,jjm+1)=dragu(l)+liftu(l)
229 ssov(i,jjm+1)=dragv(l)+liftv(l)
230 blsu(i,jjm+1)=phyu(l)-dragu(l)-liftu(l)
231 blsv(i,jjm+1)=phyv(l)-dragv(l)-liftv(l)
232 ub(i,jjm+1)=ub(1,jjm+1)
233 vb(i,jjm+1)=vb(1,jjm+1)
234 enddo
235
236 C
237 C MOMENT ANGULAIRE
238 C
239 DO j=1,jjm
240 DO i=1,iim
241
242 raam(1)=raam(1)-rea**3*dlon*dlat*0.5*
243 c (cos(zlon(i ))*sin(zlat(j ))*cos(zlat(j ))*ub(i ,j )
244 c +cos(zlon(i ))*sin(zlat(j+1))*cos(zlat(j+1))*ub(i ,j+1))
245 c +rea**3*dlon*dlat*0.5*
246 c (sin(zlon(i ))*cos(zlat(j ))*vb(i ,j )
247 c +sin(zlon(i ))*cos(zlat(j+1))*vb(i ,j+1))
248
249 oaam(1)=oaam(1)-ome*rea**4*dlon*dlat/rg*0.5*
250 c (cos(zlon(i ))*cos(zlat(j ))**2*sin(zlat(j ))*ps(i ,j )
251 c +cos(zlon(i ))*cos(zlat(j+1))**2*sin(zlat(j+1))*ps(i ,j+1))
252
253 raam(2)=raam(2)-rea**3*dlon*dlat*0.5*
254 c (sin(zlon(i ))*sin(zlat(j ))*cos(zlat(j ))*ub(i ,j )
255 c +sin(zlon(i ))*sin(zlat(j+1))*cos(zlat(j+1))*ub(i ,j+1))
256 c -rea**3*dlon*dlat*0.5*
257 c (cos(zlon(i ))*cos(zlat(j ))*vb(i ,j )
258 c +cos(zlon(i ))*cos(zlat(j+1))*vb(i ,j+1))
259
260 oaam(2)=oaam(2)-ome*rea**4*dlon*dlat/rg*0.5*
261 c (sin(zlon(i ))*cos(zlat(j ))**2*sin(zlat(j ))*ps(i ,j )
262 c +sin(zlon(i ))*cos(zlat(j+1))**2*sin(zlat(j+1))*ps(i ,j+1))
263
264 raam(3)=raam(3)+rea**3*dlon*dlat*0.5*
265 c (cos(zlat(j))**2*ub(i,j)+cos(zlat(j+1))**2*ub(i,j+1))
266
267 oaam(3)=oaam(3)+ome*rea**4*dlon*dlat/rg*0.5*
268 c (cos(zlat(j))**3*ps(i,j)+cos(zlat(j+1))**3*ps(i,j+1))
269
270 ENDDO
271 ENDDO
272
273 C
274 C COUPLE DES MONTAGNES:
275 C
276
277 DO j=1,jjm
278 DO i=1,iim
279 tmou(1)=tmou(1)-rea**2*dlon*0.5*sin(zlon(i))
280 c *(zs(i,j)-zs(i,j+1))
281 c *(cos(zlat(j+1))*ps(i,j+1)+cos(zlat(j))*ps(i,j))
282 tmou(2)=tmou(2)+rea**2*dlon*0.5*cos(zlon(i))
283 c *(zs(i,j)-zs(i,j+1))
284 c *(cos(zlat(j+1))*ps(i,j+1)+cos(zlat(j))*ps(i,j))
285 ENDDO
286 ENDDO
287
288 DO j=2,jjm
289 DO i=1,iim
290 tmou(1)=tmou(1)+rea**2*dlat*0.5*sin(zlat(j))
291 c *(zs(i+1,j)-zs(i,j))
292 c *(cos(zlon(i+1))*ps(i+1,j)+cos(zlon(i))*ps(i,j))
293 tmou(2)=tmou(2)+rea**2*dlat*0.5*sin(zlat(j))
294 c *(zs(i+1,j)-zs(i,j))
295 c *(sin(zlon(i+1))*ps(i+1,j)+sin(zlon(i))*ps(i,j))
296 tmou(3)=tmou(3)-rea**2*dlat*0.5*
297 c cos(zlat(j))*(zs(i+1,j)-zs(i,j))*(ps(i+1,j)+ps(i,j))
298 ENDDO
299 ENDDO
300
301 C
302 C COUPLES DES DIFFERENTES FRICTION AU SOL:
303 C
304 l=1
305 DO j=2,jjm
306 DO i=1,iim
307 l=l+1
308 tsso(1)=tsso(1)-rea**3*cos(zlat(j))*dlon*dlat*
309 c ssou(i,j) *sin(zlat(j))*cos(zlon(i))
310 c +rea**3*cos(zlat(j))*dlon*dlat*
311 c ssov(i,j) *sin(zlon(i))
312
313 tsso(2)=tsso(2)-rea**3*cos(zlat(j))*dlon*dlat*
314 c ssou(i,j) *sin(zlat(j))*sin(zlon(i))
315 c -rea**3*cos(zlat(j))*dlon*dlat*
316 c ssov(i,j) *cos(zlon(i))
317
318 tsso(3)=tsso(3)+rea**3*cos(zlat(j))*dlon*dlat*
319 c ssou(i,j) *cos(zlat(j))
320
321 tbls(1)=tbls(1)-rea**3*cos(zlat(j))*dlon*dlat*
322 c blsu(i,j) *sin(zlat(j))*cos(zlon(i))
323 c +rea**3*cos(zlat(j))*dlon*dlat*
324 c blsv(i,j) *sin(zlon(i))
325
326 tbls(2)=tbls(2)-rea**3*cos(zlat(j))*dlon*dlat*
327 c blsu(i,j) *sin(zlat(j))*sin(zlon(i))
328 c -rea**3*cos(zlat(j))*dlon*dlat*
329 c blsv(i,j) *cos(zlon(i))
330
331 tbls(3)=tbls(3)+rea**3*cos(zlat(j))*dlon*dlat*
332 c blsu(i,j) *cos(zlat(j))
333
334 ENDDO
335 ENDDO
336
337
338 100 format(F12.5,15(1x,F12.5))
339
340 aam=raam(3)
341 torsfc= tmou(3)+tsso(3)+tbls(3)
342 c
343 RETURN
344 END
  ViewVC Help
Powered by ViewVC 1.1.21