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

Annotation of /trunk/phylmd/aaam_bud.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 47 - (hide annotations)
Fri Jul 1 15:00:48 2011 UTC (12 years, 11 months ago) by guez
Original Path: trunk/libf/phylmd/aaam_bud.f
File size: 10737 byte(s)
Split "thermcell.f" and "cv3_routines.f".
Removed copies of files that are now in "L_util".
Moved "mva9" and "diagetpq" to their own files.
Unified variable names across procedures.

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

  ViewVC Help
Powered by ViewVC 1.1.21