mus_derivedQuantities_module2 Module

This module provides functions for calculating macroscopic quantities

! ************ ! !> Calculate the velocity in all 3 directions !! pure function getVelocity_forElemFromState_Force( state, elem, force, & & dtLB, stencil, varPos, nScalars ) result( vel ) ! --------------------------------------------------------------------------- real(kind=rk), intent(in) :: state(:) !< state array of a level type(tem_stencilHeader_type), intent(in) :: stencil real(kind=rk), intent(in) :: force(3) real(kind=rk), intent(in) :: dtLB integer, intent(in) :: elem integer, intent(in) :: varPos(:) !< varPos of current field variable integer, intent(in) :: nScalars !< number of scalars in global system real(kind=rk) :: vel(3) !< return value ! --------------------------------------------------------------------------- real(kind=rk) :: dens integer :: iDir ! --------------------------------------------------------------------------- integer :: nElems ! --------------------------------------------------------------------------- nElems = size( state ) / nScalars

vel = 0._rk
dens = 0._rk
do iDir = 1,stencil%QQ
  vel(:) = vel(:) + state( ( elem-1)* nscalars+ varpos(idir))  &
    &               * stencil%cxDirRK( :, iDir )           &
    ! Forcing Term
    & + dtLB*0.5_rk*force(:)

  dens = dens + state( ( elem-1)* nscalars+ varpos(idir))
enddo
vel = vel / dens

end function getVelocity_forElemFromState_Force ! ************ !

! ************ ! !> Calculate the shear rate tensor (strain rate) by acoustic scaling. !!

Rho  = 0._rk
Vel  = 0._rk
do iDir = 1, layout%fStencil%QQ
  fTmp( iDir ) = subset( iDir )
  vel(:) = vel(:) + fTmp( iDir ) * dble(layout%fStencil%cxDir(:,iDir))
end do
rho = sum( fTmp(:) )

if ( .not. incomp ) then
  vel  = vel / rho
  rhoD = rho
else
  rhoD = rho0
end if

! Get the second velocity moments of the source element's pdf
diagVal(:,:) = 0._rk
diagVal(1,1) = cs2*Rho
diagVal(2,2) = cs2*Rho
diagVal(3,3) = cs2*Rho
do jVal = 1, layout%fStencil%nDims
  do iVal = 1, layout%fStencil%nDims
    order = 0
    order(iVal) = order(iVal) +1
    order(jVal) = order(jVal) +1
    S(iVal,jVal) = get_moment( layout%fStencil%QQ, layout%fStencil%cxDir, order, fTmp)

    ! Now compute the S(1) from the second velocity moment
    S(iVal,jVal) = rhoD * Vel(iVal) * Vel(jVal) + diagVal(iVal,jVal) - S(iVal,jVal)
  end do
end do

S(:,:) = S(:,:) * omega * cs2inv * convPrePost(omega) * div1_2

end function getShearRateTensor_acoustic_forPdfSubset ! ************ !


