subroutine createwalls
use modglobal, only:ib, ie, jb, je, jgb, jge, kb, ke, jmax, nblocks, &
nsv, cexpnr, ifinput, libm, ih, kh, iwallmom, iwalltemp, iwallmoist, rslabs, bldT
use modsurfdata, only:thls, z0h, z0, thvs
use modfields, only:sv0, svm, thl0, thlm, qtp, qt0, IIc, IIu, IIv, IIw, IIct, IIwt, IIcs, IIus, IIvs, IIws
use modmpi, only:myid, comm3d, mpierr, MPI_INTEGER, MPI_DOUBLE_PRECISION, MY_REAL, nprocs, cmyid, &
MPI_REAL8, MPI_REAL4, MPI_SUM
use initfac, only:block
integer n, nn, pn, mn, jbeg, jend, nxn, nxs, nyn, nzn, nzs, iu, il, ju, jl, ku, kl, sc, &
i, k, dbi, dbj, dbk
integer :: IIcl(kb:ke + kh), IIul(kb:ke + kh), IIvl(kb:ke + kh), IIwl(kb:ke + kh)
integer :: IIcd(ib - ih:ie + ih, kb:ke + kh)
integer :: IIwd(ib - ih:ie + ih, kb:ke + kh)
character(80) chmess, name2
if (.not. libm) return
! check if walls are at least 2 cells in each dimension
do n = 1, nblocks
dbi = block(n, 2) - block(n, 1)
dbj = block(n, 4) - block(n, 3)
dbk = block(n, 6) - block(n, 5)
if (any((/dbi, dbj, dbk/) < 1)) then
write(6, *) "blocks not at least 2 cells in each dimension, or upper limit < lower limit"
!stop !ils13 19.07.17, don't stop for now
end if
end do
! For all blocks set the internal concentrations to zero and internal
! temperature to building temperature
do n = 1, nblocks
il = block(n, 1)
iu = block(n, 2)
kl = block(n, 5)
ku = block(n, 6)
jl = block(n, 3) - myid*jmax
ju = block(n, 4) - myid*jmax
if ((ju < jb - 1) .or. (jl > je + 1)) then ! The block is entirely out of this partition
cycle
end if
if (ju > je) ju=je !tg3315 and bss116 added 23.10.18 as bad allocation otherwise.
if (jl < jb) jl=jb
do sc = 1, nsv
!sv0(il:iu, jl:ju, kl:ku, sc) = svprof(kl:ku) !internal ! tg3315 commented to avoid flux at startup
!svm(il:iu, jl:ju, kl:ku, sc) = svprof(kl:ku) !internal
end do
thl0(il:iu, jl:ju, kl:ku) = bldT !internal ! make sure bldT is equal to init thl prof
thlm(il:iu, jl:ju, kl:ku) = bldT !internal
end do
nxwall = 0
do n = 1, nblocks ! first x and z walls
if (block(n, 4) < jb + myid*jmax) then ! no x-wall/z-wall in the range of this proc
cycle
elseif (block(n, 3) > je + myid*jmax) then ! no x-wall/z-wall in the range of this proc
cycle
else ! x-wall/z-wall found on this proc
nxwall = nxwall + 1
end if
end do
allocate (ixwall(nxwall)) !allocate the list that stores the indeces of the blocks on this cpu
k = 1
do n = 1, nblocks ! save indeces of the found x/z-walls by re-iterating
if (block(n, 4) < jb + myid*jmax) then ! no x-wall/z-wall in the range of this proc
cycle
elseif (block(n, 3) > je + myid*jmax) then ! no x-wall/z-wall in the range of this proc
cycle
else
ixwall(k) = n ! save index of block which is on this processor
k = k + 1
end if
end do
!!new approach both y walls#################################################
!!store index of block and index of the wall (since block might not be on this cpu, but is needed for x and z coords)
!!check if wall is on next cpu but not on this
!!check if wall is on last cpu but not on first (periodicity in y)
!check if wall is on first cpu but not on last (periodicity in y)
!!check if wall is on this cpu and another one on the next (i.e. both blocks end at cpu boundary, but touch each other)
nyminwall = 0
nypluswall = 0
do n = 1, nblocks
jl = block(n, 3) - myid*jmax
ju = block(n, 4) - myid*jmax
!IMPORTANT: THESE LINES OF CODE SHOULD BE HERE BUT CAUSE TROUBLE! MAKE SURE BLOCK DOES NOT TOUCH BOUNDARY (EXCECPT ALL FLOORS)
!SEE ALSO BELOW! (Like 40 lines or so)
if ((myid == 0) .and. (block(n, 4) == jge)) then ! periodicity!
nypluswall = nypluswall + 1
else if ((block(n, 3) == jgb) .and. (myid == (nprocs - 1))) then ! periodicity!
nyminwall = nyminwall + 1
end if
if ((ju < (jb - 1)) .or. (jl > (je + 1))) then
cycle
end if
if (ju == (jb - 1)) then !block on previous cpu, north wall on this
nypluswall = nypluswall + 1 !
cycle
end if
if (jl == (je + 1)) then
nyminwall = nyminwall + 1 !block on next cpu, southwall on this
cycle
end if
if ((ju < je) .and. (ju >= jb)) then !block & northwall on this cpu
nypluswall = nypluswall + 1
end if
if ((jl > jb) .and. (jl <= je)) then !block & southwall on this cpu
nyminwall = nyminwall + 1
end if
end do
allocate (iyminwall(1:nyminwall, 1:2)) !two indeces to store wall index and block index
allocate (iypluswall(1:nypluswall, 1:2))
iyminwall(:, 1) = 0
iyminwall(:, 2) = 0
iypluswall(:, 1) = 0
iypluswall(:, 2) = 0
pn = 1
mn = 1
do n = 1, nblocks
jl = block(n, 3) - myid*jmax
ju = block(n, 4) - myid*jmax
!IMPORTANT: THESE LINES OF CODE SHOULD BE HERE BUT CAUSE TROUBLE! MAKE SURE BLOCK DOES NOT TOUCH BOUNDARY (EXCECPT ALL FLOORS)
if ((myid == 0) .and. (block(n, 4) == jge)) then ! periodicity!
iypluswall(pn, 1) = n
iypluswall(pn, 2) = jb
pn = pn + 1
else if ((block(n, 3) == jgb) .and. (myid == (nprocs - 1))) then ! periodicity!
iyminwall(mn, 1) = n
iyminwall(mn, 2) = je
mn = mn + 1
end if
if ((ju < (jb - 1)) .or. (jl > (je + 1))) then
cycle
end if
if (ju == (jb - 1)) then !block on previous cpu, north wall on this
iypluswall(pn, 1) = n
iypluswall(pn, 2) = jb
pn = pn + 1
cycle
end if
if (jl == (je + 1)) then !block on next cpu, south wall on this
iyminwall(mn, 1) = n
iyminwall(mn, 2) = je
mn = mn + 1
cycle
end if
if ((ju < je) .and. (ju >= jb)) then !block & northwall on this cpu !ILS13, 5.12.17 following Tom
iypluswall(pn, 1) = n
iypluswall(pn, 2) = ju + 1
pn = pn + 1
end if
if ((jl > jb) .and. (jl <= je)) then !block & southwall on this cpu
iyminwall(mn, 1) = n
iyminwall(mn, 2) = jl - 1
mn = mn + 1
end if
end do
end subroutine createwalls