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

Last change on this file since 856 was 850, checked in by raasch, 13 years ago

last commit documented

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