mus_bc_header_module Module

This module contains boundary type definitions, different boundary treatment types, abtract interface and routine to load boundary table from config file

The boundary conditions are stored in boundary type definitions for each boundary which is defined in the main lua file.

 boundary_condition = {  }

The definitions there must equal the boundary definition within the mesh on disk, in terms of the amount of boundary conditions and the labels. A detailed description on the implementation details are given in tem_bc_module "boundary condition implementation"

This module will be used as a header to use boundary definitions in other modules.


Uses

Used by

  • module~~mus_bc_header_module~~UsedByGraph module~mus_bc_header_module mus_bc_header_module module~mus_bndforce_module mus_bndForce_module module~mus_bndforce_module->module~mus_bc_header_module module~mus_bc_passivescalar_module mus_bc_passiveScalar_module module~mus_bc_passivescalar_module->module~mus_bc_header_module module~mus_bc_general_module mus_bc_general_module module~mus_bc_general_module->module~mus_bc_header_module module~mus_bc_fluid_noneqexpol_module mus_bc_fluid_nonEqExpol_module module~mus_bc_fluid_noneqexpol_module->module~mus_bc_header_module module~mus_bc_fluid_wall_module mus_bc_fluid_wall_module module~mus_bc_fluid_wall_module->module~mus_bc_header_module module~mus_weights_module mus_weights_module module~mus_weights_module->module~mus_bc_header_module module~mus_scheme_type_module mus_scheme_type_module module~mus_scheme_type_module->module~mus_bc_header_module module~mus_bc_nernstplanck_module mus_bc_nernstPlanck_module module~mus_bc_nernstplanck_module->module~mus_bc_header_module module~mus_construction_module mus_construction_module module~mus_construction_module->module~mus_bc_header_module module~mus_field_module mus_field_module module~mus_field_module->module~mus_bc_header_module module~mus_bc_species_module mus_bc_species_module module~mus_bc_species_module->module~mus_bc_header_module module~mus_bc_fluid_turbulent_module mus_bc_fluid_turbulent_module module~mus_bc_fluid_turbulent_module->module~mus_bc_header_module module~mus_bc_var_module mus_bc_var_module module~mus_bc_var_module->module~mus_bc_header_module module~mus_bc_fluid_module mus_bc_fluid_module module~mus_bc_fluid_module->module~mus_bc_header_module module~mus_bc_fluid_experimental_module mus_bc_fluid_experimental_module module~mus_bc_fluid_experimental_module->module~mus_bc_header_module module~mus_dynloadbal_module mus_dynLoadBal_module module~mus_dynloadbal_module->module~mus_bc_header_module module~mus_bc_poisson_module mus_bc_poisson_module module~mus_bc_poisson_module->module~mus_bc_header_module module~mus_connectivity_module mus_connectivity_module module~mus_connectivity_module->module~mus_bc_header_module

Contents


Abstract Interfaces

abstract interface

Interface definition for the boundary condition subroutines

  • private subroutine boundaryRoutine(me, state, bcBuffer, globBC, levelDesc, tree, nSize, iLevel, sim_time, neigh, layout, fieldProp, varPos, nScalars, varSys, derVarPos, physics, iField, mixture)

    This abstract interface defines the interface for boundary routines. All boundary routines that need to be called via fnct function pointer need to implement this interface. In case of new boundary routines, mark them with a comment reffering to fnct in order to find the routine in case the interface needs to be changed.

    Arguments

    TypeIntentOptionalAttributesName
    class(boundary_type) :: me

    field boundary type

    real(kind=rk), intent(inout) :: state(:)

    State array

    real(kind=rk), intent(in) :: bcBuffer(:)

    state values of boundary elements of all fields of iLevel

    type(glob_boundary_type), intent(in) :: globBC

    scheme global boundary type

    type(tem_levelDesc_type), intent(in) :: levelDesc

    iLevel descriptor

    type(treelmesh_type), intent(in) :: tree

    global treelm mesh

    integer, intent(in) :: nSize

    size of state array ( in terms of elements )

    integer, intent(in) :: iLevel

    level which invokes boundary

    type(tem_time_type), intent(in) :: sim_time

    global time information

    integer, intent(in) :: neigh(:)

    global parameters

    type(mus_scheme_layout_type), intent(in) :: layout

    scheme layout

    type(mus_field_prop_type), intent(in) :: fieldProp

    Fluid property

    integer, intent(in) :: varPos(:)

    pointer to field variable in the state vector

    integer, intent(in) :: nScalars

    number of Scalars in the scheme var system

    type(tem_varSys_type), intent(in) :: varSys

    scheme variable system

    type(mus_derVarPos_type), intent(in) :: derVarPos

    position of derived quantities in varsys

    type(mus_physics_type), intent(in) :: physics

    scheme global boundary type

    integer, intent(in) :: iField

    current field

    type(mus_mixture_type), intent(in) :: mixture

    mixture info

