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

Last change on this file since 1359 was 1354, checked in by heinze, 11 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 10.1 KB
Line 
1 MODULE subsidence_mod
2
3!--------------------------------------------------------------------------------!
4! This file is part of PALM.
5!
6! PALM is free software: you can redistribute it and/or modify it under the terms
7! of the GNU General Public License as published by the Free Software Foundation,
8! either version 3 of the License, or (at your option) any later 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-2014 Leibniz Universitaet Hannover
18!--------------------------------------------------------------------------------!
19!
20! Current revisions:
21! -----------------
22!
23!
24! Former revisions:
25! -----------------
26! $Id: subsidence.f90 1354 2014-04-08 15:22:57Z hoffmann $
27!
28! 1353 2014-04-08 15:21:23Z heinze
29! REAL constants provided with KIND-attribute
30!
31! 1320 2014-03-20 08:40:49Z raasch
32! ONLY-attribute added to USE-statements,
33! kind-parameters added to all INTEGER and REAL declaration statements,
34! kinds are defined in new module kinds,
35! old module precision_kind is removed,
36! revision history before 2012 removed,
37! comment fields (!:) to be used for variable explanations added to
38! all variable declaration statements
39!
40! 1036 2012-10-22 13:43:42Z raasch
41! code put under GPL (PALM 3.9)
42!
43! Revision 3.7 2009-12-11 14:15:58Z heinze
44! Initial revision
45!
46! Description:
47! ------------
48! Impact of large-scale subsidence or ascent as tendency term for use
49! in the prognostic equation of potential temperature. This enables the
50! construction of a constant boundary layer height z_i with time.
51!-----------------------------------------------------------------------------!
52
53
54    IMPLICIT NONE
55
56    PRIVATE
57    PUBLIC  init_w_subsidence, subsidence
58
59    INTERFACE init_w_subsidence
60       MODULE PROCEDURE init_w_subsidence
61    END INTERFACE init_w_subsidence
62
63    INTERFACE subsidence
64       MODULE PROCEDURE subsidence
65       MODULE PROCEDURE subsidence_ij
66    END INTERFACE subsidence
67
68 CONTAINS
69
70    SUBROUTINE init_w_subsidence 
71
72       USE arrays_3d,                                                          &
73           ONLY:  dzu, w_subs, zu
74
75       USE control_parameters,                                                 &
76           ONLY:  message_string, ocean, subs_vertical_gradient,               &
77                  subs_vertical_gradient_level, subs_vertical_gradient_level_i
78
79       USE indices,                                                            &
80           ONLY:  nzb, nzt
81
82       USE kinds
83
84       IMPLICIT NONE
85
86       INTEGER(iwp) ::  i !:
87       INTEGER(iwp) ::  k !:
88
89       REAL(wp)     ::  gradient   !:
90       REAL(wp)     ::  ws_surface !:
91
92       IF ( .NOT. ALLOCATED( w_subs )) THEN
93          ALLOCATE( w_subs(nzb:nzt+1) )
94          w_subs = 0.0_wp
95       ENDIF
96
97      IF ( ocean )  THEN
98          message_string = 'Applying large scale vertical motion is not ' // &
99                           'allowed for ocean runs'
100          CALL message( 'init_w_subsidence', 'PA0324', 2, 2, 0, 6, 0 )
101       ENDIF
102
103!
104!--   Compute the profile of the subsidence/ascent velocity
105!--   using the given gradients
106      i = 1
107      gradient = 0.0_wp
108      ws_surface = 0.0_wp
109     
110
111      subs_vertical_gradient_level_i(1) = 0
112      DO  k = 1, nzt+1
113         IF ( i < 11 ) THEN
114            IF ( subs_vertical_gradient_level(i) < zu(k)  .AND. &
115                 subs_vertical_gradient_level(i) >= 0.0_wp )  THEN
116               gradient = subs_vertical_gradient(i) / 100.0_wp
117               subs_vertical_gradient_level_i(i) = k - 1
118               i = i + 1
119            ENDIF
120         ENDIF
121         IF ( gradient /= 0.0_wp )  THEN
122            IF ( k /= 1 )  THEN
123               w_subs(k) = w_subs(k-1) + dzu(k) * gradient
124            ELSE
125               w_subs(k) = ws_surface   + 0.5_wp * dzu(k) * gradient
126            ENDIF
127         ELSE
128            w_subs(k) = w_subs(k-1)
129         ENDIF
130      ENDDO
131
132!
133!--   In case of no given gradients for the subsidence/ascent velocity,
134!--   choose zero gradient
135      IF ( subs_vertical_gradient_level(1) == -9999999.9_wp )  THEN
136         subs_vertical_gradient_level(1) = 0.0_wp
137      ENDIF
138
139    END SUBROUTINE init_w_subsidence
140
141
142    SUBROUTINE subsidence( tendency, var, var_init ) 
143
144       USE arrays_3d,                                                          &
145           ONLY:  ddzu, w_subs
146
147       USE control_parameters,                                                 &
148           ONLY:  dt_3d
149
150       USE indices,                                                            &
151           ONLY:  nxl, nxlg, nxr, nxrg, nyn, nyng, nys, nysg, nzb, nzb_s_inner,&
152                  nzt
153
154       USE kinds
155
156       IMPLICIT NONE
157 
158       INTEGER(iwp) ::  i !:
159       INTEGER(iwp) ::  j !:
160       INTEGER(iwp) ::  k !:
161
162       REAL(wp)     ::  tmp_grad !:
163   
164       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  var      !:
165       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  tendency !:
166       REAL(wp), DIMENSION(nzb:nzt+1) ::  var_init                     !:
167       REAL(wp), DIMENSION(nzb:nzt+1) ::  var_mod                      !:
168
169       var_mod = var_init
170
171!
172!--    Influence of w_subsidence on the current tendency term
173       DO  i = nxl, nxr
174          DO  j = nys, nyn
175             DO  k = nzb_s_inner(j,i)+1, nzt 
176                IF ( w_subs(k) < 0.0_wp ) THEN    ! large-scale subsidence
177                   tendency(k,j,i) = tendency(k,j,i) - w_subs(k) *  &
178                                     ( var(k+1,j,i) - var(k,j,i) ) * ddzu(k+1)
179                ELSE                           ! large-scale ascent
180                   tendency(k,j,i) = tendency(k,j,i) - w_subs(k) *  &
181                                     ( var(k,j,i) - var(k-1,j,i) ) * ddzu(k)
182                ENDIF
183             ENDDO
184          ENDDO
185       ENDDO
186
187!
188!--    Shifting of the initial profile is especially necessary with Rayleigh
189!--    damping switched on
190
191       DO  k = nzb, nzt
192          IF ( w_subs(k) < 0.0_wp ) THEN      ! large-scale subsidence
193             var_mod(k) = var_init(k) - dt_3d * w_subs(k) *  &
194                               ( var_init(k+1) - var_init(k) ) * ddzu(k+1)
195          ENDIF
196       ENDDO
197!
198!--   At the upper boundary, the initial profile is shifted with aid of
199!--   the gradient tmp_grad. (This is ok if the gradients are linear.)
200      IF ( w_subs(nzt) < 0.0_wp ) THEN
201         tmp_grad = ( var_init(nzt+1) - var_init(nzt) ) * ddzu(nzt+1)
202         var_mod(nzt+1) = var_init(nzt+1) -  &
203                              dt_3d * w_subs(nzt+1) * tmp_grad
204      ENDIF
205       
206
207      DO  k = nzt+1, nzb+1, -1
208         IF ( w_subs(k) >= 0.0_wp ) THEN  ! large-scale ascent
209            var_mod(k) = var_init(k) - dt_3d * w_subs(k) *  &
210                               ( var_init(k) - var_init(k-1) ) * ddzu(k) 
211         ENDIF
212      ENDDO
213!
214!--   At the lower boundary shifting is not necessary because the
215!--   subsidence velocity w_subs(nzb) vanishes.
216
217
218      IF ( w_subs(nzb+1) >= 0.0_wp ) THEN
219         var_mod(nzb) = var_init(nzb)
220      ENDIF
221
222      var_init = var_mod
223
224
225 END SUBROUTINE subsidence
226
227 SUBROUTINE subsidence_ij( i, j, tendency, var, var_init ) 
228
229       USE arrays_3d,                                                          &
230           ONLY:  ddzu, w_subs
231
232       USE control_parameters,                                                 &
233           ONLY:  dt_3d
234
235       USE indices,                                                            &
236           ONLY:  nxl, nxlg, nxrg, nyng, nys, nysg, nzb_s_inner, nzb, nzt
237
238       USE kinds
239
240       IMPLICIT NONE
241 
242       INTEGER(iwp) ::  i !:
243       INTEGER(iwp) ::  j !:
244       INTEGER(iwp) ::  k !:
245
246       REAL(wp)     ::  tmp_grad !:
247   
248       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  var      !:
249       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  tendency !:
250       REAL(wp), DIMENSION(nzb:nzt+1) ::  var_init                     !:
251       REAL(wp), DIMENSION(nzb:nzt+1) ::  var_mod                      !:
252
253       var_mod = var_init
254
255!
256!--    Influence of w_subsidence on the current tendency term
257       DO  k = nzb_s_inner(j,i)+1, nzt 
258          IF ( w_subs(k) < 0.0_wp ) THEN      ! large-scale subsidence
259             tendency(k,j,i) = tendency(k,j,i) - w_subs(k) *  &
260                               ( var(k+1,j,i) - var(k,j,i) ) * ddzu(k+1)
261          ELSE                             ! large-scale ascent
262             tendency(k,j,i) = tendency(k,j,i) - w_subs(k) *  &
263                               ( var(k,j,i) - var(k-1,j,i) ) * ddzu(k)
264          ENDIF
265       ENDDO
266
267
268!
269!--    Shifting of the initial profile is especially necessary with Rayleigh
270!--    damping switched on
271       IF ( i == nxl .AND. j == nys ) THEN ! shifting only once per PE
272
273          DO  k = nzb, nzt
274             IF ( w_subs(k) < 0.0_wp ) THEN      ! large-scale subsidence
275                var_mod(k) = var_init(k) - dt_3d * w_subs(k) *  &
276                                  ( var_init(k+1) - var_init(k) ) * ddzu(k+1)
277             ENDIF
278          ENDDO
279!
280!--       At the upper boundary, the initial profile is shifted with aid of
281!--       the gradient tmp_grad. (This is ok if the gradients are linear.)
282          IF ( w_subs(nzt) < 0.0_wp ) THEN
283             tmp_grad = ( var_init(nzt+1) - var_init(nzt) ) * ddzu(nzt+1)
284             var_mod(nzt+1) = var_init(nzt+1) -  &
285                                  dt_3d * w_subs(nzt+1) * tmp_grad
286          ENDIF
287       
288
289          DO  k = nzt+1, nzb+1, -1
290             IF ( w_subs(k) >= 0.0_wp ) THEN  ! large-scale ascent
291                var_mod(k) = var_init(k) - dt_3d * w_subs(k) *  &
292                                   ( var_init(k) - var_init(k-1) ) * ddzu(k)
293             ENDIF
294          ENDDO
295!
296!--       At the lower boundary shifting is not necessary because the
297!--       subsidence velocity w_subs(nzb) vanishes.
298
299
300          IF ( w_subs(nzb+1) >= 0.0_wp ) THEN
301             var_mod(nzb) = var_init(nzb)
302          ENDIF
303
304          var_init = var_mod 
305
306       ENDIF
307
308 END SUBROUTINE subsidence_ij
309
310
311 END MODULE subsidence_mod
Note: See TracBrowser for help on using the repository browser.