Uses

  • module~~mus_derivedquantities_module2~~UsesGraph module~mus_derivedquantities_module2 mus_derivedQuantities_module2 module~tem_stencil_module tem_stencil_module module~mus_derivedquantities_module2->module~tem_stencil_module module~tem_param_module tem_param_module module~mus_derivedquantities_module2->module~tem_param_module module~mus_scheme_layout_module mus_scheme_layout_module module~mus_derivedquantities_module2->module~mus_scheme_layout_module module~mus_moments_module mus_moments_module module~mus_derivedquantities_module2->module~mus_moments_module module~tem_logging_module tem_logging_module module~mus_derivedquantities_module2->module~tem_logging_module module~tem_float_module tem_float_module module~mus_derivedquantities_module2->module~tem_float_module module~mus_graddata_module mus_gradData_module module~mus_derivedquantities_module2->module~mus_graddata_module module~env_module env_module module~mus_derivedquantities_module2->module~env_module module~mus_scheme_layout_module->module~tem_stencil_module module~mus_scheme_layout_module->module~tem_param_module module~mus_scheme_layout_module->module~tem_logging_module module~mus_scheme_layout_module->module~env_module module~tem_comm_env_module tem_comm_env_module module~mus_scheme_layout_module->module~tem_comm_env_module module~aot_table_module aot_table_module module~mus_scheme_layout_module->module~aot_table_module mpi mpi module~mus_scheme_layout_module->mpi module~aotus_module aotus_module module~mus_scheme_layout_module->module~aotus_module module~mus_moments_type_module mus_moments_type_module module~mus_scheme_layout_module->module~mus_moments_type_module module~tem_grow_array_module tem_grow_array_module module~mus_scheme_layout_module->module~tem_grow_array_module module~tem_tools_module tem_tools_module module~mus_scheme_layout_module->module~tem_tools_module module~tem_dyn_array_module tem_dyn_array_module module~mus_scheme_layout_module->module~tem_dyn_array_module module~aot_out_module aot_out_module module~mus_scheme_layout_module->module~aot_out_module module~tem_aux_module tem_aux_module module~mus_scheme_layout_module->module~tem_aux_module module~mus_moments_module->module~tem_logging_module module~mus_moments_module->module~env_module module~tem_debug_module tem_debug_module module~mus_moments_module->module~tem_debug_module module~mus_scheme_header_module mus_scheme_header_module module~mus_moments_module->module~mus_scheme_header_module module~mus_moments_module->module~mus_moments_type_module module~tem_matrix_module tem_matrix_module module~mus_moments_module->module~tem_matrix_module module~tem_math_module tem_math_module module~mus_moments_module->module~tem_math_module module~mus_moments_module->module~tem_aux_module module~mus_graddata_module->module~tem_stencil_module module~mus_graddata_module->module~tem_param_module module~mus_graddata_module->module~tem_logging_module module~mus_graddata_module->module~env_module module~mus_graddata_module->module~tem_debug_module module~tem_construction_module tem_construction_module module~mus_graddata_module->module~tem_construction_module module~mus_graddata_module->module~tem_aux_module module~mus_scheme_header_module->module~tem_logging_module module~mus_scheme_header_module->module~env_module module~mus_scheme_header_module->module~aot_table_module module~mus_scheme_header_module->module~aotus_module module~mus_scheme_header_module->module~tem_tools_module module~mus_scheme_header_module->module~aot_out_module module~mus_scheme_header_module->module~tem_aux_module module~mus_moments_type_module->module~env_module module~mus_moments_type_module->module~tem_matrix_module

Used by

  • module~~mus_derivedquantities_module2~~UsedByGraph module~mus_derivedquantities_module2 mus_derivedQuantities_module2 module~mus_interpolate_d3q27_module mus_interpolate_d3q27_module module~mus_interpolate_d3q27_module->module~mus_derivedquantities_module2 module~mus_interpolate_d2q9_module mus_interpolate_d2q9_module module~mus_interpolate_d2q9_module->module~mus_derivedquantities_module2 module~mus_bc_fluid_noneqexpol_module mus_bc_fluid_nonEqExpol_module module~mus_bc_fluid_noneqexpol_module->module~mus_derivedquantities_module2 module~mus_smagorinsky_module mus_Smagorinsky_module module~mus_smagorinsky_module->module~mus_derivedquantities_module2 module~mus_ibm_module mus_IBM_module module~mus_ibm_module->module~mus_derivedquantities_module2 module~mus_interpolate_debug_module mus_interpolate_debug_module module~mus_interpolate_debug_module->module~mus_derivedquantities_module2 module~mus_derquanisothermaceq_module mus_derQuanIsothermAcEq_module module~mus_derquanisothermaceq_module->module~mus_derivedquantities_module2 module~mus_derquanpoisson_module mus_derQuanPoisson_module module~mus_derquanpoisson_module->module~mus_derivedquantities_module2 module~mus_flow_module mus_flow_module module~mus_flow_module->module~mus_derivedquantities_module2 module~mus_derquanps_module mus_derQuanPS_module module~mus_derquanps_module->module~mus_derivedquantities_module2 module~mus_vreman_module mus_Vreman_module module~mus_vreman_module->module~mus_derivedquantities_module2 module~mus_derquan_module mus_derQuan_module module~mus_derquan_module->module~mus_derivedquantities_module2 module~mus_operation_var_module mus_operation_var_module module~mus_operation_var_module->module~mus_derivedquantities_module2 module~mus_derquanincomp_module mus_derQuanIncomp_module module~mus_derquanincomp_module->module~mus_derivedquantities_module2 module~mus_wale_module mus_WALE_module module~mus_wale_module->module~mus_derivedquantities_module2