abstract interface

Interface definition for the boundary condition subroutines

  • private subroutine mus_proc_calcBndForce(me, bndForce, posInBndID, globBC, currState, levelDesc, nSize, iLevel, neigh, layout, nScalars)

    This abstract interface defines the interface for bndForce calculation. calcBndForce in order to find the routine in case the interface needs to be changed.

    Arguments

    TypeIntentOptionalAttributesName
    class(boundary_type), intent(in) :: me

    field boundary type

    real(kind=rk), intent(inout) :: bndForce(:,:)

    bndForce to fill

    integer, intent(in) :: posInBndID(:)
    type(glob_boundary_type), intent(in) :: globBC

    scheme global boundary type

    real(kind=rk), intent(in) :: currState(:)

    current state array to access post-collision values

    type(tem_levelDesc_type), intent(in) :: levelDesc

    iLevel descriptor

    integer, intent(in) :: nSize

    size of state array ( in terms of elements )

    integer, intent(in) :: iLevel

    level which invokes boundary

    integer, intent(in) :: neigh(:)

    global parameters

    type(mus_scheme_layout_type), intent(in) :: layout

    scheme layout

    integer, intent(in) :: nScalars

    number of Scalars in the scheme var system


Derived Types

type, public :: glob_boundary_type

Information about boundary elements

Read more…

Components

TypeVisibilityAttributesNameInitial
character(len=LabelLen), private :: label

name of this boundary

type(bc_elems_type), private, allocatable:: elemLvl(:)

elements list of this boundary in different refinement level size: minLevel -> maxLevel

logical, private :: treat_corner

Collect corner BC i.e elements which are interesected by multiple boundaries only for moments based BC

type(mus_bc_corner_type), private :: cornerBC

list corner boundary elements

integer, private, allocatable:: nElems(:)

number of local (this process) boundary elements per level including ghostFromCoarser. GhostFromFiner are interpolated but ghostFromCoarser needs correct boundary value in fine level sub-iteration.

integer, private, allocatable:: nElems_Fluid(:)

number of local (this process) boundary elements per level wihtout Ghost elements

integer, private, allocatable:: nElems_totalLevel(:)

number of total (total processes) boundary elements per level

integer, private :: nElems_local =0

number of local boundary elements

integer, private :: nElems_total =0

number of total boundary elements

logical, private :: hasQVal =.false.

whether this boundary requires qVal

logical, private :: qValInitialized =.false.

has qVal initialized with default qVal = 0.5

logical, private :: isWall =.false.

if this is a wall or symmetry BC (implicit treated during streaming)

integer, private :: normal(3)

Average normal direction of this boundary along layout%prevailDir Normal in elemLvl is normal direction per element and this is normal per boundary

integer, private :: normalInd

which is the index in the stencil corresponding to the normal vector?

type, public :: boundary_type

Boundary treatment data type. It include boundary treatment for only ONE field. For the same boundary elements, each field may have a different treatment. It includes function pointers and neighbor elements required by this boundary treatment

Components

TypeVisibilityAttributesNameInitial
integer, private :: bc_id

