source: lmdz_wrf/trunk/WRFV3/lmdz/aaam_bud.F90 @ 409

Last change on this file since 409 was 1, checked in by lfita, 10 years ago
  • -- --- Opening of the WRF+LMDZ coupling repository --- -- -

WRF: version v3.3
LMDZ: version v1818

More details in:

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