Contents


Interfaces

public interface getShearStressTensor

  • private pure function getShearStressTensor_forElemFromState(state, neigh, elem, omega, layout, iField, varPos, nScalars) result(tau)

    Author
    Jiaxing Qi

    Calculate the viscous shear stress (exclude pressure) This function returns a one-dimensional array with 6 entries: tau(1:6) = [ Sxx, Syy, Szz, Sxy, Syz, Sxz ]

    Read more…

    Arguments

    TypeIntentOptionalAttributesName
    real(kind=rk), intent(in) :: state(:)
    integer, intent(in) :: neigh(:)
    integer, intent(in) :: elem
    real(kind=rk), intent(in) :: omega
    type(mus_scheme_layout_type), intent(in) :: layout
    integer, intent(in) :: iField
    integer, intent(in) :: varPos(:)
    integer, intent(in) :: nScalars

    Return Value real(kind=rk)(6)

  • private pure function getShearStressTensor_forPdfSubset(subset, omega, layout, varPos) result(tau)

    Author
    Manuel Hasert

    Calculate the viscous shear stress (exclude pressure) This function returns a one-dimensional array with 6 entries: tau(1:6) = [ Sxx, Syy, Szz, Sxy, Syz, Sxz ]

    Read more…

    Arguments

    TypeIntentOptionalAttributesName
    real(kind=rk), intent(in) :: subset(:)
    real(kind=rk), intent(in) :: omega
    type(mus_scheme_layout_type), intent(in) :: layout
    integer, intent(in) :: varPos(:)

    Return Value real(kind=rk)(6)

  • private pure function getShearStressTensorIncomp_forPdfSubset(subset, omega, layout, rho0) result(tau)

    Author
    Manuel Hasert

    Calculate the viscous shear stress (exclude pressure) This function returns a one-dimensional array with 6 entries: tau(1:6) = [ Sxx, Syy, Szz, Sxy, Syz, Sxz ]

    Read more…

    Arguments

    TypeIntentOptionalAttributesName
    real(kind=rk), intent(in) :: subset(:)
    real(kind=rk), intent(in) :: omega
    type(mus_scheme_layout_type), intent(in) :: layout
    real(kind=rk), intent(in) :: rho0

    Return Value real(kind=rk)(6)

  • private pure function getShearRateTensor_diffusive_forPdfSubset(f, omega, layout) result(S)

    Author
    Jiaxing Qi

    Calculate the strain rate tensor through 2nd moment. This function returns a two-dimensional 3 x 3 symmetirc array:

    Read more…

    Arguments

    TypeIntentOptionalAttributesName
    real(kind=rk), intent(in) :: f(layout%fStencil%QQ)

    pdf array ( post-collision value )

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

    relaxation parameter

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

    stencil layout

    Return Value real(kind=rk)(3,3)

    output array: strain rate tensor

public interface getShearRateTensor_acoustic

  • private pure function getShearRateTensor_acoustic_lbm(subset, omega, layout) result(S)

    This routine calculates shear rate tensor (i.e. strain rate tensor) by acoustic scaling (i.e. CE analysis)

    Arguments

    TypeIntentOptionalAttributesName
    real(kind=rk), intent(in) :: subset(:)
    real(kind=rk), intent(in) :: omega
    type(mus_scheme_layout_type), intent(in) :: layout

    Return Value real(kind=rk)(3,3)

  • private pure function getShearRateTensor_acoustic_incomp(subset, omega, layout, rho0) result(S)

    This routine calculates shear rate tensor (i.e. strain rate tensor) by acoustic scaling (i.e. CE analysis)

    Arguments

    TypeIntentOptionalAttributesName
    real(kind=rk), intent(in) :: subset(:)
    real(kind=rk), intent(in) :: omega
    type(mus_scheme_layout_type), intent(in) :: layout
    real(kind=rk), intent(in) :: rho0

    Return Value real(kind=rk)(3,3)