Boundary ID

character(len=LabelLen), private :: BC_kind

kind of this boundary

character(len=LabelLen), private :: label

name of this boundary

integer, private :: nNeighs =0

number of neighbors fluid element needed for extrapolation to get data at boundary element

logical, private :: curved =.false.

Curved boundaries are treated with qValues and extrapolation of nonEq from 2nd neighbor

type(bc_neigh_type), private, allocatable:: neigh(:)

Level wise neighbor information of boundary elements Each field may have a different handlement on the same boundary Thus requires a different neighborhood allocated in mus_alloc_fieldBC

type(bc_field_elems_type), private, allocatable:: elemLvl(:)

Level wise element information of boundary elements like stencil position in array of stencils, moments type for momentBC and LODI array for NRBC allocated in mus_alloc_fieldBC

procedure(boundaryRoutine), private, pointer:: fnct=> null()

function pointer to boundary routine to call

procedure(mus_proc_calcBndForce), private, pointer:: calcBndForce=> null()

function pointer to compute bndForce on wall

type(bc_states_type), private :: BC_states

musubi bc state variables

type(grw_stringkeyvaluepairarray_type), private :: varDict

Dictionary of boundary variables with varDict%val()%key is the name of boundary variable and varDict%val()%value is the name of spacetime function variable

integer, private :: order

order of boundary which determines number of neighbor need for this boundary

real(kind=rk), private :: qVal

standard q-value for higher order boundaries without any definition in seeder

type(bc_nrbc_type), private :: nrbc

Non-reflective boundary type

real(kind=rk), private :: slip_fac

slip wall factor

type(bc_blackbox_membrane_type), private :: blackbox_mem

black box membrane boundary

logical, private :: useComputeNeigh =.false.

whether this boundary needs neighbors on compute stencil

logical, private :: evalBcVar_link =.false.

Is true if this boundary requires boundary variables on boundary surface linkwise. Required to compute real coordinates points on boundary surface

logical, private :: requireNeighBufPre_nNext =.false.

Is true only if boundary requires neighbuf Pre-collision values from next time step

logical, private :: requireNeighBufPre =.false.

Is true only if boundary requires neighbuf Pre-collision values from current time step

logical, private :: requireNeighBufPost =.false.

Is true only if boundary requires neighbuf Post-collision values

logical, private :: treat_corner =.false.

Collect corner BC i.e elements which are interesected by multiple boundaries only for moments based BC

type(array_type), private, allocatable:: links(:)

Level-wise links that are to be updated size: minLevel:maxLevel allocated in mus_alloc_fieldBC

type(bc_wall_bouzidi_type), private, allocatable:: bouzidi(:)

Level wise wall bouzidi data size: minLevel:maxLevel

type(bc_inlet_bouzidi_type), private, allocatable:: inletUbbQVal(:)

Level wise inletUbbQVal data size: minLevel:maxLevel

type(bc_inlet_bouzidi_type), private, allocatable:: inletBFL(:)

Level wise inletBfl data size: minLevel:maxLevel

type(bc_nonEqExpol_type), private, allocatable:: nonEqExpol(:)

Level wise nonEquilibrium extrapolation data for link-wise implementation size: minLevel:maxLevel

type(bc_outlet_type), private, allocatable:: outletExpol(:)

Level wise outletExpol data size: minLevel:maxLevel

type(mus_turb_wallFunc_type), private :: turbwallFunc

Turbulent wall model type contains function pointers to compute friction and stream-wise velocity. Also stores friction velocity and turbulent viscosity on boundary elements

type, private :: bc_moments_type

information needed for moments based boundary condition

Components

TypeVisibilityAttributesNameInitial
integer, private :: nUnKnownPdfs

Number of unknown state links to update

integer, private, allocatable:: knownMom_pos(:)

known moments position in moments vector

real(kind=rk), private, allocatable:: unKnownPdfs_MatInv(:,:)

inverse matrix of unknown pdf matrix

