source: palm/trunk/SOURCE/subsidence.f90 @ 1660

Last change on this file since 1660 was 1490, checked in by raasch, 10 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 12.2 KB
Line 
1 MODULE subsidence_mod
2
3!--------------------------------------------------------------------------------!
4! This file is part of PALM.
5!
6! PALM is free software: you can redistribute it and/or modify it under the terms
7! of the GNU General Public License as published by the Free Software Foundation,
8! either version 3 of the License, or (at your option) any later version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 1997-2014 Leibniz Universitaet Hannover
18!--------------------------------------------------------------------------------!
19!
20! Current revisions:
21! -----------------
22!
23!
24! Former revisions:
25! -----------------
26! $Id: subsidence.f90 1490 2014-10-30 08:10:16Z gronemeier $
27!
28! 1489 2014-10-30 08:09:12Z raasch
29! bugfix: sums_ls_l can only be used if large_scale_forcing is switched on
30!
31! 1382 2014-04-30 12:15:41Z boeske
32! Changed the weighting factor that is used in the summation of subsidence
33! tendencies for profile data output from weight_pres to weight_substep
34! added Neumann boundary conditions for profile data output of subsidence terms
35! at nzt+1
36!
37! 1380 2014-04-28 12:40:45Z heinze
38! Shifting only necessary in case of scalar Rayleigh damping
39!
40! 1365 2014-04-22 15:03:56Z boeske
41! Summation of subsidence tendencies for data output added
42! +ls_index, sums_ls_l, tmp_tend
43!
44! 1353 2014-04-08 15:21:23Z heinze
45! REAL constants provided with KIND-attribute
46!
47! 1320 2014-03-20 08:40:49Z raasch
48! ONLY-attribute added to USE-statements,
49! kind-parameters added to all INTEGER and REAL declaration statements,
50! kinds are defined in new module kinds,
51! old module precision_kind is removed,
52! revision history before 2012 removed,
53! comment fields (!:) to be used for variable explanations added to
54! all variable declaration statements
55!
56! 1036 2012-10-22 13:43:42Z raasch
57! code put under GPL (PALM 3.9)
58!
59! Revision 3.7 2009-12-11 14:15:58Z heinze
60! Initial revision
61!
62! Description:
63! ------------
64! Impact of large-scale subsidence or ascent as tendency term for use
65! in the prognostic equation of potential temperature. This enables the
66! construction of a constant boundary layer height z_i with time.
67!-----------------------------------------------------------------------------!
68
69
70    IMPLICIT NONE
71
72    PRIVATE
73    PUBLIC  init_w_subsidence, subsidence
74
75    INTERFACE init_w_subsidence
76       MODULE PROCEDURE init_w_subsidence
77    END INTERFACE init_w_subsidence
78
79    INTERFACE subsidence
80       MODULE PROCEDURE subsidence
81       MODULE PROCEDURE subsidence_ij
82    END INTERFACE subsidence
83
84 CONTAINS
85
86    SUBROUTINE init_w_subsidence 
87
88       USE arrays_3d,                                                          &
89           ONLY:  dzu, w_subs, zu
90
91       USE control_parameters,                                                 &
92           ONLY:  message_string, ocean, subs_vertical_gradient,               &
93                  subs_vertical_gradient_level, subs_vertical_gradient_level_i
94
95       USE indices,                                                            &
96           ONLY:  nzb, nzt
97
98       USE kinds
99
100       IMPLICIT NONE
101
102       INTEGER(iwp) ::  i !:
103       INTEGER(iwp) ::  k !:
104
105       REAL(wp)     ::  gradient   !:
106       REAL(wp)     ::  ws_surface !:
107
108       IF ( .NOT. ALLOCATED( w_subs ))  THEN
109          ALLOCATE( w_subs(nzb:nzt+1) )
110          w_subs = 0.0_wp
111       ENDIF
112
113       IF ( ocean )  THEN
114          message_string = 'Applying large scale vertical motion is not ' // &
115                           'allowed for ocean runs'
116          CALL message( 'init_w_subsidence', 'PA0324', 2, 2, 0, 6, 0 )
117       ENDIF
118
119!
120!--   Compute the profile of the subsidence/ascent velocity
121!--   using the given gradients
122      i = 1
123      gradient = 0.0_wp
124      ws_surface = 0.0_wp
125     
126
127      subs_vertical_gradient_level_i(1) = 0
128      DO  k = 1, nzt+1
129         IF ( i < 11 )  THEN
130            IF ( subs_vertical_gradient_level(i) < zu(k)  .AND. &
131                 subs_vertical_gradient_level(i) >= 0.0_wp )  THEN
132               gradient = subs_vertical_gradient(i) / 100.0_wp
133               subs_vertical_gradient_level_i(i) = k - 1
134               i = i + 1
135            ENDIF
136         ENDIF
137         IF ( gradient /= 0.0_wp )  THEN
138            IF ( k /= 1 )  THEN
139               w_subs(k) = w_subs(k-1) + dzu(k) * gradient
140            ELSE
141               w_subs(k) = ws_surface   + 0.5_wp * dzu(k) * gradient
142            ENDIF
143         ELSE
144            w_subs(k) = w_subs(k-1)
145         ENDIF
146      ENDDO
147
148!
149!--   In case of no given gradients for the subsidence/ascent velocity,
150!--   choose zero gradient
151      IF ( subs_vertical_gradient_level(1) == -9999999.9_wp )  THEN
152         subs_vertical_gradient_level(1) = 0.0_wp
153      ENDIF
154
155    END SUBROUTINE init_w_subsidence
156
157
158    SUBROUTINE subsidence( tendency, var, var_init, ls_index ) 
159
160       USE arrays_3d,                                                          &
161           ONLY:  ddzu, w_subs
162
163       USE control_parameters,                                                 &
164           ONLY:  dt_3d, intermediate_timestep_count, large_scale_forcing,     &
165                  scalar_rayleigh_damping
166
167       USE indices,                                                            &
168           ONLY:  nxl, nxlg, nxr, nxrg, nyn, nyng, nys, nysg, nzb, nzb_s_inner,&
169                  nzt
170
171       USE kinds
172
173       USE statistics,                                                         &
174           ONLY:  sums_ls_l, weight_substep
175
176       IMPLICIT NONE
177 
178       INTEGER(iwp) ::  i !:
179       INTEGER(iwp) ::  j !:
180       INTEGER(iwp) ::  k !:
181       INTEGER(iwp) ::  ls_index !:
182
183       REAL(wp)     ::  tmp_tend !:
184       REAL(wp)     ::  tmp_grad !:
185   
186       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  var      !:
187       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  tendency !:
188       REAL(wp), DIMENSION(nzb:nzt+1) ::  var_init                     !:
189       REAL(wp), DIMENSION(nzb:nzt+1) ::  var_mod                      !:
190
191       var_mod = var_init
192
193!
194!--    Influence of w_subsidence on the current tendency term
195       DO  i = nxl, nxr
196          DO  j = nys, nyn
197
198             DO  k = nzb_s_inner(j,i)+1, nzt 
199                IF ( w_subs(k) < 0.0_wp )  THEN    ! large-scale subsidence
200                   tmp_tend = - w_subs(k) *                                    &
201                              ( var(k+1,j,i) - var(k,j,i) ) * ddzu(k+1)
202                ELSE                               ! large-scale ascent
203                   tmp_tend = - w_subs(k) *                                    &
204                              ( var(k,j,i) - var(k-1,j,i) ) * ddzu(k)
205                ENDIF
206
207                tendency(k,j,i) = tendency(k,j,i) + tmp_tend
208
209                IF ( large_scale_forcing )  THEN
210                   sums_ls_l(k,ls_index) = sums_ls_l(k,ls_index) + tmp_tend    &
211                                   * weight_substep(intermediate_timestep_count)
212                ENDIF
213             ENDDO
214
215             sums_ls_l(nzt+1,ls_index) = sums_ls_l(nzt,ls_index)
216
217          ENDDO
218       ENDDO
219
220!
221!--    Shifting of the initial profile is especially necessary with Rayleigh
222!--    damping switched on
223       IF ( scalar_rayleigh_damping ) THEN
224          DO  k = nzb, nzt
225             IF ( w_subs(k) < 0.0_wp )  THEN      ! large-scale subsidence
226                var_mod(k) = var_init(k) - dt_3d * w_subs(k) *  &
227                                  ( var_init(k+1) - var_init(k) ) * ddzu(k+1)
228             ENDIF
229          ENDDO
230!
231!--      At the upper boundary, the initial profile is shifted with aid of
232!--      the gradient tmp_grad. (This is ok if the gradients are linear.)
233         IF ( w_subs(nzt) < 0.0_wp )  THEN
234            tmp_grad = ( var_init(nzt+1) - var_init(nzt) ) * ddzu(nzt+1)
235            var_mod(nzt+1) = var_init(nzt+1) -  &
236                                 dt_3d * w_subs(nzt+1) * tmp_grad
237         ENDIF
238       
239
240         DO  k = nzt+1, nzb+1, -1
241            IF ( w_subs(k) >= 0.0_wp )  THEN  ! large-scale ascent
242               var_mod(k) = var_init(k) - dt_3d * w_subs(k) *  &
243                                  ( var_init(k) - var_init(k-1) ) * ddzu(k) 
244            ENDIF
245         ENDDO
246!
247!--      At the lower boundary shifting is not necessary because the
248!--      subsidence velocity w_subs(nzb) vanishes.
249         IF ( w_subs(nzb+1) >= 0.0_wp )  THEN
250            var_mod(nzb) = var_init(nzb)
251         ENDIF
252
253         var_init = var_mod
254      ENDIF
255
256
257 END SUBROUTINE subsidence
258
259 SUBROUTINE subsidence_ij( i, j, tendency, var, var_init, ls_index ) 
260
261       USE arrays_3d,                                                          &
262           ONLY:  ddzu, w_subs
263
264       USE control_parameters,                                                 &
265           ONLY:  dt_3d, intermediate_timestep_count, large_scale_forcing,     &
266                  scalar_rayleigh_damping
267
268       USE indices,                                                            &
269           ONLY:  nxl, nxlg, nxrg, nyng, nys, nysg, nzb_s_inner, nzb, nzt
270
271       USE kinds
272
273       USE statistics,                                                         &
274           ONLY:  sums_ls_l, weight_substep
275
276       IMPLICIT NONE
277 
278       INTEGER(iwp) ::  i !:
279       INTEGER(iwp) ::  j !:
280       INTEGER(iwp) ::  k !:
281       INTEGER(iwp) ::  ls_index !:
282
283       REAL(wp)     ::  tmp_tend !:
284       REAL(wp)     ::  tmp_grad !:
285   
286       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  var      !:
287       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  tendency !:
288       REAL(wp), DIMENSION(nzb:nzt+1) ::  var_init                     !:
289       REAL(wp), DIMENSION(nzb:nzt+1) ::  var_mod                      !:
290
291       var_mod = var_init
292
293!
294!--    Influence of w_subsidence on the current tendency term
295       DO  k = nzb_s_inner(j,i)+1, nzt 
296          IF ( w_subs(k) < 0.0_wp )  THEN      ! large-scale subsidence
297             tmp_tend = - w_subs(k) * ( var(k+1,j,i) - var(k,j,i) ) * ddzu(k+1)
298          ELSE                                 ! large-scale ascent
299             tmp_tend = - w_subs(k) * ( var(k,j,i) - var(k-1,j,i) ) * ddzu(k)
300          ENDIF
301
302          tendency(k,j,i) = tendency(k,j,i) + tmp_tend
303
304          IF ( large_scale_forcing )  THEN
305             sums_ls_l(k,ls_index) = sums_ls_l(k,ls_index) + tmp_tend          &
306                                   * weight_substep(intermediate_timestep_count)
307          ENDIF
308       ENDDO
309
310       IF ( large_scale_forcing )  THEN
311          sums_ls_l(nzt+1,ls_index) = sums_ls_l(nzt,ls_index)
312       ENDIF
313
314!
315!--    Shifting of the initial profile is especially necessary with Rayleigh
316!--    damping switched on
317       IF ( scalar_rayleigh_damping ) THEN
318          IF ( i == nxl .AND. j == nys )  THEN ! shifting only once per PE
319
320             DO  k = nzb, nzt
321                IF ( w_subs(k) < 0.0_wp )  THEN      ! large-scale subsidence
322                   var_mod(k) = var_init(k) - dt_3d * w_subs(k) *  &
323                                     ( var_init(k+1) - var_init(k) ) * ddzu(k+1)
324                ENDIF
325             ENDDO
326!
327!--          At the upper boundary, the initial profile is shifted with aid of
328!--          the gradient tmp_grad. (This is ok if the gradients are linear.)
329             IF ( w_subs(nzt) < 0.0_wp )  THEN
330                tmp_grad = ( var_init(nzt+1) - var_init(nzt) ) * ddzu(nzt+1)
331                var_mod(nzt+1) = var_init(nzt+1) -  &
332                                     dt_3d * w_subs(nzt+1) * tmp_grad
333             ENDIF
334       
335
336             DO  k = nzt+1, nzb+1, -1
337                IF ( w_subs(k) >= 0.0_wp )  THEN  ! large-scale ascent
338                   var_mod(k) = var_init(k) - dt_3d * w_subs(k) *  &
339                                      ( var_init(k) - var_init(k-1) ) * ddzu(k)
340                ENDIF
341             ENDDO
342!
343!--          At the lower boundary shifting is not necessary because the
344!--          subsidence velocity w_subs(nzb) vanishes.
345             IF ( w_subs(nzb+1) >= 0.0_wp )  THEN
346                var_mod(nzb) = var_init(nzb)
347             ENDIF
348
349             var_init = var_mod 
350
351          ENDIF
352       ENDIF
353
354 END SUBROUTINE subsidence_ij
355
356
357 END MODULE subsidence_mod
Note: See TracBrowser for help on using the repository browser.