My Project
 All Classes Files Functions Variables Macros
sw_case_williamson91_6_loc.F
Go to the documentation of this file.
1 !
2 ! $Id $
3 !
4  SUBROUTINE sw_case_williamson91_6_loc(vcov,ucov,teta,masse,ps)
5 
6 c=======================================================================
7 c
8 c Author: Thomas Dubos original: 26/01/2010
9 c -------
10 c
11 c Subject:
12 c ------
13 c Realise le cas-test 6 de Williamson et al. (1991) : onde de Rossby-Haurwitz
14 c
15 c Method:
16 c --------
17 c
18 c Interface:
19 c ----------
20 c
21 c Input:
22 c ------
23 c
24 c Output:
25 c -------
26 c
27 c=======================================================================
28  USE parallel
29 
30  IMPLICIT NONE
31 c-----------------------------------------------------------------------
32 c Declararations:
33 c ---------------
34 
35 #include "dimensions.h"
36 #include "paramet.h"
37 #include "comvert.h"
38 #include "comconst.h"
39 #include "comgeom.h"
40 #include "iniprint.h"
41 
42 c Arguments:
43 c ----------
44 
45 c variables dynamiques
46  REAL vcov(ijb_v:ije_v,llm),ucov(ijb_u:ije_u,llm) ! vents covariants
47  REAL teta(ijb_u:ije_u,llm) ! temperature potentielle
48  REAL ps(ijb_u:ije_u) ! pression au sol
49  REAL masse(ijb_u:ije_u,llm) ! masse d'air
50  REAL phis(ijb_u:ije_u) ! geopotentiel au sol
51 
52 c Local:
53 c ------
54 
55  real,allocatable :: ucov_glo(:,:)
56  real,allocatable :: vcov_glo(:,:)
57  real,allocatable :: teta_glo(:,:)
58  real,allocatable :: masse_glo(:,:)
59  real,allocatable :: ps_glo(:)
60 
61 ! REAL p (ip1jmp1,llmp1 ) ! pression aux interfac.des couches
62 ! REAL pks(ip1jmp1) ! exner au sol
63 ! REAL pk(ip1jmp1,llm) ! exner au milieu des couches
64 ! REAL pkf(ip1jmp1,llm) ! exner filt.au milieu des couches
65 ! REAL alpha(ip1jmp1,llm),beta(ip1jmp1,llm)
66 
67  real,allocatable :: p(:,:)
68  real,allocatable :: pks(:)
69  real,allocatable :: pk(:,:)
70  real,allocatable :: pkf(:,:)
71  real,allocatable :: alpha(:,:),beta(:,:)
72 
73  REAL :: sinth,costh,costh2, ath,bth,cth, lon,dps
74  INTEGER i,j,ij
75 
76  REAL, PARAMETER :: rho=1 ! masse volumique de l'air (arbitraire)
77  REAL, PARAMETER :: k = 7.848e-6 ! K = \omega
78  REAL, PARAMETER :: gh0 = 9.80616 * 8e3
79  INTEGER, PARAMETER :: r0=4, r1=r0+1, r2=r0+2 ! mode 4
80 c NB : rad = 6371220 dans W91 (6371229 dans LMDZ)
81 c omeg = 7.292e-5 dans W91 (7.2722e-5 dans LMDZ)
82 
83 
84  ! allocate (global) arrays
85  allocate(vcov_glo(ip1jm,llm))
86  allocate(ucov_glo(ip1jmp1,llm))
87  allocate(teta_glo(ip1jmp1,llm))
88  allocate(ps_glo(ip1jmp1))
89  allocate(masse_glo(ip1jmp1,llm))
90 
91  allocate(p(ip1jmp1,llmp1))
92  allocate(pks(ip1jmp1))
93  allocate(pk(ip1jmp1,llm))
94  allocate(pkf(ip1jmp1,llm))
95  allocate(alpha(ip1jmp1,llm))
96  allocate(beta(ip1jmp1,llm))
97 
98  IF(0==0) THEN
99 !c Williamson et al. (1991) : onde de Rossby-Haurwitz
100  teta_glo(:,:) = preff/rho/cpp
101 !c geopotentiel (pression de surface)
102  do j=1,jjp1
103  costh2 = cos(rlatu(j))**2
104  ath = (r0+1)*(costh2**2) + (2*r0*r0-r0-2)*costh2 - 2*r0*r0
105  ath = .25*(k**2)*(costh2**(r0-1))*ath
106  ath = .5*k*(2*omeg+k)*costh2 + ath
107  bth = (r1*r1+1)-r1*r1*costh2
108  bth = 2*(omeg+k)*k/(r1*r2) * (costh2**(r0/2))*bth
109  cth = r1*costh2 - r2
110  cth = .25*k*k*(costh2**r0)*cth
111  do i=1,iip1
112  ij=(j-1)*iip1+i
113  lon = rlonv(i)
114  dps = ath + bth*cos(r0*lon) + cth*cos(2*r0*lon)
115  ps_glo(ij) = rho*(gh0 + (rad**2)*dps)
116  enddo
117  enddo
118 ! write(lunout,*) 'W91 ps', MAXVAL(ps), MINVAL(ps)
119 c vitesse zonale ucov
120  do j=1,jjp1
121  costh = cos(rlatu(j))
122  costh2 = costh**2
123  ath = rad*k*costh
124  bth = r0*(1-costh2)-costh2
125  bth = rad*k*bth*(costh**(r0-1))
126  do i=1,iip1
127  ij=(j-1)*iip1+i
128  lon = rlonu(i)
129  ucov_glo(ij,1) = (ath + bth*cos(r0*lon))
130  enddo
131  enddo
132 ! write(lunout,*) 'W91 u', MAXVAL(ucov(:,1)), MINVAL(ucov(:,1))
133  ucov_glo(:,1)=ucov_glo(:,1)*cu
134 c vitesse meridienne vcov
135  do j=1,jjm
136  sinth = sin(rlatv(j))
137  costh = cos(rlatv(j))
138  ath = -rad*k*r0*sinth*(costh**(r0-1))
139  do i=1,iip1
140  ij=(j-1)*iip1+i
141  lon = rlonv(i)
142  vcov_glo(ij,1) = ath*sin(r0*lon)
143  enddo
144  enddo
145  write(lunout,*) 'W91 v', maxval(vcov(:,1)), minval(vcov(:,1))
146  vcov_glo(:,1)=vcov_glo(:,1)*cv
147 
148 c ucov_glo=0
149 c vcov_glo=0
150  ELSE
151 c test non-tournant, onde se propageant en latitude
152  do j=1,jjp1
153  do i=1,iip1
154  ij=(j-1)*iip1+i
155  ps_glo(ij) = 1e5*(1 + .1*exp(-100*(1+sin(rlatu(j)))**2))
156  enddo
157  enddo
158 
159 c rho = preff/(cpp*teta)
160  teta_glo(:,:) = .01*preff/cpp ! rho = 100 ; phi = ps/rho = 1e3 ; c=30 m/s = 2600 km/j = 23 degres / j
161  ucov_glo(:,:)=0.
162  vcov_glo(:,:)=0.
163  END IF
164 
165  CALL pression( ip1jmp1, ap, bp, ps_glo, p )
166  CALL massdair(p,masse_glo)
167 
168  ! copy data from global array to local array:
169  teta(ijb_u:ije_u,:)=teta_glo(ijb_u:ije_u,:)
170  ucov(ijb_u:ije_u,:)=ucov_glo(ijb_u:ije_u,:)
171  vcov(ijb_v:ije_v,:)=vcov_glo(ijb_v:ije_v,:)
172  masse(ijb_u:ije_u,:)=masse_glo(ijb_u:ije_u,:)
173  ps(ijb_u:ije_u)=ps_glo(ijb_u:ije_u)
174 
175  ! cleanup
176  deallocate(teta_glo)
177  deallocate(ucov_glo)
178  deallocate(vcov_glo)
179  deallocate(masse_glo)
180  deallocate(ps_glo)
181  deallocate(p)
182  deallocate(pks)
183  deallocate(pk)
184  deallocate(pkf)
185  deallocate(alpha)
186  deallocate(beta)
187 
188  END
189 c-----------------------------------------------------------------------