type, private :: bc_elems_type

Level wise boundary elements information

Components

TypeVisibilityAttributesNameInitial
type(dyn_intarray_type), private :: elem

Positions in levelDesc total list Its purpose is to get the treeID of BC elements size: globBC%nElems to use: levelDesc(iLevel)%total( globBC%elemLvl(iLevel)%elem%val(iElem) )

type(grw_intarray_type), private :: posInBcElemBuf

Position of this boundary in the bc_elemBuffer bc_elemBuffer is growing array initial size: globBC%nElems It is initiated in mus_init_bc_elems of mus_bc_header_module Only non-wall BC needs elemBuffer

type(grw_int2darray_type), private :: normal

Normal vector for each boundary element pointing into the domain. normal is a growing array with size: normal%val(3, nElems)

type(grw_intarray_type), private :: normalInd

which is the index in the stencil corresponding to the normal vector?

type(grw_logical2darray_type), private :: bitmask

bit mask for each node holding directions which have to be updated. The bitmask points into the incoming direction into the flow domain, which actually we want to update * For PUSH, we write to the incoming position, as the kernel reads it from there without propagation. * For PULL, we need to write to the inverse direction, as the kernel performs a bounce back before reading it. However, this bounced back direction actually comes from the non-existent boundary element and would point into the incoming direction, so the value has to be treated and set as if it points really into the incoming direction. bitmask is a growing array, the values are in bitmask%val(:,:) 1st index size is QQN since center is not treated for bcElems

type(grw_real2darray_type), private :: qVal

The q-Values for the exact wall positions Its size: QQN, nElemsBC It is allocated in routine: allocateBCList assigned in routine: assignBCList

type, private :: mus_bc_corner_type

contains information needed to treat corner nodes or nodes intersected by more than one boundary type

Components

TypeVisibilityAttributesNameInitial
integer, private, allocatable:: nElems(:)

number of local boundary elements per level

type(bc_elems_type), private, allocatable:: elemLvl(:)

boundary neighbor in different refinement level

integer, private, allocatable:: bcid(:,:)

boundary id in each direction for each corner node size layout%stencil%QQ, nElems

type, private :: bc_wall_bouzidi_type

Wall bouzidi data for one level

Components

TypeVisibilityAttributesNameInitial
real(kind=rk), private, allocatable:: qVal(:)

size: links(iLevel)%nVals

real(kind=rk), private, allocatable:: cIn(:)
real(kind=rk), private, allocatable:: cOut(:)
real(kind=rk), private, allocatable:: cNgh(:)
integer, private, allocatable:: inPos(:)

set in routine: mus_set_bouzidi

integer, private, allocatable:: outPos(:)
integer, private, allocatable:: nghPos(:)

type, private :: bc_inlet_bouzidi_type

Provides bouzidi coefficients and q-values needed for link-wise implementation of certain inlet boundary conditions.

Components

TypeVisibilityAttributesNameInitial
real(kind=rk), private, allocatable:: qVal(:)

size: links(iLevel)%nVals

real(kind=rk), private, allocatable:: cIn(:)
real(kind=rk), private, allocatable:: cOut(:)
real(kind=rk), private, allocatable:: cNgh(:)
integer, private, allocatable:: inPos(:)

set in routine: mus_set_bouzidi

integer, private, allocatable:: outPos(:)
integer, private, allocatable:: nghPos(:)
real(kind=rk), private, allocatable:: cVel(:)

size: links(iLevel)%nVals

integer, private, allocatable:: iDir(:)

set in routine: mus_set_bouzidi

integer, private, allocatable:: posInBuffer(:)

type, private :: bc_nonEqExpol_type

Provides coefficients needed for link-wise implementation of nonEquilibrium extrapolation scheme for boundary conditions.

Components

TypeVisibilityAttributesNameInitial
real(kind=rk), private, allocatable:: c_w(:)

size: links(iLevel)%nVals

