turbulent_wall_noneq_expol_incomp Subroutine

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

BC routine for turbulent wall based on Guo's nonequilibrium extrapolation. The implementation is based on the following paper: Haussmann, M. et al. (2019) ‘Large-eddy simulation coupled with wall models for turbulent channel flows at high Reynolds numbers with a lattice Boltzmann method — Application to Coriolis mass flowmeter’, Computers & Mathematics with Applications. Elsevier Ltd, 78(10), pp. 3285–3302.

It uses wall model to compute velocity on the boundary node. All directions of PDF in the boundary elements are updated with Equilibrium plus non-equilibrium. Density is computed using Zho-He approach for straight walls. "On pressure and velocity boundary conditions for the lattice Boltzmann BGK model", Physics of Fluids 9, 1591-1598 (1997) https://doi.org/10.1063/1.869307

non-equilibrium are computed from PDF on neighbor and extrapolated to boundary. This routine is used for straight wall boundaries.

Usage


```lua boundary_condition = { { label = 'wall', kind = 'turbulent_wall_noneq_expol', wall_model = 'musker', nonlinear_solver = 'fixed_point' } }

This subroutine's interface must match the abstract interface definition boundaryRoutine in bc/mus_bc_header_module.f90 in order to be callable via fnct function pointer.

Arguments

TypeIntentOptionalAttributesName
class(boundary_type) :: me

global boundary type

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

Current state vector of iLevel

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

Treelm Mesh

integer, intent(in) :: nSize

size of state array ( in terms of elements )

integer, intent(in) :: iLevel

the level On which this boundary was invoked

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

global time information

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

connectivity array corresponding to state vector

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

stencil layout information

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

fluid parameters and properties

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


Calls

proc~~turbulent_wall_noneq_expol_incomp~~CallsGraph proc~turbulent_wall_noneq_expol_incomp turbulent_wall_noneq_expol_incomp proc~calcvelsw_unitsw_veltau_tvisc_incomp calcVelSW_unitSW_velTau_tVisc_incomp proc~turbulent_wall_noneq_expol_incomp->proc~calcvelsw_unitsw_veltau_tvisc_incomp

Contents


Variables

TypeVisibilityAttributesNameInitial
integer, private :: iElem
integer, private :: iDir
integer, private :: QQ
integer, private :: elemPos
integer, private :: normalInd_inv
real(kind=rk), private :: f_pre(layout%fStencil%QQ)
real(kind=rk), private :: f_neigh(layout%fStencil%QQ)
real(kind=rk), private :: dens_neigh
real(kind=rk), private :: velSW(3)
real(kind=rk), private :: vel_neigh(3)
real(kind=rk), private :: dens_bnd
real(kind=rk), private :: vel_normal
real(kind=rk), private :: normal(3)
real(kind=rk), private :: fEq
real(kind=rk), private :: fEq_bnd
real(kind=rk), private :: unitSW(3,globBC%nElems(iLevel))
real(kind=rk), private :: velSW_bnd(globBC%nElems(iLevel))