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

Last change on this file since 1846 was 1818, checked in by maronga, 9 years ago

last commit documented / copyright update

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