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

Last change on this file since 1036 was 1036, checked in by raasch, 9 years ago

code has been put under the GNU General Public License (v3)

  • Property svn:keywords set to Id
File size: 8.5 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-2012  Leibniz University Hannover
18!--------------------------------------------------------------------------------!
19!
20! Current revisions:
21! -----------------
22!
23! Former revisions:
24! -----------------
25! $Id: subsidence.f90 1036 2012-10-22 13:43:42Z raasch $
26!
27! 671 2011-01-11 12:04:00Z heinze $
28! bugfix: access to ddzu(nzt+2) which is not defined
29!
30! 667 2010-12-23 12:06:00Z suehring
31! nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng
32!
33! 580 2010-10-05 13:59:11Z heinze
34! Renaming of ws_vertical_gradient to subs_vertical_gradient,
35! ws_vertical_gradient_level to subs_vertical_gradient_level and
36! ws_vertical_gradient_level_ind to subs_vertical_gradient_level_i
37!
38! Revision 3.7 2009-12-11 14:15:58Z heinze
39! Initial revision
40!
41! Description:
42! ------------
43! Impact of large-scale subsidence or ascent as tendency term for use
44! in the prognostic equation of potential temperature. This enables the
45! construction of a constant boundary layer height z_i with time.
46!-----------------------------------------------------------------------------!
47
48
49    IMPLICIT NONE
50
51    PRIVATE
52    PUBLIC  init_w_subsidence, subsidence
53
54    INTERFACE init_w_subsidence
55       MODULE PROCEDURE init_w_subsidence
56    END INTERFACE init_w_subsidence
57
58    INTERFACE subsidence
59       MODULE PROCEDURE subsidence
60       MODULE PROCEDURE subsidence_ij
61    END INTERFACE subsidence
62
63 CONTAINS
64
65    SUBROUTINE init_w_subsidence
66
67       USE arrays_3d
68       USE control_parameters
69       USE grid_variables
70       USE indices
71       USE pegrid
72       USE statistics
73
74       IMPLICIT NONE
75
76       INTEGER :: i, k
77       REAL    :: gradient, ws_surface
78
79       IF ( .NOT. ALLOCATED( w_subs )) THEN
80          ALLOCATE( w_subs(nzb:nzt+1) )
81          w_subs = 0.0
82       ENDIF
83
84      IF ( ocean )  THEN
85          message_string = 'Applying large scale vertical motion is not ' // &
86                           'allowed for ocean runs'
87          CALL message( 'init_w_subsidence', 'PA0324', 2, 2, 0, 6, 0 )
88       ENDIF
89
90!
91!--   Compute the profile of the subsidence/ascent velocity
92!--   using the given gradients
93      i = 1
94      gradient = 0.0
95      ws_surface = 0.0
96     
97
98      subs_vertical_gradient_level_i(1) = 0
99      DO  k = 1, nzt+1
100         IF ( i < 11 ) THEN
101            IF ( subs_vertical_gradient_level(i) < zu(k)  .AND. &
102                 subs_vertical_gradient_level(i) >= 0.0 )  THEN
103               gradient = subs_vertical_gradient(i) / 100.0
104               subs_vertical_gradient_level_i(i) = k - 1
105               i = i + 1
106            ENDIF
107         ENDIF
108         IF ( gradient /= 0.0 )  THEN
109            IF ( k /= 1 )  THEN
110               w_subs(k) = w_subs(k-1) + dzu(k) * gradient
111            ELSE
112               w_subs(k) = ws_surface   + 0.5 * dzu(k) * gradient
113            ENDIF
114         ELSE
115            w_subs(k) = w_subs(k-1)
116         ENDIF
117      ENDDO
118
119!
120!--   In case of no given gradients for the subsidence/ascent velocity,
121!--   choose zero gradient
122      IF ( subs_vertical_gradient_level(1) == -9999999.9 )  THEN
123         subs_vertical_gradient_level(1) = 0.0
124      ENDIF
125
126    END SUBROUTINE init_w_subsidence
127
128
129    SUBROUTINE subsidence( tendency, var, var_init ) 
130
131       USE arrays_3d
132       USE control_parameters
133       USE grid_variables
134       USE indices
135       USE pegrid
136       USE statistics
137
138       IMPLICIT NONE
139 
140       INTEGER :: i, j, k
141
142       REAL :: tmp_grad
143   
144       REAL, DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) :: var, tendency
145       REAL, DIMENSION(nzb:nzt+1) :: var_init, var_mod
146
147       var_mod = var_init
148
149!
150!--    Influence of w_subsidence on the current tendency term
151       DO  i = nxl, nxr
152          DO  j = nys, nyn
153             DO  k = nzb_s_inner(j,i)+1, nzt
154                IF ( w_subs(k) < 0.0 ) THEN    ! large-scale subsidence
155                   tendency(k,j,i) = tendency(k,j,i) - w_subs(k) *  &
156                                     ( var(k+1,j,i) - var(k,j,i) ) * ddzu(k+1)
157                ELSE                           ! large-scale ascent
158                   tendency(k,j,i) = tendency(k,j,i) - w_subs(k) *  &
159                                     ( var(k,j,i) - var(k-1,j,i) ) * ddzu(k)
160                ENDIF
161             ENDDO
162          ENDDO
163       ENDDO
164
165!
166!--    Shifting of the initial profile is especially necessary with Rayleigh
167!--    damping switched on
168
169       DO  k = nzb, nzt
170          IF ( w_subs(k) < 0.0 ) THEN      ! large-scale subsidence
171             var_mod(k) = var_init(k) - dt_3d * w_subs(k) *  &
172                               ( var_init(k+1) - var_init(k) ) * ddzu(k+1)
173          ENDIF
174       ENDDO
175!
176!--   At the upper boundary, the initial profile is shifted with aid of
177!--   the gradient tmp_grad. (This is ok if the gradients are linear.)
178      IF ( w_subs(nzt) < 0.0 ) THEN
179         tmp_grad = ( var_init(nzt+1) - var_init(nzt) ) * ddzu(nzt+1)
180         var_mod(nzt+1) = var_init(nzt+1) -  &
181                              dt_3d * w_subs(nzt+1) * tmp_grad
182      ENDIF
183       
184
185      DO  k = nzt+1, nzb+1, -1
186         IF ( w_subs(k) >= 0.0 ) THEN  ! large-scale ascent
187            var_mod(k) = var_init(k) - dt_3d * w_subs(k) *  &
188                               ( var_init(k) - var_init(k-1) ) * ddzu(k) 
189         ENDIF
190      ENDDO
191!
192!--   At the lower boundary shifting is not necessary because the
193!--   subsidence velocity w_subs(nzb) vanishes.
194
195
196      IF ( w_subs(nzb+1) >= 0.0 ) THEN
197         var_mod(nzb) = var_init(nzb)
198      ENDIF
199
200      var_init = var_mod
201
202
203 END SUBROUTINE subsidence
204
205 SUBROUTINE subsidence_ij( i, j, tendency, var, var_init ) 
206
207       USE arrays_3d
208       USE control_parameters
209       USE grid_variables
210       USE indices
211       USE pegrid
212       USE statistics
213
214       IMPLICIT NONE
215 
216       INTEGER :: i, j, k
217
218       REAL :: tmp_grad
219   
220       REAL, DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) :: var, tendency
221       REAL, DIMENSION(nzb:nzt+1) :: var_init, var_mod
222
223       var_mod = var_init
224
225!
226!--    Influence of w_subsidence on the current tendency term
227       DO  k = nzb_s_inner(j,i)+1, nzt
228          IF ( w_subs(k) < 0.0 ) THEN      ! large-scale subsidence
229             tendency(k,j,i) = tendency(k,j,i) - w_subs(k) *  &
230                               ( var(k+1,j,i) - var(k,j,i) ) * ddzu(k+1)
231          ELSE                             ! large-scale ascent
232             tendency(k,j,i) = tendency(k,j,i) - w_subs(k) *  &
233                               ( var(k,j,i) - var(k-1,j,i) ) * ddzu(k)
234          ENDIF
235       ENDDO
236
237
238!
239!--    Shifting of the initial profile is especially necessary with Rayleigh
240!--    damping switched on
241       IF ( i == nxl .AND. j == nys ) THEN ! shifting only once per PE
242
243          DO  k = nzb, nzt
244             IF ( w_subs(k) < 0.0 ) THEN      ! large-scale subsidence
245                var_mod(k) = var_init(k) - dt_3d * w_subs(k) *  &
246                                  ( var_init(k+1) - var_init(k) ) * ddzu(k+1)
247             ENDIF
248          ENDDO
249!
250!--       At the upper boundary, the initial profile is shifted with aid of
251!--       the gradient tmp_grad. (This is ok if the gradients are linear.)
252          IF ( w_subs(nzt) < 0.0 ) THEN
253             tmp_grad = ( var_init(nzt+1) - var_init(nzt) ) * ddzu(nzt+1)
254             var_mod(nzt+1) = var_init(nzt+1) -  &
255                                  dt_3d * w_subs(nzt+1) * tmp_grad
256          ENDIF
257       
258
259          DO  k = nzt+1, nzb+1, -1
260             IF ( w_subs(k) >= 0.0 ) THEN  ! large-scale ascent
261                var_mod(k) = var_init(k) - dt_3d * w_subs(k) *  &
262                                   ( var_init(k) - var_init(k-1) ) * ddzu(k)
263             ENDIF
264          ENDDO
265!
266!--       At the lower boundary shifting is not necessary because the
267!--       subsidence velocity w_subs(nzb) vanishes.
268
269
270          IF ( w_subs(nzb+1) >= 0.0 ) THEN
271             var_mod(nzb) = var_init(nzb)
272          ENDIF
273
274          var_init = var_mod
275
276       ENDIF
277
278 END SUBROUTINE subsidence_ij
279
280
281 END MODULE subsidence_mod
Note: See TracBrowser for help on using the repository browser.