real(kind=rk), private, allocatable:: c_f(:)
real(kind=rk), private, allocatable:: c_ff(:)
real(kind=rk), private, allocatable:: c_neq_f(:)
real(kind=rk), private, allocatable:: c_neq_ff(:)
integer, private, allocatable:: iDir(:)

set in routine: mus_set_bouzidi

integer, private, allocatable:: posInBuffer(:)
integer, private, allocatable:: posInNeighBuf(:)
integer, private, allocatable:: posInBCelems(:)

type, private :: bc_outlet_type

Provides needed positions for link-wise implementation of outlet boundary condition outlet_expol

Components

TypeVisibilityAttributesNameInitial
integer, private, allocatable:: statePos(:)

size: links(iLevel)%nVals

integer, private, allocatable:: iElem(:)

type, private :: bc_nrbc_type

variable definition for non-reflective type of boundary conditions\n These boundary condition is taken from the paper:

Read more…

Components

TypeVisibilityAttributesNameInitial
real(kind=rk), private :: kappa

specific heat ratio appears in eq 3: cs = sqrt( kappa * R * T )

real(kind=rk), private :: K_mod

from eq 17: k = sigma * ( 1 - Ma^2 ) * cs / L

real(kind=rk), private :: cs_mod

used in eq 17

real(kind=rk), private :: sigma

used in eq 17

real(kind=rk), private :: lodi_length

length between inlet and outlet, represented by L in the paper

real(kind=rk), private :: Ma_L

Lattice Mach number characterised by lattice flow velocity

type, private :: bc_blackbox_membrane_type

Contains variables for black-box membrane boundary condition

Components

TypeVisibilityAttributesNameInitial
character(len=labelLen), private :: label

type of membrance (AEM or CEM)

real(kind=rk), private :: transNr

transference number defining transfer property of ionic species on the membrance KM

Read more…
real(kind=rk), private :: solvTrans_coeff

osmotic solvent transport coefficient

type, private :: bc_states_type

contains all possible bc state variables needed for boundary

Components

TypeVisibilityAttributesNameInitial
type(tem_bc_state_type), private :: velocity

Velocity at boundary

type(tem_bc_state_type), private :: massFlowRate

mass flow rate at boundary node for mfr (mass flow rate) boundary condition

type(tem_bc_state_type), private :: pressure

Pressure at boundary

type(tem_bc_state_type), private :: moleFrac

mole fraction for species boundary

type(tem_bc_state_type), private :: moleFlux

molar flux for species (c_i*u_i)

type(tem_bc_state_type), private :: moleDens

mole density for species (c_i)

type(tem_bc_state_type), private :: moleDiff_flux

molar diffusion flux for species

type(tem_bc_state_type), private :: pdf

probability density function for boundary

type(tem_bc_state_type), private :: potential

potential for Poisson boundary

type(tem_bc_state_type), private :: surChargeDens

surface charge density for Poisson Neumann boundary

type, private :: array_type

Simple array type which include array and its size

Components

TypeVisibilityAttributesNameInitial
integer, private, allocatable:: val(:)

array

integer, private :: nVals

size

type, private :: bc_neigh_type

Level wise information about neighbor of boundary elements for a single field

Components

TypeVisibilityAttributesNameInitial
integer, private, allocatable:: posInState(:,:)

Neighbor position in level-wise state array size: nNeighs, nElems allocated in routine: setFieldBCNeigh, mus_construction_module

real(kind=rk), private, allocatable:: neighBufferPre_nNext(:,:)

Pre-collision state values of neighbors on normal direction on next time step size: nNeighs, nElems*stencil%QQ It is allocated in routine update_BClists

real(kind=rk), private, allocatable:: neighBufferPre(:,:)

Pre-collision state values of neighbors on normal direction on previous time step size: nNeighs, nElems*stencil%QQ It is allocated in routine update_BClists

real(kind=rk), private, allocatable:: neighBufferPost(:,:)

Post-collision state values of neighbors on normal direction size: nNeighs, nElems*stencil%QQ It is allocated in routine update_BClists It is filled in fill_neighBuffer

