ADflow  v1.0
ADflow is a finite volume RANS solver tailored for gradient-based aerodynamic design optimization.
flagNearWall.F90
Go to the documentation of this file.
1 subroutine flagnearwallcells(level, sps)
2 
3  ! This routine is vastly more complex that it really should
4  ! be. Essentially what we want to do is flag cells that are
5  ! approximately within a distance "nearWallDistance" of a
6  ! wall. Nominally, this just means looking at planes parallel to a
7  ! wall boundary and determining the distance to the point on the
8  ! surface. Very quick (no searching). However, if a block is split
9  ! parallel to a wall BC which is not desirable, but will happen
10  ! occasionally for block partitioning purposes, we can have "near
11  ! wall" cells that are on a split block and have no idea they are
12  ! near a wall. So what we have to do is similar to the flooding
13  ! algorithm: Starting at the wall, you proceed outwards. If you get
14  ! to the other side of the block and it is still near a wall, you
15  ! communicate and then continue on the other side. For the fringe
16  ! search, (which is built from the dual mesh), we want to know of
17  ! the *cells* of the dual mesh are near a wall. That means we
18  ! actually want check the distnaces of the dual of the
19  ! dual. Basically that means we're back to the primal.
20 
21  use constants
22  use oversetdata
23  use inputoverset
24  use blockpointers
25  use adtapi
26  use cgnsgrid
27  use communication
28  use utils
29  implicit none
30 
31  ! Input Params
32  integer(kind=intType), intent(in) :: level, sps
33 
34  ! Working paramters
35  integer(kind=intType) :: i, j, k, nn, mm, ierr, nFlaggedLocal
36  integer(kind=intType) :: iStart, iEnd, jStart, jEnd, kStart, kEnd
37  integer(kind=intType) :: loopIter, nAtBoundaryLocal, nAtBoundary
38  logical :: tostop
39  type(xplane), dimension(:), allocatable :: planes
40  real(Kind=realtype), dimension(:, :, :, :), pointer :: xx, xseed
41 
42  real(kind=realtype) :: dist
43  ! Since we need to have the 'x' and 'nearWall' arrays allocated at
44  ! once, these need to go into block. Allocate the two x arrays the
45  ! near wall arrays
46  natboundarylocal = 0
47  nflaggedlocal = 0
48  do nn = 1, ndom
49  call setpointers(nn, level, sps)
50  if (.not. associated(flowdoms(nn, level, sps)%nearWall)) then
51  allocate (flowdoms(nn, level, sps)%nearWall(1:il, 1:jl, 1:kl))
52  end if
53 
54  allocate (flowdoms(nn, level, sps)%XSeed(0:ie, 0:je, 0:ke, 3))
55 
56  ! Manaully set the three pointers
57  xseed => flowdoms(nn, level, sps)%xSeed
58  nearwall => flowdoms(nn, level, sps)%nearWall
59  nearwall = 0
60  ! Use large to denote no seed presnet
61  xseed = large
62 
63  ! Flag all nodes that are within nearWallDist as being nearWall
64  do mm = 1, nbocos
65  if (iswalltype(bctype(mm))) then
66 
67  call setboundarypointers(mm, bcfaceid(mm), .false.)
68  ! Loop over the generalized plane
69  do j = jstart, jend
70  do i = istart, iend
71 
72  ! Loop over the 'k' ie offwall direction
73  do k = kstart, kend
74 
75  dist = norm2(planes(k)%x(i, j, :) - planes(1)%x(i, j, :))
76  if (dist < nearwalldist) then
77  planes(k)%nearWall(i, j) = 1
78  nflaggedlocal = nflaggedlocal + 1
79 
80  ! If we made it all the way to the other side
81  ! of the block, set the wall point into the
82  ! xSeed array. Note that we need the *SECOND
83  ! LAST NODE*. This is the value that is
84  ! actually transfered as it is the halo node
85  ! for the other block.
86 
87  if (k == kend - 1) then
88  planes(k)%xSeed(i, j, :) = planes(1)%x(i, j, :)
89  natboundarylocal = natboundarylocal + 1
90  end if
91  end if
92  end do
93  end do
94  end do
95  deallocate (planes)
96  end if
97  end do ! BocoLoop
98  end do
99 
100  ! Determine if any cells made it to the other side of a face. If so, we have to keep going:
101  call mpi_allreduce(natboundarylocal, natboundary, 1, adflow_integer, mpi_sum, &
102  adflow_comm_world, ierr)
103  call echk(ierr, __file__, __line__)
104 
105  ! Iterative loop
106  loopiter = 1
107  parallelsyncloop: do while (natboundary > 0)
108  if (myid == 0) then
109  print *, 'Flag Near Wall Iteration:', loopiter, 'nAtBoundary', natboundary
110  end if
111 
112  ! Reset the counter
113  natboundarylocal = 0
114 
115  ! Exchange the xSeeds after the initial flooding. Set the
116  ! pointers for the comm:
117  do nn = 1, ndom
118  xseed => flowdoms(nn, level, sps)%xSeed
119  flowdoms(nn, level, sps)%realCommVars(1)%var => xseed(:, :, :, 1)
120  flowdoms(nn, level, sps)%realCommVars(2)%var => xseed(:, :, :, 2)
121  flowdoms(nn, level, sps)%realCommVars(3)%var => xseed(:, :, :, 3)
122  end do
123 
124  ! Run the generic integer exchange
125  call whalo1to1realgeneric(3, level, sps, commpatternnode_1st, internalnode_1st)
126 
127  do nn = 1, ndom
128  call setpointers(nn, level, sps)
129 
130  ! Manaully set the additional pointers
131  xseed => flowdoms(nn, level, sps)%xSeed
132  nearwall => flowdoms(nn, level, sps)%nearWall
133 
134  ! Loop over all halos. We won't be selective here, since the
135  ! seeds cound show up anywhere! :-)
136 
137  ! Generic loop over all boundaries
138  do mm = 1, 6
139  call setboundarypointers(mm, mm, .true.)
140 
141  ! Loop over the size of the generalized plane
142  do j = jstart, jend
143  do i = istart, iend
144 
145  ! Determine if we need to do the generalized 'k'
146  ! direction at all. We only need to do if a valid
147  ! seed has shown up in the xSeed:
148  if (planes(0)%xSeed(i, j, 1) < large) then
149 
150  ! Loop over the 'k' ie offwall direction
151  do k = kstart, kend
152 
153  dist = norm2(planes(k)%x(i, j, :) - planes(0)%xSeed(i, j, :))
154  if (dist < nearwalldist .and. planes(k)%nearWall(i, j) == 0) then
155  planes(k)%nearWall(i, j) = 1
156 
157  ! If we made it all the way to the other side
158  ! of the block, copy over the seed for the next exchange.
159  if (k == kend - 1) then
160  planes(k)%xSeed(i, j, :) = planes(0)%xSeed(i, j, :)
161  natboundarylocal = natboundarylocal + 1
162  end if
163  end if
164  end do
165  end if
166  end do
167  end do
168  deallocate (planes)
169  end do
170  end do
171 
172  ! Determine if any cells made it to the other side of a face. If so, we have to keep going:
173  call mpi_allreduce(natboundarylocal, natboundary, 1, adflow_integer, mpi_sum, &
174  adflow_comm_world, ierr)
175  call echk(ierr, __file__, __line__)
176 
177  loopiter = loopiter + 1
178 
179  end do parallelsyncloop
180 
181  ! Deallocate X and XSeed since they are no longer needed. We
182  ! have to hold onto nearWall a little longer since we need to use it
183  ! in initializeOBlock. It will deallocate it.
184  do nn = 1, ndom
185  deallocate (flowdoms(nn, level, sps)%XSeed)
186  end do
187 
188 contains
189  subroutine setboundarypointers(mm, faceID, fullFaces)
190  implicit none
191 
192  integer(kind=intType), intent(in) :: mm, faceID
193  logical, intent(in) :: fullFaces
194 
195  if (.not. fullfaces) then
196  istart = bcdata(mm)%inBeg; iend = bcdata(mm)%inEnd
197  jstart = bcdata(mm)%jnBeg; jend = bcdata(mm)%jnEnd
198  else
199  istart = 1; jstart = 1
200  select case (mm)
201  case (imin, imax)
202  iend = jl; jend = kl
203  case (jmin, jmax)
204  iend = il; jend = kl
205  case (kmin, kmax)
206  iend = il; jend = jl
207  end select
208  end if
209 
210  select case (faceid)
211  case (imin)
212  kstart = 1; kend = il
213  allocate (planes(0:il))
214  planes(0)%xSeed => xseed(0, 1:jl, 1:kl, :)
215  do i = 1, il
216  planes(i)%x => x(i, 1:jl, 1:kl, :)
217  planes(i)%xSeed => xseed(i, 1:jl, 1:kl, :)
218  planes(i)%nearWall => nearwall(i, :, :)
219  end do
220  case (imax)
221  kstart = 1; kend = il
222  allocate (planes(0:il))
223  planes(0)%xSeed => xseed(ie, 1:jl, 1:kl, :)
224  do i = 1, il
225  planes(i)%x => x(il - i + 1, 1:jl, 1:kl, :)
226  planes(i)%xSeed => xseed(il - i + 1, 1:jl, 1:kl, :)
227  planes(i)%nearWall => nearwall(il - i + 1, :, :)
228  end do
229  case (jmin)
230  kstart = 1; kend = jl
231  allocate (planes(0:jl))
232  planes(0)%xSeed => xseed(1:il, 0, 1:kl, :)
233  do j = 1, jl
234  planes(j)%x => x(1:il, j, 1:kl, :)
235  planes(j)%xSeed => xseed(1:il, j, 1:kl, :)
236  planes(j)%nearWall => nearwall(:, j, :)
237  end do
238  case (jmax)
239  kstart = 1; kend = jl
240  allocate (planes(0:jl))
241  planes(0)%xSeed => xseed(1:il, je, 1:kl, :)
242  do j = 1, jl
243  planes(j)%x => x(1:il, jl - j + 1, 1:kl, :)
244  planes(j)%xSeed => xseed(1:il, jl - j + 1, 1:kl, :)
245  planes(j)%nearWall => nearwall(:, jl - j + 1, :)
246  end do
247  case (kmin)
248  kstart = 1; kend = kl
249  allocate (planes(0:kl))
250  planes(0)%xSeed => xseed(1:il, 1:jl, 0, :)
251  do k = 1, kl
252  planes(k)%x => x(1:il, 1:jl, k, :)
253  planes(k)%xSeed => xseed(1:il, 1:jl, k, :)
254  planes(k)%nearWall => nearwall(:, :, k)
255  end do
256  case (kmax)
257  kstart = 1; kend = kl
258  allocate (planes(0:kl))
259  planes(0)%xSeed => xseed(1:il, 1:jl, ke, :)
260  do k = 1, kl
261  planes(k)%x => x(1:il, 1:jl, kl - k + 1, :)
262  planes(k)%xSeed => xseed(1:il, 1:jl, kl - k + 1, :)
263  planes(k)%nearWall => nearwall(:, :, kl - k + 1)
264  end do
265  end select
266  end subroutine setboundarypointers
267 end subroutine flagnearwallcells
subroutine flagnearwallcells(level, sps)
Definition: flagNearWall.F90:2
subroutine setboundarypointers(mm, faceID, fullFaces)
Definition: adtAPI.F90:1
Definition: BCData.F90:1
integer(kind=inttype) jl
integer(kind=inttype) ie
integer(kind=inttype), dimension(:), pointer bcfaceid
integer(kind=inttype) nbocos
integer(kind=inttype), dimension(:), pointer bctype
integer(kind=inttype) je
integer(kind=inttype) ke
real(kind=realtype), dimension(:, :, :, :), pointer x
real(kind=realtype), dimension(:, :, :, :), pointer xseed
integer(kind=inttype) kl
integer(kind=inttype) il
type(commtype), dimension(:), allocatable commpatternnode_1st
type(internalcommtype), dimension(:), allocatable internalnode_1st
integer adflow_comm_world
integer(kind=inttype), parameter imax
Definition: constants.F90:293
integer(kind=inttype), parameter kmin
Definition: constants.F90:296
integer(kind=inttype), parameter jmax
Definition: constants.F90:295
integer(kind=inttype), parameter imin
Definition: constants.F90:292
real(kind=realtype), parameter large
Definition: constants.F90:24
integer(kind=inttype), parameter kmax
Definition: constants.F90:297
integer(kind=inttype), parameter jmin
Definition: constants.F90:294
real(kind=realtype) nearwalldist
Definition: inputParam.F90:883
Definition: utils.F90:1
logical function iswalltype(bType)
Definition: utils.F90:1705
subroutine echk(errorcode, file, line)
Definition: utils.F90:5722
subroutine setpointers(nn, mm, ll)
Definition: utils.F90:3237