public interface getEquilibrium

  • private pure function getEquilibrium_forElemfromState(state, elem, layout, varPos, nScalars, neigh) result(equil)

    Calculate the equilibrium distribution function in all directions

    Read more…

    Arguments

    TypeIntentOptionalAttributesName
    real(kind=rk), intent(in) :: state(:)
    integer, intent(in) :: elem
    type(mus_scheme_layout_type), intent(in) :: layout
    integer, intent(in) :: varPos(:)
    integer, intent(in) :: nScalars
    integer, intent(in) :: neigh(:)

    Return Value real(kind=rk)(layout%fStencil%QQ)

  • public pure function getEqByDensVel(dens, vel, layout) result(equil)

    Calculate the equilibrium distribution function in all directions

    Read more…

    Arguments

    TypeIntentOptionalAttributesName
    real(kind=rk), intent(in) :: dens
    real(kind=rk), intent(in) :: vel(3)
    type(mus_scheme_layout_type), intent(in) :: layout

    Return Value real(kind=rk)(layout%fStencil%QQ)

  • private pure function getEquilibrium_forPdfSubset(subset, layout, varPos) result(equil)

    Calculate the equilibrium distribution function in all directions

    Read more…

    Arguments

    TypeIntentOptionalAttributesName
    real(kind=rk), intent(in) :: subset(:)
    type(mus_scheme_layout_type), intent(in) :: layout
    integer, intent(in) :: varPos(:)

    Return Value real(kind=rk)(layout%fStencil%QQ)

public interface getDensity

  • private pure function getDensity_forElemFromState(state, elem, stencil, varPos, nScalars) result(res)

    Calculate the density of a given element number with the given state vector (sum up all values)

    Arguments

    TypeIntentOptionalAttributesName
    real(kind=rk), intent(in) :: state(:)
    integer, intent(in) :: elem
    type(tem_stencilHeader_type), intent(in) :: stencil
    integer, intent(in) :: varPos(:)
    integer, intent(in) :: nScalars

    Return Value real(kind=rk)

  • private pure function getDensity_forPdfSubset(subset, stencil, varPos) result(res)

    Calculate the density of a given subset of pdfs vector (sum up all values)

    Arguments

    TypeIntentOptionalAttributesName
    real(kind=rk), intent(in) :: subset(:)
    type(tem_stencilHeader_type), intent(in) :: stencil
    integer, intent(in) :: varPos(:)

    Return Value real(kind=rk)

public interface getVelocity

  • private pure function getVelocity_forElemFromState_noForce(state, elem, stencil, varPos, nScalars) result(vel)

    Calculate the velocity in all 3 directions from the element indicated (elem), reading the pdf (state information) from the state array. state array includes all the pdfs of all elements. The access to the state array has to be done via the generic access macro IDX, as we want to access post-collision values.

    Arguments

    TypeIntentOptionalAttributesName
    real(kind=rk), intent(in) :: state(:)
    integer, intent(in) :: elem

    element index, for which to calc velocity

    type(tem_stencilHeader_type), intent(in) :: stencil
    integer, intent(in) :: varPos(:)
    integer, intent(in) :: nScalars

    Return Value real(kind=rk)(3)

  • private pure function getVelocity_forPdfSubset(subset, stencil, varPos) result(vel)

    Calculate the velocity in all 3 directions from a subset given, ordered according to the stencil

    Arguments

    TypeIntentOptionalAttributesName
    real(kind=rk), intent(in) :: subset(:)
    type(tem_stencilHeader_type), intent(in) :: stencil
    integer, intent(in) :: varPos(:)

    Return Value real(kind=rk)(3)

public interface getVelocity_incomp

  • private pure function getVelocity_forPdfSubset_incomp(subset, stencil, varPos) result(vel)

    Calculate the velocity in all 3 directions from a subset given, ordered according to the stencil

    Arguments

    TypeIntentOptionalAttributesName
    real(kind=rk), intent(in) :: subset(:)
    type(tem_stencilHeader_type), intent(in) :: stencil
    integer, intent(in) :: varPos(:)

    Return Value real(kind=rk)(3)


Functions

public pure function getEqDistr(iDir, rho, vel, layout) result(fEq)

Get the equilibrium distribution in the specified direction iDir

Arguments

TypeIntentOptionalAttributesName
integer, intent(in) :: iDir
real(kind=rk), intent(in) :: rho
real(kind=rk), intent(in) :: vel(3)
type(mus_scheme_layout_type), intent(in) :: layout

Return Value real(kind=rk)

public pure function getEqByDensVel(dens, vel, layout) result(equil)

