Boundary Cells and Boundary Conditions

  • Mohamed Houssem Kasbaoui

Boundary Cells

Specification of boundaries (other than periodic boundary conditions) are handeled by the derived type bc_set described in /path/to/leap/kernel/leapbc.f90. Here, we will describe briefly how to specify boundary cells and types using the example illustrated below.

In this example, we have four distinct regions on each side of the domain.

To define these regions in LEAP, we can use a bc_set. For example, defining the inlet region goes as follows:

use leapKinds
use leapBlock
type(parallel_obj) :: parallel
type(block_obj)    :: block
type(block_obj)    :: bcs

! Initialize the parallel machinery
call parallel%Initialize()

! Define a block as shown in the picture
call block%Initialize(ngc=1,parallel)
! ... (setup the block)

! Initialize the BC utility on this block
call bcs%Initialize(block,parallel)

! Add a new region for the inlet
xlo = [block%pmin(1),block%pmin(2),block%pmin(3)]
xhi = [block%pmin(1),block%pmax(2),block%pmax(3)]
call bcs%Add('Inlet', xlo,xhi, normal = '+x1')

Note that xlo and xhi define bounds of a 2D surface along the inlet. An error will be returned, if xlo and xhi result in a volumetric region. The normal here is '+x1' meaning that it is in the -direction towards increasing values (left to right). This code will result in the first layer of cells tagged as boundary cells with their left side being the boundary.

Note that if we had specified the normal as '-x1', this will result in the right side of the highlighted cells being tagged as boundary, which is incorrect.

Next, we specify the remaining 3 regions:

! Add a new region for the outlet
xlo = [block%pmax(1),block%pmin(2),block%pmin(3)]
xhi = [block%pmax(1),block%pmax(2),block%pmax(3)]
call bcs%Add('Outlet', xlo,xhi, normal = '-x1')

! Add a new region for the bottom wall
xlo = [block%pmin(1),block%pmin(2),block%pmin(3)]
xhi = [block%pmax(1),block%pmin(2),block%pmax(3)]
call bcs%Add('Top', xlo,xhi, normal = '+x2')

! Add a new region for the top wall
xlo = [block%pmin(1),block%pmax(2),block%pmin(3)]
xhi = [block%pmax(1),block%pmax(2),block%pmax(3)]
call bcs%Add('Bottom', xlo,xhi, normal = '-x2')

Boundary Conditions

Now that we have defined the four different regions, we can define boundary conditions on these regions. The main two boundary conditions supported currently in LEAP are the Dirichlet and 0-gradient Neumann boundary conditions.

Starting with the inlet, we can define a Dirichlet boundary condition for a variable U as follows:

real(wp), pointer :: U_BC(:,:,:)
type(extents_obj) :: extents
integer           :: i,j,k

! Add a boundary condition
call bcs%SetBC('Inlet',BC_DIRICHLET, 'U')

! Get extents and a pointer to the BC data
extents = bcs%GetExtents('Inlet')
call bcs%GetBCPointer('Inlet','U',U_BC)

! Set the BC values
do k=extents%lo(3),extents%hi(3)
  do j=extents%lo(2),extents%hi(2)
    do i=extents%lo(1),extents%hi(1)
      U_BC(i,j,k) =  block%ym(j)
    end do
  end do
end do

To add a 0-gradient Neumann BC on the outlet, we can use the following:

! Add a boundary condition
call bcs%SetBC('Outlet',BC_NEUMANN, 'U')

Writing/Reading Boundary Conditions

Reading and wring to disk is straightforward with the following commands:

  ! Load the regions and different boundary conditions from file
  call bcs%Read('bcs.hdf5')

  ! Write regions and boundary conditions to file
  ! then finalize.
  call bcs%Write('bcs.hdf5')
  call bcs%Finalize()