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

Last change on this file since 4083 was 3655, checked in by knoop, 6 years ago

Bugfix: made "unit" and "found" intend INOUT in module interface subroutines + automatic copyright update

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