Calculate the equilibrium distribution function in all directions

Read more…

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: dens
real(kind=rk), intent(in) :: vel(3)
type(mus_scheme_layout_type), intent(in) :: layout

Return Value real(kind=rk)(layout%fStencil%QQ)

public pure function getShearRate(strain) result(shear)

Author
Jiaxing Qi

Calculate the Shear Rate

Read more…

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: strain(:)

strain rate tensor: xx, yy, zz, xy, yz, zx

Return Value real(kind=rk)

public function getWSS(state, neigh, elem, omega, layout, iField, varPos, nScalars) result(wss)

Author
Jiaxing Qi

Calculate the wall shear stress (WSS)

Read more…

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: state(:)
integer, intent(in) :: neigh(:)
integer, intent(in) :: elem
real(kind=rk), intent(in) :: omega
type(mus_scheme_layout_type), intent(in) :: layout
integer, intent(in) :: iField
integer, intent(in) :: varPos(:)
integer, intent(in) :: nScalars

Return Value real(kind=rk)

public function getWSS2D(state, neigh, elem, omega, layout, iField, varPos, nScalars) result(wss)

Author
Jiaxing Qi

Calculate the wall shear stress (WSS) 2D

Read more…

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: state(:)
integer, intent(in) :: neigh(:)
integer, intent(in) :: elem
real(kind=rk), intent(in) :: omega
type(mus_scheme_layout_type), intent(in) :: layout
integer, intent(in) :: iField
integer, intent(in) :: varPos(:)
integer, intent(in) :: nScalars

Return Value real(kind=rk)

public function getNEq_acoustic(layout, omega, Sxx) result(nEq)

Author
Jiaxing Qi

Setting the non-equilibrium part based on the acoustic scaling

Read more…

Arguments

TypeIntentOptionalAttributesName
type(mus_scheme_layout_type), intent(in) :: layout
real(kind=rk), intent(in) :: omega
real(kind=rk), intent(in) :: Sxx(3,3)

strain rate tensor

Return Value real(kind=rk)(layout%fStencil%QQ)

public function getNEq_diffusive(layout, omega, Sxx) result(nEq)

Author
Jiaxing Qi

Calculate the non-equilibrium part of pdf from strain rate tensor based on the diffusive scaling

Read more…

Arguments

TypeIntentOptionalAttributesName
type(mus_scheme_layout_type), intent(in) :: layout
real(kind=rk), intent(in) :: omega
real(kind=rk), intent(in) :: Sxx(3,3)

Strain rate tensor. It is a symmetric 3x3 matrix

Return Value real(kind=rk)(layout%fStencil%QQ)

public function setThirdOrder_diffusive(layout, omega, x, dx, dt) result(f3)

Set third order for diffusive scaling (for push algorithm?!) Have a look at the initial condition file for TGV in Hasert/2013/multilevel_diffusive/ic Note that this is only for a 2d taylor green test case in the region x = y = 0:2pi with reference velocity u0 = 1

Arguments

TypeIntentOptionalAttributesName
type(mus_scheme_layout_type), intent(in) :: layout
real(kind=rk), intent(in) :: omega
real(kind=rk), intent(in) :: x(3)
real(kind=rk), intent(in) :: dx
real(kind=rk), intent(in) :: dt

Return Value real(kind=rk)(layout%fStencil%QQ)

public pure function convPrePost(omega) result(conv)

Author
Jiaxing Qi

Conversion factor betwen the pre- and post-collision quantity for the shear stress.

Read more…

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: omega

Return Value real(kind=rk)

public pure function getNonEqFac(omegaS, omegaT) result(fac)

Calculate the conversion factor for nonEq in difference levels

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: omegaS
real(kind=rk), intent(in) :: omegaT

Return Value real(kind=rk)

public pure function getNonEqFac_intp(omegaS, omegaT) result(fac)

Calculate the conversion factor to convert nonEq moments between fine and coarser.

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: omegaS
real(kind=rk), intent(in) :: omegaT

Return Value real(kind=rk)

public pure function getEquilibriumIncomp(dens, vel, layout, rho0) result(equil)

Calculate the equilibrium distribution function in all directions This is the incompressible formulation with reference density rho0

Read more…

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: dens
real(kind=rk), intent(in) :: vel(3)
type(mus_scheme_layout_type), intent(in) :: layout
real(kind=rk), intent(in) :: rho0