real(kind=rk), private, allocatable:: computeNeighBuf(:)

Post-collision state values of neighbors on all directions at current time step. Its also a Pre-collision state values of an element next to boundary at next time step size: nElems * computeStencil%QQ It is allocated in routine update_BClists (mus_construction_module), filled in routine: fill_computeNeighBuf (mus_bc_general_module). it use AOS layout always!

type, private :: bc_field_elems_type

Contains field specific element information for each level

Components

TypeVisibilityAttributesNameInitial
integer, private, allocatable:: stencilPos(:)

Position of stencil of each element in array of stencils in layout. Unique stencil label for boundary stencils are created with boundary label and stencil%cxDir therefore each stencil is limited to one boundary type Dimension: globBC%nElems(iLevel) allocated in mus_build_BCStencils (mus_construction_module)

real(kind=rk), private, allocatable:: lodi(:,:,:)

LODI array for the NRBC Dimension: 4, 3, globBC%nElems(iLevel) allocated in init_nrbc routine

type(bc_moments_type), private, allocatable:: moments(:)

defined moments position in moments matrix and inverse of unknown pdf matrix Dimension: bc%nElems(iLevel)

integer, private, allocatable:: posInNghElems(:)

Position in LevelDesc%neigh%NghElems size: globBC%nElems(iLevel) allocated in mus_build_BCStencils routine, mus_construction_module


Subroutines

public subroutine mus_load_bc(me, bc_prop, conf, parent, varSys, stencil)

Read in the boundary conditions from the LUA parameter file.\n

Arguments

TypeIntentOptionalAttributesName
type(boundary_type), allocatable:: me(:)

field boundary type

type(tem_BC_prop_type), intent(in) :: bc_prop

boundary data from mesh

type(flu_state) :: conf

lua flu state

integer, intent(in), optional :: parent

bc parent handle

type(tem_varSys_type), intent(inout) :: varSys

Global variable system required to append annoymous boundary variables

type(tem_stencilHeader_type) :: stencil

stencil information

public subroutine mus_init_bc_elems(me, nElems, QQN, hasQval)

Arguments

TypeIntentOptionalAttributesName
type(bc_elems_type), intent(inout) :: me
integer, intent(in) :: nElems
integer, intent(in) :: QQN

Number of direction excluding center since center is never treated for bc elems

logical, intent(in) :: hasQval

public subroutine rearrange_bc_elems(minLevel, maxLevel, nValid, posInBCElem, nElems, elemLvl)

It removes non-valid elements while still maintaining the origianl order. The given bc elements (elems) contains both valid and non-valid elements. Position of valid elements are given by posInBCElem. Valid elements are moved towards the start of the elems so that they become continuous in the elems.

Arguments

TypeIntentOptionalAttributesName
integer, intent(in) :: minLevel
integer, intent(in) :: maxLevel
integer, intent(in) :: nValid(minLevel:maxLevel)

number of valid BC elements

integer, intent(in) :: posInBCElem(maxval(nElems),minLevel:maxLevel)

their position in bc_elems_type

integer, intent(in) :: nElems(minLevel:maxLevel)

BC elements information

type(bc_elems_type), intent(inout) :: elemLvl(minLevel:maxLevel)

BC elements information

public subroutine check_solid_in_bc(minLevel, maxLevel, levelPointer, levelDesc, nElems, elemLvl, nValid, posInBCElem)

It count valid (non-solid) elements in BC elements list. Input: minLevel, maxLevel LevelPointer LevelDesc nElems - number of BC elements elems - positions of BC elements in tree or levelPointer Output: nValid - number of valid BC elements posInBCElem - positions of valid elements in BC elements list

Arguments

TypeIntentOptionalAttributesName
integer, intent(in) :: minLevel
integer, intent(in) :: maxLevel
integer, intent(in) :: levelPointer(:)

Level Pointer ( position of a element in level desc )

type(tem_levelDesc_type), intent(in) :: levelDesc(minLevel:maxLevel)

