LMDZ
top_bound_p.F
Go to the documentation of this file.
1 !
2 ! $Id: top_bound_p.F 1907 2013-11-26 13:10:46Z lguez $
3 !
4  SUBROUTINE top_bound_p(vcov,ucov,teta,masse,dt)
5  USE parallel_lmdz
6  IMPLICIT NONE
7 c
8 #include "dimensions.h"
9 #include "paramet.h"
10 #include "comconst.h"
11 #include "comvert.h"
12 #include "comgeom2.h"
13 
14 
15 c .. DISSIPATION LINEAIRE A HAUT NIVEAU, RUN MESO,
16 C F. LOTT DEC. 2006
17 c ( 10/12/06 )
18 
19 c=======================================================================
20 c
21 c Auteur: F. LOTT
22 c -------
23 c
24 c Objet:
25 c ------
26 c
27 c Dissipation linéaire (ex top_bound de la physique)
28 c
29 c=======================================================================
30 
31 ! top_bound sponge layer model:
32 ! Quenching is modeled as: A(t)=Am+A0*exp(-lambda*t)
33 ! where Am is the zonal average of the field (or zero), and lambda the inverse
34 ! of the characteristic quenching/relaxation time scale
35 ! Thus, assuming Am to be time-independent, field at time t+dt is given by:
36 ! A(t+dt)=A(t)-(A(t)-Am)*(1-exp(-lambda*t))
37 ! Moreover lambda can be a function of model level (see below), and relaxation
38 ! can be toward the average zonal field or just zero (see below).
39 
40 ! NB: top_bound sponge is only called from leapfrog if ok_strato=.true.
41 
42 ! sponge parameters: (loaded/set in conf_gcm.F ; stored in comconst.h)
43 ! iflag_top_bound=0 for no sponge
44 ! iflag_top_bound=1 for sponge over 4 topmost layers
45 ! iflag_top_bound=2 for sponge from top to ~1% of top layer pressure
46 ! mode_top_bound=0: no relaxation
47 ! mode_top_bound=1: u and v relax towards 0
48 ! mode_top_bound=2: u and v relax towards their zonal mean
49 ! mode_top_bound=3: u,v and pot. temp. relax towards their zonal mean
50 ! tau_top_bound : inverse of charactericstic relaxation time scale at
51 ! the topmost layer (Hz)
52 
53 
54 #include "comdissipn.h"
55 #include "iniprint.h"
56 
57 c Arguments:
58 c ----------
59 
60  real,intent(inout) :: ucov(iip1,jjp1,llm) ! covariant zonal wind
61  real,intent(inout) :: vcov(iip1,jjm,llm) ! covariant meridional wind
62  real,intent(inout) :: teta(iip1,jjp1,llm) ! potential temperature
63  real,intent(in) :: masse(iip1,jjp1,llm) ! mass of atmosphere
64  real,intent(in) :: dt ! time step (s) of sponge model
65 
66 c Local:
67 c ------
68  REAL massebx(iip1,jjp1,llm),masseby(iip1,jjm,llm),zm
69  REAL uzon(jjp1,llm),vzon(jjm,llm),tzon(jjp1,llm)
70 
71  integer i
72  REAL,SAVE :: rdamp(llm) ! quenching coefficient
73  real,save :: lambda(llm) ! inverse or quenching time scale (Hz)
74  LOGICAL,SAVE :: first=.true.
75  INTEGER j,l,jjb,jje
76 
77 
78  if (iflag_top_bound == 0) return
79 
80  if (first) then
81 c$OMP BARRIER
82 c$OMP MASTER
83  if (iflag_top_bound == 1) then
84 ! sponge quenching over the topmost 4 atmospheric layers
85  lambda(:)=0.
86  lambda(llm)=tau_top_bound
87  lambda(llm-1)=tau_top_bound/2.
88  lambda(llm-2)=tau_top_bound/4.
89  lambda(llm-3)=tau_top_bound/8.
90  else if (iflag_top_bound == 2) then
91 ! sponge quenching over topmost layers down to pressures which are
92 ! higher than 100 times the topmost layer pressure
93  lambda(:)=tau_top_bound
94  s *max(presnivs(llm)/presnivs(:)-0.01,0.)
95  endif
96 
97 ! quenching coefficient rdamp(:)
98 ! rdamp(:)=dt*lambda(:) ! Explicit Euler approx.
99  rdamp(:)=1.-exp(-lambda(:)*dt)
100 
101  write(lunout,*)'TOP_BOUND mode',mode_top_bound
102  write(lunout,*)'Sponge layer coefficients'
103  write(lunout,*)'p (Pa) z(km) tau(s) 1./tau (Hz)'
104  do l=1,llm
105  if (rdamp(l).ne.0.) then
106  write(lunout,'(6(1pe12.4,1x))')
107  & presnivs(l),log(preff/presnivs(l))*scaleheight,
108  & 1./lambda(l),lambda(l)
109  endif
110  enddo
111  first=.false.
112 c$OMP END MASTER
113 c$OMP BARRIER
114  endif ! of if (first)
115 
116 
117  CALL massbar_p(masse,massebx,masseby)
118 
119  ! compute zonal average of vcov (or set it to zero)
120  if (mode_top_bound.ge.2) then
121  jjb=jj_begin
122  jje=jj_end
123  IF (pole_sud) jje=jj_end-1
124 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
125  do l=1,llm
126  do j=jjb,jje
127  zm=0.
128  vzon(j,l)=0
129  do i=1,iim
130 ! NB: we can work using vcov zonal mean rather than v since the
131 ! cv coefficient (which relates the two) only varies with latitudes
132  vzon(j,l)=vzon(j,l)+vcov(i,j,l)*masseby(i,j,l)
133  zm=zm+masseby(i,j,l)
134  enddo
135  vzon(j,l)=vzon(j,l)/zm
136  enddo
137  enddo
138 c$OMP END DO NOWAIT
139  else
140 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
141  do l=1,llm
142  vzon(:,l)=0.
143  enddo
144 c$OMP END DO NOWAIT
145  endif ! of if (mode_top_bound.ge.2)
146 
147  ! compute zonal average of u (or set it to zero)
148  if (mode_top_bound.ge.2) then
149  jjb=jj_begin
150  jje=jj_end
151  IF (pole_nord) jjb=jj_begin+1
152  IF (pole_sud) jje=jj_end-1
153 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
154  do l=1,llm
155  do j=jjb,jje
156  uzon(j,l)=0.
157  zm=0.
158  do i=1,iim
159  uzon(j,l)=uzon(j,l)+massebx(i,j,l)*ucov(i,j,l)/cu(i,j)
160  zm=zm+massebx(i,j,l)
161  enddo
162  uzon(j,l)=uzon(j,l)/zm
163  enddo
164  enddo
165 c$OMP END DO NOWAIT
166  else
167 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
168  do l=1,llm
169  uzon(:,l)=0.
170  enddo
171 c$OMP END DO NOWAIT
172  endif ! of if (mode_top_bound.ge.2)
173 
174  ! compute zonal average of potential temperature, if necessary
175  if (mode_top_bound.ge.3) then
176  jjb=jj_begin
177  jje=jj_end
178  IF (pole_nord) jjb=jj_begin+1
179  IF (pole_sud) jje=jj_end-1
180 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
181  do l=1,llm
182  do j=jjb,jje
183  zm=0.
184  tzon(j,l)=0.
185  do i=1,iim
186  tzon(j,l)=tzon(j,l)+teta(i,j,l)*masse(i,j,l)
187  zm=zm+masse(i,j,l)
188  enddo
189  tzon(j,l)=tzon(j,l)/zm
190  enddo
191  enddo
192 c$OMP END DO NOWAIT
193  endif ! of if (mode_top_bound.ge.3)
194 
195  if (mode_top_bound.ge.1) then
196  ! Apply sponge quenching on vcov:
197  jjb=jj_begin
198  jje=jj_end
199  IF (pole_sud) jje=jj_end-1
200 
201 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
202  do l=1,llm
203  do j=jjb,jje
204  do i=1,iip1
205  vcov(i,j,l)=vcov(i,j,l)
206  & -rdamp(l)*(vcov(i,j,l)-vzon(j,l))
207  enddo
208  enddo
209  enddo
210 c$OMP END DO NOWAIT
211 
212  ! Apply sponge quenching on ucov:
213  jjb=jj_begin
214  jje=jj_end
215  IF (pole_nord) jjb=jj_begin+1
216  IF (pole_sud) jje=jj_end-1
217 
218 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
219  do l=1,llm
220  do j=jjb,jje
221  do i=1,iip1
222  ucov(i,j,l)=ucov(i,j,l)
223  & -rdamp(l)*(ucov(i,j,l)-cu(i,j)*uzon(j,l))
224  enddo
225  enddo
226  enddo
227 c$OMP END DO NOWAIT
228  endif ! of if (mode_top_bound.ge.1)
229 
230  if (mode_top_bound.ge.3) then
231  ! Apply sponge quenching on teta:
232  jjb=jj_begin
233  jje=jj_end
234  IF (pole_nord) jjb=jj_begin+1
235  IF (pole_sud) jje=jj_end-1
236 
237 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
238  do l=1,llm
239  do j=jjb,jje
240  do i=1,iip1
241  teta(i,j,l)=teta(i,j,l)
242  & -rdamp(l)*(teta(i,j,l)-tzon(j,l))
243  enddo
244  enddo
245  enddo
246 c$OMP END DO NOWAIT
247  endif ! of if (mode_top_bond.ge.3)
248 
249  END
preff
!$Id preff
Definition: comvert.h:8
parallel_lmdz::jj_end
integer, save jj_end
Definition: parallel_lmdz.F90:20
iflag_top_bound
!$Id && iflag_top_bound
Definition: comconst.h:7
parallel_lmdz::jj_begin
integer, save jj_begin
Definition: parallel_lmdz.F90:19
parallel_lmdz::pole_sud
logical, save pole_sud
Definition: parallel_lmdz.F90:25
llm
!$Id Turb_fcg_gcssold get_uvd hqturb_gcssold endif!large scale llm day day1 day day1 *dt_toga endif!time annee_ref dt_toga u_toga vq_toga w_prof vq_prof llm day day1 day day1 *dt_dice endif!time annee_ref dt_dice swup_dice vg_dice omega_dice tg_prof vg_profd w_profd omega_profd!do llm!print llm l llm
Definition: 1D_interp_cases.h:103
presnivs
!$Id presnivs(llm)
tau_top_bound
!$Id mode_top_bound COMMON comconstr omeg dissip_zref tau_top_bound
Definition: comconst.h:7
false
!$Id itapm1 ENDIF!IM on interpole les champs sur les niveaux STD de pression!IM a chaque pas de temps de la physique c!positionnement de l argument logique a false c!pour ne pas recalculer deux fois la meme chose!c!a cet effet un appel a plevel_new a ete deplace c!a la fin de la serie d appels c!la boucle DO nlevSTD a ete internalisee c!dans d ou la creation de cette routine c c!CALL false
Definition: calcul_STDlev.h:26
massbar_p
subroutine massbar_p(masse, massebx, masseby)
Definition: massbar_p.F:2
parallel_lmdz::pole_nord
logical, save pole_nord
Definition: parallel_lmdz.F90:24
parallel_lmdz
Definition: parallel_lmdz.F90:4
jjp1
!$Header jjp1
Definition: paramet.h:14
top_bound_p
subroutine top_bound_p(vcov, ucov, teta, masse, dt)
Definition: top_bound_p.F:5
true
!$Id itapm1 ENDIF!IM on interpole les champs sur les niveaux STD de pression!IM a chaque pas de temps de la physique c!positionnement de l argument logique a false c!pour ne pas recalculer deux fois la meme chose!c!a cet effet un appel a plevel_new a ete deplace c!a la fin de la serie d appels c!la boucle DO nlevSTD a ete internalisee c!dans d ou la creation de cette routine c c!CALL ulevSTD CALL &zphi philevSTD CALL &zx_rh rhlevSTD!DO klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon du jour ou toutes les read_climoz CALL true
Definition: calcul_STDlev.h:139
scaleheight
!$Id scaleheight
Definition: comvert.h:9
cu
!$Id!Parameters for parameters that control the rate of approach!to quasi equilibrium noff nlm real tlcrit real entp real sigd real coeffs real dtmax real cu real betad real damp real delta COMMON cvparam nlm tlcrit sigd coeffs cu
Definition: cvparam.h:12
iim
c c zjulian c cym CALL iim cym klev iim
Definition: ini_bilKP_ave.h:24
lunout
!$Header!gestion des impressions de sorties et de débogage la sortie standard prt_level COMMON comprint lunout
Definition: iniprint.h:7