source: palm/trunk/SOURCE/lpm_init_sgs_tke.f90 @ 1310

Last change on this file since 1310 was 1310, checked in by raasch, 7 years ago

update of GPL copyright

  • Property svn:keywords set to Id
File size: 7.9 KB
Line 
1 SUBROUTINE lpm_init_sgs_tke
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: lpm_init_sgs_tke.f90 1310 2014-03-14 08:01:56Z raasch $
27!
28! 1036 2012-10-22 13:43:42Z raasch
29! code put under GPL (PALM 3.9)
30!
31! 849 2012-03-15 10:35:09Z raasch
32! initial revision (former part of advec_particles)
33!
34!
35! Description:
36! ------------
37! Calculates quantities required for considering the SGS velocity fluctuations
38! in the particle transport by a stochastic approach. The respective
39! quantities are: SGS-TKE gradients and horizontally averaged profiles of the
40! SGS TKE and the resolved-scale velocity variances.
41!------------------------------------------------------------------------------!
42
43    USE arrays_3d
44    USE control_parameters
45    USE grid_variables
46    USE indices
47    USE particle_attributes
48    USE pegrid
49    USE statistics
50
51    IMPLICIT NONE
52
53    INTEGER ::  i, j, k
54
55
56!
57!-- TKE gradient along x and y
58    DO  i = nxl, nxr
59       DO  j = nys, nyn
60          DO  k = nzb, nzt+1
61
62             IF ( k <= nzb_s_inner(j,i-1)  .AND.  k > nzb_s_inner(j,i)  .AND.  &
63                  k  > nzb_s_inner(j,i+1) )                                    &
64             THEN
65                de_dx(k,j,i) = 2.0 * sgs_wfu_part * ( e(k,j,i+1) - e(k,j,i) )  &
66                               * ddx
67             ELSEIF ( k  > nzb_s_inner(j,i-1)  .AND.  k > nzb_s_inner(j,i)     &
68                      .AND.  k <= nzb_s_inner(j,i+1) )                         &
69             THEN
70                de_dx(k,j,i) = 2.0 * sgs_wfu_part * ( e(k,j,i) - e(k,j,i-1) )  &
71                               * ddx
72             ELSEIF ( k < nzb_s_inner(j,i)  .AND.  k < nzb_s_inner(j,i+1) )    &
73             THEN
74                de_dx(k,j,i) = 0.0
75             ELSEIF ( k < nzb_s_inner(j,i-1)  .AND.  k < nzb_s_inner(j,i) )    &
76             THEN
77                de_dx(k,j,i) = 0.0
78             ELSE
79                de_dx(k,j,i) = sgs_wfu_part * ( e(k,j,i+1) - e(k,j,i-1) ) * ddx
80             ENDIF
81
82             IF ( k <= nzb_s_inner(j-1,i)  .AND.  k > nzb_s_inner(j,i)  .AND.  &
83                  k  > nzb_s_inner(j+1,i) )                                    &
84             THEN
85                de_dy(k,j,i) = 2.0 * sgs_wfv_part * ( e(k,j+1,i) - e(k,j,i) )  &
86                               * ddy
87             ELSEIF ( k  > nzb_s_inner(j-1,i)  .AND.  k  > nzb_s_inner(j,i)    &
88                      .AND.  k <= nzb_s_inner(j+1,i) )                         &
89             THEN
90                de_dy(k,j,i) = 2.0 * sgs_wfv_part * ( e(k,j,i) - e(k,j-1,i) )  &
91                               * ddy
92             ELSEIF ( k < nzb_s_inner(j,i)  .AND.  k < nzb_s_inner(j+1,i) )    &
93             THEN
94                de_dy(k,j,i) = 0.0
95             ELSEIF ( k < nzb_s_inner(j-1,i)  .AND.  k < nzb_s_inner(j,i) )    &
96             THEN
97                de_dy(k,j,i) = 0.0
98             ELSE
99                de_dy(k,j,i) = sgs_wfv_part * ( e(k,j+1,i) - e(k,j-1,i) ) * ddy
100             ENDIF
101
102          ENDDO
103       ENDDO
104    ENDDO
105
106!
107!-- TKE gradient along z, including bottom and top boundary conditions
108    DO  i = nxl, nxr
109       DO  j = nys, nyn
110
111          DO  k = nzb_s_inner(j,i)+2, nzt-1
112             de_dz(k,j,i)  = 2.0 * sgs_wfw_part * &
113                             ( e(k+1,j,i) - e(k-1,j,i) ) / ( zu(k+1)-zu(k-1) )
114          ENDDO
115
116          k = nzb_s_inner(j,i)
117          de_dz(nzb:k,j,i) = 0.0
118          de_dz(k+1,j,i)   = 2.0 * sgs_wfw_part * ( e(k+2,j,i) - e(k+1,j,i) ) &
119                                                / ( zu(k+2) - zu(k+1) )
120          de_dz(nzt,j,i)   = 0.0
121          de_dz(nzt+1,j,i) = 0.0
122       ENDDO
123    ENDDO
124
125
126!
127!-- Lateral boundary conditions
128    CALL exchange_horiz( de_dx, nbgp )
129    CALL exchange_horiz( de_dy, nbgp )
130    CALL exchange_horiz( de_dz, nbgp )
131    CALL exchange_horiz( diss, nbgp  )
132
133
134!
135!-- Calculate the horizontally averaged profiles of SGS TKE and resolved
136!-- velocity variances (they may have been already calculated in routine
137!-- flow_statistics).
138    IF ( .NOT. flow_statistics_called )  THEN
139
140!
141!--    First calculate horizontally averaged profiles of the horizontal
142!--    velocities.
143       sums_l(:,1,0) = 0.0
144       sums_l(:,2,0) = 0.0
145
146       DO  i = nxl, nxr
147          DO  j =  nys, nyn
148             DO  k = nzb_s_outer(j,i), nzt+1
149                sums_l(k,1,0)  = sums_l(k,1,0)  + u(k,j,i)
150                sums_l(k,2,0)  = sums_l(k,2,0)  + v(k,j,i)
151             ENDDO
152          ENDDO
153       ENDDO
154
155#if defined( __parallel )
156!
157!--    Compute total sum from local sums
158       IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
159       CALL MPI_ALLREDUCE( sums_l(nzb,1,0), sums(nzb,1), nzt+2-nzb, &
160                           MPI_REAL, MPI_SUM, comm2d, ierr )
161       IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
162       CALL MPI_ALLREDUCE( sums_l(nzb,2,0), sums(nzb,2), nzt+2-nzb, &
163                              MPI_REAL, MPI_SUM, comm2d, ierr )
164#else
165       sums(:,1) = sums_l(:,1,0)
166       sums(:,2) = sums_l(:,2,0)
167#endif
168
169!
170!--    Final values are obtained by division by the total number of grid
171!--    points used for the summation.
172       hom(:,1,1,0) = sums(:,1) / ngp_2dh_outer(:,0)   ! u
173       hom(:,1,2,0) = sums(:,2) / ngp_2dh_outer(:,0)   ! v
174
175!
176!--    Now calculate the profiles of SGS TKE and the resolved-scale
177!--    velocity variances
178       sums_l(:,8,0)  = 0.0
179       sums_l(:,30,0) = 0.0
180       sums_l(:,31,0) = 0.0
181       sums_l(:,32,0) = 0.0
182       DO  i = nxl, nxr
183          DO  j = nys, nyn
184             DO  k = nzb_s_outer(j,i), nzt+1
185                sums_l(k,8,0)  = sums_l(k,8,0)  + e(k,j,i)
186                sums_l(k,30,0) = sums_l(k,30,0) + ( u(k,j,i) - hom(k,1,1,0) )**2
187                sums_l(k,31,0) = sums_l(k,31,0) + ( v(k,j,i) - hom(k,1,2,0) )**2
188                sums_l(k,32,0) = sums_l(k,32,0) + w(k,j,i)**2
189             ENDDO
190          ENDDO
191       ENDDO
192
193#if defined( __parallel )
194!
195!--    Compute total sum from local sums
196       IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
197       CALL MPI_ALLREDUCE( sums_l(nzb,8,0), sums(nzb,8), nzt+2-nzb, &
198                           MPI_REAL, MPI_SUM, comm2d, ierr )
199       IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
200       CALL MPI_ALLREDUCE( sums_l(nzb,30,0), sums(nzb,30), nzt+2-nzb, &
201                           MPI_REAL, MPI_SUM, comm2d, ierr )
202       IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
203       CALL MPI_ALLREDUCE( sums_l(nzb,31,0), sums(nzb,31), nzt+2-nzb, &
204                           MPI_REAL, MPI_SUM, comm2d, ierr )
205       IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
206       CALL MPI_ALLREDUCE( sums_l(nzb,32,0), sums(nzb,32), nzt+2-nzb, &
207                           MPI_REAL, MPI_SUM, comm2d, ierr )
208
209#else
210       sums(:,8)  = sums_l(:,8,0)
211       sums(:,30) = sums_l(:,30,0)
212       sums(:,31) = sums_l(:,31,0)
213       sums(:,32) = sums_l(:,32,0)
214#endif
215
216!
217!--    Final values are obtained by division by the total number of grid
218!--    points used for the summation.
219       hom(:,1,8,0)  = sums(:,8)  / ngp_2dh_outer(:,0)   ! e
220       hom(:,1,30,0) = sums(:,30) / ngp_2dh_outer(:,0)   ! u*2
221       hom(:,1,31,0) = sums(:,31) / ngp_2dh_outer(:,0)   ! v*2
222       hom(:,1,32,0) = sums(:,32) / ngp_2dh_outer(:,0)   ! w*2
223
224    ENDIF
225
226 END SUBROUTINE lpm_init_sgs_tke
Note: See TracBrowser for help on using the repository browser.