Level Descriptor

integer, intent(in) :: nElems(minLevel:maxLevel)

number of BC elements

type(bc_elems_type), intent(in) :: elemLvl(minLevel:maxLevel)

BC elements list that contains their position in levelPointer

integer, intent(out) :: nValid(minLevel:maxLevel)

number of valid (non-solid) elements

integer, intent(out) :: posInBCElem(maxval(nElems),minLevel:maxLevel)

positions of valid elements in globBC elements list

public subroutine debug_glob_boundary_type(nBCs, minLevel, maxLevel, me)

Arguments

TypeIntentOptionalAttributesName
integer, intent(in) :: nBCs
integer, intent(in) :: minLevel
integer, intent(in) :: maxLevel
type(glob_boundary_type), intent(in) :: me(nBCs)

public subroutine mus_set_posInNghElems(minLevel, maxLevel, nStencils, globBC, BC)

Set BC elements positions in LevelDesc%neigh%nghElems

Arguments

TypeIntentOptionalAttributesName
integer, intent(in) :: minLevel
integer, intent(in) :: maxLevel
integer, intent(in) :: nStencils
type(glob_boundary_type), intent(in) :: globBC
type(boundary_type), intent(inout) :: BC

public subroutine mus_set_bcLinks(iField, QQ, QQN, nScalars, nElems, elemLvl, nSize, neigh, links)

Arguments

TypeIntentOptionalAttributesName
integer, intent(in) :: iField
integer, intent(in) :: QQ
integer, intent(in) :: QQN
integer, intent(in) :: nScalars
integer, intent(in) :: nElems

number of BC elements

type(bc_elems_type), intent(in) :: elemLvl
integer, intent(in) :: nSize

number of BC elements

integer, intent(in) :: neigh(:)
type(array_type), intent(out) :: links

public subroutine mus_set_bouzidi(iLevel, QQ, QQN, nScalars, globBC, cxDirInv, varPos, bouzidi)

Set necessary data for Wall Bouzidi BC bouzidi should be allocated beforehand

Arguments

TypeIntentOptionalAttributesName
integer, intent(in) :: iLevel
integer, intent(in) :: QQ
integer, intent(in) :: QQN
integer, intent(in) :: nScalars
type(glob_boundary_type), intent(in) :: globBC
integer, intent(in) :: cxDirInv(QQ)
integer, intent(in) :: varPos(:)
type(bc_wall_bouzidi_type), intent(inout) :: bouzidi

cIn, cOut, cNgh should have size of nLinks

public subroutine mus_alloc_bouzidi(me, nVals, minLevel, maxLevel)

Arguments

TypeIntentOptionalAttributesName
type(bc_wall_bouzidi_type), allocatable:: me(:)
integer, intent(in) :: nVals(minLevel:maxLevel)
integer, intent(in) :: minLevel
integer, intent(in) :: maxLevel

public subroutine mus_alloc_fieldBC(me, minLevel, maxLevel)

Arguments

TypeIntentOptionalAttributesName
type(boundary_type) :: me
integer, intent(in) :: minLevel
integer, intent(in) :: maxLevel

public subroutine mus_set_inletUbb(me, tree, stencil, nScalars, globBC, levelDesc, varPos, nLinks, minLevel, maxLevel)

Set necessary data for BC velocity_bounceback_qval

Read more…

Arguments

TypeIntentOptionalAttributesName
type(bc_inlet_bouzidi_type), intent(out), allocatable:: me(:)

setting type for bc

type(treelmesh_type), intent(in) :: tree

using mesh information

type(tem_stencilHeader_type), intent(in) :: stencil

for directions

integer, intent(in) :: nScalars

number of scalars

type(glob_boundary_type), intent(in) :: globBC

for number of elements in boundary and position in buffer

type(tem_levelDesc_type), intent(in) :: levelDesc(minLevel:maxLevel)

level descriptor

integer, intent(in) :: varPos(:)

for position of outgoing elements

integer, intent(in) :: nLinks(minLevel:maxLevel)

