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

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

last commit documented

  • Property svn:keywords set to Id
File size: 5.3 KB
RevLine 
[151]1 SUBROUTINE inflow_turbulence
2
3!------------------------------------------------------------------------------!
[484]4! Current revisions:
[151]5! -----------------
[710]6!
[151]7!
8! Former revisions:
9! -----------------
10! $Id: inflow_turbulence.f90 710 2011-03-30 09:45:27Z suehring $
11!
[710]12! 709 2011-03-30 09:31:40Z raasch
13! formatting adjustments
14!
[668]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!
[623]19! 622 2010-12-10 08:08:13Z raasch
20! optional barriers included in order to speed up collective operations
21!
[226]22! 222 2009-01-12 16:04:16Z letzel
23! Bugfix for nonparallel execution
24!
[198]25! Initial version (2008/03/07)
[151]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
[667]44    INTEGER ::  i, imax, j, k, l, ngp_ifd, ngp_pr
[151]45
46    REAL, DIMENSION(1:2) ::  volume_flow_l, volume_flow_offset
[667]47    REAL, DIMENSION(nzb:nzt+1,5,nbgp) ::  avpr, avpr_l
48    REAL, DIMENSION(nzb:nzt+1,nysg:nyng,5,nbgp) ::  inflow_dist
[151]49
50    CALL cpu_log( log_point(40), 'inflow_turbulence', 'start' )
51
52!
[667]53!-- Carry out spanwise averaging in the recycling plane
[151]54    avpr_l = 0.0
[667]55    ngp_pr = ( nzt - nzb + 2 ) * 5 * nbgp
56    ngp_ifd = ngp_pr * ( nyn - nys + 1 + 2 * nbgp )
[151]57
58!
59!-- First, local averaging within the recycling domain
[667]60    i = recycling_plane
[151]61
[667]62#if defined( __parallel )
63    IF ( myidx == id_recycling )  THEN
64       
65       DO  l = 1, nbgp
[151]66          DO  j = nys, nyn
[667]67             DO  k = nzb, nzt + 1
[151]68
[667]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)
[151]74
75             ENDDO
76          ENDDO
[667]77          i = i + 1
[151]78       ENDDO
79
80    ENDIF
81!
82!-- Now, averaging over all PEs
[622]83    IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
[709]84    CALL MPI_ALLREDUCE( avpr_l(nzb,1,1), avpr(nzb,1,1), ngp_pr, MPI_REAL, &
85                        MPI_SUM, comm2d, ierr )
[667]86
[151]87#else
[667]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   
[151]103    avpr = avpr_l
104#endif
105
[667]106    avpr = avpr / ( ny + 1 )
[151]107!
108!-- Calculate the disturbances at the recycling plane
109    i = recycling_plane
110
[222]111#if defined( __parallel )
[163]112    IF ( myidx == id_recycling )  THEN
[667]113       DO  l = 1, nbgp
114          DO  j = nysg, nyng
115             DO  k = nzb, nzt + 1
[151]116
[667]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
[151]124          ENDDO
[667]125          i = i + 1
[151]126       ENDDO
127
128    ENDIF
[222]129#else
[667]130    DO  l = 1, nbgp
131       DO  j = nysg, nyng
132          DO  k = nzb, nzt+1
[151]133
[667]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
[222]141       ENDDO
[667]142       i = i + 1
[222]143    ENDDO
144#endif
145
[151]146!
147!-- For parallel runs, send the disturbances to the respective inflow PE
148#if defined( __parallel )
[163]149    IF ( myidx == id_recycling  .AND.  myidx /= id_inflow )  THEN
[151]150
[667]151       CALL MPI_SEND( inflow_dist(nzb,nysg,1,1), ngp_ifd, MPI_REAL, &
[151]152                      id_inflow, 1, comm1dx, ierr )
153
[163]154    ELSEIF ( myidx /= id_recycling  .AND.  myidx == id_inflow )  THEN
[151]155
[163]156       inflow_dist = 0.0
[667]157       CALL MPI_RECV( inflow_dist(nzb,nysg,1,1), ngp_ifd, MPI_REAL, &
[163]158                      id_recycling, 1, comm1dx, status, ierr )
[151]159
160    ENDIF
161#endif
162
163!
164!-- Add the disturbance at the inflow
165    IF ( nxl == 0 )  THEN
166
[667]167       DO  j = nysg, nyng
168          DO  k = nzb, nzt + 1
[151]169
[709]170              u(k,j,-nbgp+1:0) = mean_inflow_profiles(k,1) + &
[667]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)
[709]174              w(k,j,-nbgp:-1)  =                             &
175                           inflow_dist(k,j,3,1:nbgp) * inflow_damping_factor(k)
[667]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 )
[151]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.