source: palm/trunk/SOURCE/inflow_turbulence.f90 @ 821

Last change on this file since 821 was 710, checked in by raasch, 14 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 5.3 KB
Line 
1 SUBROUTINE inflow_turbulence
2
3!------------------------------------------------------------------------------!
4! Current revisions:
5! -----------------
6!
7!
8! Former revisions:
9! -----------------
10! $Id: inflow_turbulence.f90 710 2011-03-30 09:45:27Z maronga $
11!
12! 709 2011-03-30 09:31:40Z raasch
13! formatting adjustments
14!
15! 667 2010-12-23 12:06:00Z suehring/gryschka
16! Using nbgp recycling planes for a better resolution of the turbulent flow
17! near the inflow.
18!
19! 622 2010-12-10 08:08:13Z raasch
20! optional barriers included in order to speed up collective operations
21!
22! 222 2009-01-12 16:04:16Z letzel
23! Bugfix for nonparallel execution
24!
25! Initial version (2008/03/07)
26!
27! Description:
28! ------------
29! Imposing turbulence at the respective inflow using the turbulence
30! recycling method of Kataoka and Mizuno (2002).
31!------------------------------------------------------------------------------!
32
33    USE arrays_3d
34    USE control_parameters
35    USE cpulog
36    USE grid_variables
37    USE indices
38    USE interfaces
39    USE pegrid
40
41
42    IMPLICIT NONE
43
44    INTEGER ::  i, imax, j, k, l, ngp_ifd, ngp_pr
45
46    REAL, DIMENSION(1:2) ::  volume_flow_l, volume_flow_offset
47    REAL, DIMENSION(nzb:nzt+1,5,nbgp) ::  avpr, avpr_l
48    REAL, DIMENSION(nzb:nzt+1,nysg:nyng,5,nbgp) ::  inflow_dist
49
50    CALL cpu_log( log_point(40), 'inflow_turbulence', 'start' )
51
52!
53!-- Carry out spanwise averaging in the recycling plane
54    avpr_l = 0.0
55    ngp_pr = ( nzt - nzb + 2 ) * 5 * nbgp
56    ngp_ifd = ngp_pr * ( nyn - nys + 1 + 2 * nbgp )
57
58!
59!-- First, local averaging within the recycling domain
60    i = recycling_plane
61
62#if defined( __parallel )
63    IF ( myidx == id_recycling )  THEN
64       
65       DO  l = 1, nbgp
66          DO  j = nys, nyn
67             DO  k = nzb, nzt + 1
68
69                avpr_l(k,1,l) = avpr_l(k,1,l) + u(k,j,i)
70                avpr_l(k,2,l) = avpr_l(k,2,l) + v(k,j,i)
71                avpr_l(k,3,l) = avpr_l(k,3,l) + w(k,j,i)
72                avpr_l(k,4,l) = avpr_l(k,4,l) + pt(k,j,i)
73                avpr_l(k,5,l) = avpr_l(k,5,l) + e(k,j,i)
74
75             ENDDO
76          ENDDO
77          i = i + 1
78       ENDDO
79
80    ENDIF
81!
82!-- Now, averaging over all PEs
83    IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
84    CALL MPI_ALLREDUCE( avpr_l(nzb,1,1), avpr(nzb,1,1), ngp_pr, MPI_REAL, &
85                        MPI_SUM, comm2d, ierr )
86
87#else
88    DO  l = 1, nbgp
89       DO  j = nys, nyn
90          DO  k = nzb, nzt + 1
91
92             avpr_l(k,1,l) = avpr_l(k,1,l) + u(k,j,i)
93             avpr_l(k,2,l) = avpr_l(k,2,l) + v(k,j,i)
94             avpr_l(k,3,l) = avpr_l(k,3,l) + w(k,j,i)
95             avpr_l(k,4,l) = avpr_l(k,4,l) + pt(k,j,i)
96             avpr_l(k,5,l) = avpr_l(k,5,l) + e(k,j,i)
97
98          ENDDO
99       ENDDO
100       i = i + 1 
101    ENDDO
102   
103    avpr = avpr_l
104#endif
105
106    avpr = avpr / ( ny + 1 )
107!
108!-- Calculate the disturbances at the recycling plane
109    i = recycling_plane
110
111#if defined( __parallel )
112    IF ( myidx == id_recycling )  THEN
113       DO  l = 1, nbgp
114          DO  j = nysg, nyng
115             DO  k = nzb, nzt + 1
116
117                inflow_dist(k,j,1,l) = u(k,j,i+1) - avpr(k,1,l)
118                inflow_dist(k,j,2,l) = v(k,j,i)   - avpr(k,2,l)
119                inflow_dist(k,j,3,l) = w(k,j,i)   - avpr(k,3,l)
120                inflow_dist(k,j,4,l) = pt(k,j,i)  - avpr(k,4,l)
121                inflow_dist(k,j,5,l) = e(k,j,i)   - avpr(k,5,l)
122             
123            ENDDO
124          ENDDO
125          i = i + 1
126       ENDDO
127
128    ENDIF
129#else
130    DO  l = 1, nbgp
131       DO  j = nysg, nyng
132          DO  k = nzb, nzt+1
133
134             inflow_dist(k,j,1,l) = u(k,j,i+1) - avpr(k,1,l)
135             inflow_dist(k,j,2,l) = v(k,j,i)   - avpr(k,2,l)
136             inflow_dist(k,j,3,l) = w(k,j,i)   - avpr(k,3,l)
137             inflow_dist(k,j,4,l) = pt(k,j,i)  - avpr(k,4,l)
138             inflow_dist(k,j,5,l) = e(k,j,i)   - avpr(k,5,l)
139             
140          ENDDO
141       ENDDO
142       i = i + 1
143    ENDDO
144#endif
145
146!
147!-- For parallel runs, send the disturbances to the respective inflow PE
148#if defined( __parallel )
149    IF ( myidx == id_recycling  .AND.  myidx /= id_inflow )  THEN
150
151       CALL MPI_SEND( inflow_dist(nzb,nysg,1,1), ngp_ifd, MPI_REAL, &
152                      id_inflow, 1, comm1dx, ierr )
153
154    ELSEIF ( myidx /= id_recycling  .AND.  myidx == id_inflow )  THEN
155
156       inflow_dist = 0.0
157       CALL MPI_RECV( inflow_dist(nzb,nysg,1,1), ngp_ifd, MPI_REAL, &
158                      id_recycling, 1, comm1dx, status, ierr )
159
160    ENDIF
161#endif
162
163!
164!-- Add the disturbance at the inflow
165    IF ( nxl == 0 )  THEN
166
167       DO  j = nysg, nyng
168          DO  k = nzb, nzt + 1
169
170              u(k,j,-nbgp+1:0) = mean_inflow_profiles(k,1) + &
171                           inflow_dist(k,j,1,1:nbgp) * inflow_damping_factor(k)
172              v(k,j,-nbgp:-1)  = mean_inflow_profiles(k,2) + &
173                           inflow_dist(k,j,2,1:nbgp) * inflow_damping_factor(k)
174              w(k,j,-nbgp:-1)  =                             &
175                           inflow_dist(k,j,3,1:nbgp) * inflow_damping_factor(k)
176              pt(k,j,-nbgp:-1) = mean_inflow_profiles(k,4) + &
177                           inflow_dist(k,j,4,1:nbgp) * inflow_damping_factor(k)
178              e(k,j,-nbgp:-1)  = mean_inflow_profiles(k,5) + &
179                           inflow_dist(k,j,5,1:nbgp) * inflow_damping_factor(k)
180              e(k,j,-nbgp:-1)  = MAX( e(k,j,-nbgp:-1), 0.0 )
181
182          ENDDO
183       ENDDO
184
185    ENDIF
186
187    CALL cpu_log( log_point(40), 'inflow_turbulence', 'stop' )
188
189
190 END SUBROUTINE inflow_turbulence
Note: See TracBrowser for help on using the repository browser.