for linkwise treatment

integer, intent(in) :: minLevel

minimum and maximum level

integer, intent(in) :: maxLevel

minimum and maximum level

public subroutine mus_set_inletBfl(me, tree, stencil, nScalars, globBC, levelDesc, varPos, nLinks, minLevel, maxLevel)

Bitmask is true for incoming direction

Arguments

TypeIntentOptionalAttributesName
type(bc_inlet_bouzidi_type), intent(out), allocatable:: me(:)

setting type for bc

type(treelmesh_type), intent(in) :: tree

using mesh information

type(tem_stencilHeader_type), intent(in) :: stencil

for directions

integer, intent(in) :: nScalars

number of scalars

type(glob_boundary_type), intent(in) :: globBC

for number of elements in boundary and position in buffer

type(tem_levelDesc_type), intent(in) :: levelDesc(minLevel:maxLevel)

level descriptor

integer, intent(in) :: varPos(:)

for position of elements

integer, intent(in) :: nLinks(minLevel:maxLevel)

for linkwise treatment

integer, intent(in) :: minLevel

minimum and maximum level

integer, intent(in) :: maxLevel

minimum and maximum level

public subroutine mus_set_nonEqExpol(me, curved, tree, stencil, nScalars, globBC, bc_neigh, pdf, levelDesc, varPos, nLinks, minLevel, maxLevel)

Linkwise non-equilibrium extrapolation (can handle curved walls)

Read more…

Arguments

TypeIntentOptionalAttributesName
type(bc_nonEqExpol_type), intent(out), allocatable:: me(:)

setting type for bc

logical, intent(in) :: curved

Curved or straight boundary

type(treelmesh_type), intent(in) :: tree

using mesh information

type(tem_stencilHeader_type), intent(in) :: stencil

for directions

integer, intent(in) :: nScalars

number of scalars

type(glob_boundary_type), intent(in) :: globBC

scheme global boundary type

type(bc_neigh_type), intent(in) :: bc_neigh(minLevel:maxLevel)

boundary neighbor

type(pdf_data_type), intent(in) :: pdf(tree%global%minlevel:tree%global%maxlevel)

contains global state vector

type(tem_levelDesc_type), intent(in) :: levelDesc(minLevel:maxLevel)

level descriptor

integer, intent(in) :: varPos(:)

for position of elements

integer, intent(in) :: nLinks(minLevel:maxLevel)

for linkwise treatment

integer, intent(in) :: minLevel

minimum and maximum level

integer, intent(in) :: maxLevel

minimum and maximum level

public subroutine mus_set_outletExpol(me, stencil, globBC, nLinks, minLevel, maxLevel)

Arguments

TypeIntentOptionalAttributesName
type(bc_outlet_type), intent(out), allocatable:: me(:)

setting type for bc

type(tem_stencilHeader_type), intent(in) :: stencil

for directions

type(glob_boundary_type), intent(in) :: globBC

number of scalars for number of elements in boundary and position in buffer

integer, intent(in) :: nLinks(minLevel:maxLevel)

for position of elements Number of links on each level

integer, intent(in) :: minLevel

minimum and maximum level

integer, intent(in) :: maxLevel

minimum and maximum level

public subroutine mus_fieldBC_cleanup(me, nBCs, minLevel, maxLevel)

This routines deallocates allocatables in field%bc boundary_type for dynamic load balancing

Arguments

TypeIntentOptionalAttributesName
type(boundary_type), intent(inout) :: me(nBCs)
integer, intent(in) :: nBCs
integer, intent(in) :: minLevel
integer, intent(in) :: maxLevel

private subroutine set_bouzidi_coeff(qVal, cIn, cOut, cNgh, cVel)

Set the coefficients of bouzidi linear interpolation boundary condition.

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: qVal
real(kind=rk), intent(out) :: cIn
real(kind=rk), intent(out) :: cOut
real(kind=rk), intent(out) :: cNgh
real(kind=rk), intent(out) :: cVel