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

Last change on this file since 3294 was 3294, checked in by raasch, 3 years ago

modularization of the ocean code

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