Return Value real(kind=rk)(layout%fStencil%QQ)

public function set_pdfDiffusive(layout, omega, rho0, mom) result(pdf)

Calculate the distribution function in all directions by using the fEq + fNeq

Read more…

Arguments

TypeIntentOptionalAttributesName
type(mus_scheme_layout_type), intent(in) :: layout
real(kind=rk), intent(in) :: omega
real(kind=rk), intent(in) :: rho0
real(kind=rk), intent(in) :: mom(layout%fStencil%QQ)

Return Value real(kind=rk)(layout%fStencil%QQ)

public function set_pdfAcoustic(layout, omega, rho0, mom, incompressible) result(pdf)

Calculate the distribution function in all directions by using the fEq + fNeq

Read more…

Arguments

TypeIntentOptionalAttributesName
type(mus_scheme_layout_type), intent(in) :: layout
real(kind=rk), intent(in) :: omega
real(kind=rk), intent(in) :: rho0
real(kind=rk), intent(in) :: mom(layout%fStencil%QQ)
logical, intent(in), optional :: incompressible

Return Value real(kind=rk)(layout%fStencil%QQ)

public pure function secondMom(cxcx, f, QQ) result(m)

Calculate second moments of some quantity where Q is number of discrete velocity.\n The output is 1 dimentional array which has 6 componenents.\n Specifically, This function is used by shear stress and strain rate.

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: cxcx(6,QQ)
real(kind=rk), intent(in) :: f(QQ)

quantity to which second moment is applied

integer, intent(in) :: QQ

Return Value real(kind=rk)(6)

public pure function getStrainFacDffs(omegaS, omegaT) result(fac)

Calculate the conversion factor for nonEq in difference levels

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: omegaS
real(kind=rk), intent(in) :: omegaT

Return Value real(kind=rk)

public pure function getGradU(auxField, gradData, velPos, nAuxScalars, nDims, nSolve, elemOffset) result(gradU)

This function computes gradient of velocity from gradient and veleocity data. Gradient is computed using central difference. if an element has an boundary then neighbor refers to current element then forward difference is used

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: auxField(:)

auxField

type(mus_gradData_type), intent(in) :: gradData

gradient data

integer, intent(in) :: velPos(3)

Position of velocity field in auxField

integer, intent(in) :: nAuxScalars

Number of scalars in auxField array

integer, intent(in) :: nDims

Dimensions

integer, intent(in) :: nSolve

Number of element to solve in this level

integer, intent(in) :: elemOffset

Offset for elements when computing chunkwise

Return Value real(kind=rk)(nDims,nDims,nSolve)

output: gradient of velocity

private pure function getDensity_forPdfSubset(subset, stencil, varPos) result(res)

Calculate the density of a given subset of pdfs vector (sum up all values)

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: subset(:)
type(tem_stencilHeader_type), intent(in) :: stencil
integer, intent(in) :: varPos(:)

Return Value real(kind=rk)

private pure function getDensity_forElemFromState(state, elem, stencil, varPos, nScalars) result(res)

Calculate the density of a given element number with the given state vector (sum up all values)

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: state(:)
integer, intent(in) :: elem
type(tem_stencilHeader_type), intent(in) :: stencil
integer, intent(in) :: varPos(:)
integer, intent(in) :: nScalars

Return Value real(kind=rk)

private pure function getVelocity_forPdfSubset(subset, stencil, varPos) result(vel)

Calculate the velocity in all 3 directions from a subset given, ordered according to the stencil

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: subset(:)
type(tem_stencilHeader_type), intent(in) :: stencil
integer, intent(in) :: varPos(:)

Return Value real(kind=rk)(3)

private pure function getVelocity_forPdfSubset_incomp(subset, stencil, varPos) result(vel)

Calculate the velocity in all 3 directions from a subset given, ordered according to the stencil

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: subset(:)
type(tem_stencilHeader_type), intent(in) :: stencil
integer, intent(in) :: varPos(:)

Return Value real(kind=rk)(3)

private pure function getVelocity_forElemFromState_noForce(state, elem, stencil, varPos, nScalars) result(vel)

