source: palm/trunk/SOURCE/subsidence_mod.f90 @ 2696

Last change on this file since 2696 was 2696, checked in by kanani, 4 years ago

Merge of branch palm4u into trunk

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