bcp Subroutine

public subroutine bcp(p)

Uses

  • proc~~bcp~~UsesGraph proc~bcp modboundary::bcp module~modfields modfields proc~bcp->module~modfields module~modglobal modglobal proc~bcp->module~modglobal module~modmpi modmpi proc~bcp->module~modmpi mpi mpi module~modmpi->mpi

Arguments

Type IntentOptional Attributes Name
real, intent(inout), dimension(ib - ih:ie + ih, jb - jh:je + jh, kb - kh:ke + kh) :: p

Calls

proc~~bcp~~CallsGraph proc~bcp modboundary::bcp proc~excj modmpi::excj proc~bcp->proc~excj mpi_sendrecv mpi_sendrecv proc~excj->mpi_sendrecv

Called by

proc~~bcp~~CalledByGraph proc~bcp modboundary::bcp proc~tderive modpois::tderive proc~tderive->proc~bcp proc~poisson modpois::poisson proc~poisson->proc~tderive program~dalesurban DALESURBAN program~dalesurban->proc~poisson

Contents

Source Code

bcp

Source Code

   subroutine bcp(p)

      use modglobal, only:ib, ie, jb, je, ih, jh, kb, ke, kh, linoutflow, dxfi
      use modfields, only:pres0, up, u0, um, uouttot
      use modmpi, only:excj

      real, dimension(ib - ih:ie + ih, jb - jh:je + jh, kb - kh:ke + kh), intent(inout) :: p !< pressure
      integer i, j, k

      if (linoutflow) then
         do k = kb, ke
            do j = jb, je
               p(ib - 1, j, k) = p(ib, j, k) ! inflow:  dp/dn=0
               pres0(ib - 1, j, k) = pres0(ib, j, k) ! inflow:  dp/dn=0
               p(ie + 1, j, k) = -p(ie, j, k) ! outflow: p=0
               pres0(ie + 1, j, k) = -pres0(ie, j, k) ! outflow: p=0
               up(ie + 1, j, k) = -(u0(ie + 1, j, k) - u0(ie, j, k))*dxfi(ie)*uouttot
            enddo
         enddo
      else
         do k = kb, ke
            do j = jb, je
               p(ib - 1, j, k) = p(ie, j, k)
               p(ie + 1, j, k) = p(ib, j, k)
            enddo
         enddo
      endif

      call excj(p, ib - 1, ie + 1, jb - 1, je + 1, kb - 1, ke + 1) ! cyclic
      call excj(pres0, ib - 1, ie + 1, jb - 1, je + 1, kb - 1, ke + 1) ! cyclic

   end subroutine bcp