Calculate the velocity in all 3 directions from the element indicated (elem), reading the pdf (state information) from the state array. state array includes all the pdfs of all elements. The access to the state array has to be done via the generic access macro IDX, as we want to access post-collision values.

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: state(:)
integer, intent(in) :: elem

element index, for which to calc velocity

type(tem_stencilHeader_type), intent(in) :: stencil
integer, intent(in) :: varPos(:)
integer, intent(in) :: nScalars

Return Value real(kind=rk)(3)

private pure function getEquilibrium_forPdfSubset(subset, layout, varPos) result(equil)

Calculate the equilibrium distribution function in all directions

Read more…

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: subset(:)
type(mus_scheme_layout_type), intent(in) :: layout
integer, intent(in) :: varPos(:)

Return Value real(kind=rk)(layout%fStencil%QQ)

private pure function getEquilibrium_forElemfromState(state, elem, layout, varPos, nScalars, neigh) result(equil)

Calculate the equilibrium distribution function in all directions

Read more…

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: state(:)
integer, intent(in) :: elem
type(mus_scheme_layout_type), intent(in) :: layout
integer, intent(in) :: varPos(:)
integer, intent(in) :: nScalars
integer, intent(in) :: neigh(:)

Return Value real(kind=rk)(layout%fStencil%QQ)

private pure function getShearStressTensorIncomp_forPdfSubset(subset, omega, layout, rho0) result(tau)

Author
Manuel Hasert

Calculate the viscous shear stress (exclude pressure) This function returns a one-dimensional array with 6 entries: tau(1:6) = [ Sxx, Syy, Szz, Sxy, Syz, Sxz ]

Read more…

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: subset(:)
real(kind=rk), intent(in) :: omega
type(mus_scheme_layout_type), intent(in) :: layout
real(kind=rk), intent(in) :: rho0

Return Value real(kind=rk)(6)

private pure function getShearStressTensor_forPdfSubset(subset, omega, layout, varPos) result(tau)

Author
Manuel Hasert

Calculate the viscous shear stress (exclude pressure) This function returns a one-dimensional array with 6 entries: tau(1:6) = [ Sxx, Syy, Szz, Sxy, Syz, Sxz ]

Read more…

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: subset(:)
real(kind=rk), intent(in) :: omega
type(mus_scheme_layout_type), intent(in) :: layout
integer, intent(in) :: varPos(:)

Return Value real(kind=rk)(6)

private pure function getShearStressTensor_forElemFromState(state, neigh, elem, omega, layout, iField, varPos, nScalars) result(tau)

Author
Jiaxing Qi

Calculate the viscous shear stress (exclude pressure) This function returns a one-dimensional array with 6 entries: tau(1:6) = [ Sxx, Syy, Szz, Sxy, Syz, Sxz ]

Read more…

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: state(:)
integer, intent(in) :: neigh(:)
integer, intent(in) :: elem
real(kind=rk), intent(in) :: omega
type(mus_scheme_layout_type), intent(in) :: layout
integer, intent(in) :: iField
integer, intent(in) :: varPos(:)
integer, intent(in) :: nScalars

Return Value real(kind=rk)(6)

private pure function getShearRateTensor_acoustic_incomp(subset, omega, layout, rho0) result(S)

This routine calculates shear rate tensor (i.e. strain rate tensor) by acoustic scaling (i.e. CE analysis)

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: subset(:)
real(kind=rk), intent(in) :: omega
type(mus_scheme_layout_type), intent(in) :: layout
real(kind=rk), intent(in) :: rho0

Return Value real(kind=rk)(3,3)

private pure function getShearRateTensor_acoustic_lbm(subset, omega, layout) result(S)

This routine calculates shear rate tensor (i.e. strain rate tensor) by acoustic scaling (i.e. CE analysis)

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: subset(:)
real(kind=rk), intent(in) :: omega
type(mus_scheme_layout_type), intent(in) :: layout

Return Value real(kind=rk)(3,3)

private pure function getShearRateTensor_diffusive_forPdfSubset(f, omega, layout) result(S)

Author
Jiaxing Qi

Calculate the strain rate tensor through 2nd moment. This function returns a two-dimensional 3 x 3 symmetirc array:

Read more…

Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: f(layout%fStencil%QQ)

pdf array ( post-collision value )

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

relaxation parameter

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

stencil layout

Return Value real(kind=rk)(3,3)

output array: strain rate tensor