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

Last change on this file since 2232 was 2232, checked in by suehring, 4 years ago

Adjustments according new topography and surface-modelling concept implemented

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