LBMComputeVelocity

Compute object for macroscopic velocity reconstruction.

LBMComputeVelocity computes macroscopic velocity from distribution functions. Forces can be added to velocity.

Overview

Computes velocity from the distribution f and density rho using the stencil direction vectors e_x,e_y,e_z:

u_x = (1/rho) sum_q f e_{x,q}, with analogous expressions for u_y and u_z in 2D/3D.

Optional force contributions are added as (forces)/(2 rho); constant body forces can be enabled on each axis via add_body_force = true and body_force_{x,y,z} constants.

Example Input File Syntax

[TensorComputes]
  [Solve]
    [rho]
      type = LBMComputeDensity
      buffer = rho
      f = f
    []
    [vel]
      type = LBMComputeVelocity
      buffer = u
      f = f
      rho = rho
      enable_forces = true
      forces = F
      add_body_force = true
      body_force_x = '1e-6'
      body_force_y = '0.0'
      body_force_z = '0.0'
    []
  []
[]
!listing-end

Input Parameters

bufferThe buffer this compute is writing to
C++ Type:std::string
Controllable:No
Description:The buffer this compute is writing to
fDistribution function
C++ Type:std::string
Controllable:No
Description:Distribution function
rhoDensity
C++ Type:std::string
Controllable:No
Description:Density

Required Parameters

add_body_forceFalseWhether to enable forces or no
Default:False
C++ Type:bool
Controllable:No
Description:Whether to enable forces or no
body_force_x0.0Body force to be added in x-dir
Default:0.0
C++ Type:std::string
Controllable:No
Description:Body force to be added in x-dir
body_force_y0.0Body force to be added in y-dir
Default:0.0
C++ Type:std::string
Controllable:No
Description:Body force to be added in y-dir
body_force_z0.0Body force to be added in z-dir
Default:0.0
C++ Type:std::string
Controllable:No
Description:Body force to be added in z-dir
enable_forcesFalseWhether to enable forces or no
Default:False
C++ Type:bool
Controllable:No
Description:Whether to enable forces or no
forcesforcesForce tensor
Default:forces
C++ Type:std::string
Controllable:No
Description:Force tensor

Optional Parameters

control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Controllable:No
Description:Adds user-defined labels for accessing object parameters via control logic.
enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Controllable:No
Description:Set the enabled status of the MooseObject.

Advanced Parameters

Input Files

(examples/lbm/Phase-field/static_bubble_3d.i)
(test/tests/lbm/vertical_density_bcs.i)
(examples/lbm/natural_convection_3D/convection.i)
(test/tests/lbm/convective_outflow_2d_right.i)
(test/tests/lbm/nee_3d_box_first.i)
(test/tests/lbm/convective_outflow_2d_left.i)
(test/tests/lbm/mixed_bcs_d3q19_reverse.i)
(examples/lbm/rarefied_gas/channel.i)
(test/tests/lbm/convective_outflow_3d_left.i)
(test/tests/lbm/convective_outflow_3d_front.i)
(test/tests/lbm/channel3D.i)
(test/tests/lbm/convective_outflow_3d_bottom.i)
(test/tests/lbm/nee_3d_box_second.i)
(test/tests/lbm/mixed_bcs_d3q19.i)
(test/tests/lbm/horizontal_mixed_bcs_d2q9_reverse.i)
(examples/lbm/Rayleigh-Benard/solve.i)
(test/tests/lbm/horizontal_mixed_bcs_d2q9.i)
(test/tests/lbm/nee_2d_channel.i)
(test/tests/lbm/convective_outflow_2d_bottom.i)
(test/tests/lbm/convective_outflow_2d_top.i)
(examples/lbm/Phase-field/droplet_impact.i)
(examples/lbm/Phase-field/spinodal_decomposition.i)
(test/tests/lbm/channel2D.i)
(examples/lbm/Karman-vortex/cylinder.i)
(test/tests/lbm/convective_outflow_3d_right.i)
(test/tests/lbm/mixed_bcs_d3q27.i)
(test/tests/lbm/rarefied_gas.i)
(test/tests/lbm/obstacle.i)
(test/tests/lbm/nee_2d_all_walls.i)
(test/tests/lbm/isotropic_stencil_mrt.i)
(test/tests/lbm/vertical_velocity_bcs.i)
(examples/lbm/Phase-field/static_bubble.i)
(test/tests/lbm/phase.i)
(test/tests/lbm/convective_outflow_3d_top.i)
(test/tests/lbm/smagorinsky_mrt.i)
(examples/lbm/Formula1-aerodynamics/f1.i)
(test/tests/lbm/mixed_bcs_d3q27_reverse.i)
(test/tests/lbm/phase_3D.i)
(test/tests/lbm/convective_outflow_3d_back.i)
(examples/lbm/Phase-field/layered_poiseuille.i)
(examples/lbm/Pebbles/pebbles.i)

(examples/lbm/Phase-field/static_bubble_3d.i)

#
# 3D Static bubble test case
# Based on the 2D static bubble (PHYSICAL REVIEW E 97, 033309 - Section III.A)
# Small domain to verify 3D gradient/laplacian correctness
#

# Domain (small for quick testing)
Nx = 40
Ny = 40
Nz = 40

# Bubble parameters
Cx = '${Nx}/2.0'
Cy = '${Ny}/2.0'
Cz = '${Nz}/2.0'
R = 12

# Fluid properties
rho_l = 1000.0
rho_g = 1.0
nu_l = 0.1
nu_g = 0.1
sigma = 0.001

# Phase field parameters
tau_h = 0.7
D = 4

[Domain]
  dim = 3
  nx = '${Nx}'
  ny = '${Ny}'
  nz = '${Nz}'
  xmax = '${Nx}'
  ymax = '${Ny}'
  zmax = '${Nz}'
  device_names = 'cpu'
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [d3q27]
    type = LBMD3Q27
  []
[]

[TensorBuffers]
  # Macroscopic phase field variables
  [phi]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [grad_phi]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [laplacian_phi]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [mu]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [forces]
    type = LBMTensorBuffer
    buffer_type = mv
  []

  # Macroscopic hydrodynamic variables
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [pressure]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [rho]
    type = LBMTensorBuffer
    buffer_type = ms
  []

  # LBM phase field variables
  [h]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [h_post_collision]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [h_eq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [relaxation_tensor]
    type = LBMTensorBuffer
    buffer_type = ms
  []

  # LBM hydrodynamic variables
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_post_collision]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_eq]
    type = LBMTensorBuffer
    buffer_type = df
  []
[]

[TensorComputes/Initialize]
  [phi_init]
    type = ParsedCompute
    buffer = phi
    expression = '0.5 + 0.5 * tanh(2*(R - sqrt((x - Cx)^2 + (y - Cy)^2 + (z - Cz)^2)) / D)'
    constant_names = 'Cx Cy Cz R D'
    constant_expressions = '${Cx} ${Cy} ${Cz} ${R} ${D}'
    extra_symbols = true
  []
  [grad_phi_init]
    type = LBMIsotropicGradient
    buffer = grad_phi
    scalar_field = phi
  []
  [rho_init]
    type = ParsedCompute
    buffer = rho
    extra_symbols = true
    expression = 'phi*(rho_l - rho_g) + rho_g'
    constant_names = 'rho_l rho_g'
    constant_expressions = '${rho_l} ${rho_g}'
    inputs = phi
  []
  [pressure_init]
    type = LBMConstantTensor
    buffer = pressure
    constants = 0.3
  []
  [h_eq_init]
    type = LBMPhaseEquilibrium
    buffer = h_eq
    phi = phi
    velocity = velocity
  []
  [h_post_collision_init]
    type = LBMPhaseEquilibrium
    buffer = h_post_collision
    phi = phi
    velocity = velocity
  []
  [h_init]
    type = LBMPhaseEquilibrium
    buffer = h
    phi = phi
    velocity = velocity
  []
  [f_eq_init]
    type = LBMPressureCorrectedEquilibrium
    buffer = f_eq
    rho = rho
    velocity = velocity
    pressure = pressure
  []
  [f_post_collision_init]
    type = ParsedCompute
    buffer = f_post_collision
    expression = 'f_eq'
    inputs = f_eq
  []
  [f_init]
    type = ParsedCompute
    buffer = f
    expression = 'f_eq'
    inputs = f_eq
  []
[]

[TensorComputes/Solve]
  # Phase Field
  [compute_phi]
    type = LBMComputeDensity
    buffer = phi
    f = h
  []
  [grad_phi]
    type = LBMIsotropicGradient
    buffer = grad_phi
    scalar_field = phi
  []
  [laplacian_phi]
    type = LBMIsotropicLaplacian
    buffer = laplacian_phi
    scalar_field = phi
  []
  [potential]
    type = LBMComputeChemicalPotential
    buffer = mu
    phi = phi
    laplacian_phi = laplacian_phi
    thickness = D
    sigma = sigma
  []
  [forces]
    type = LBMComputeSurfaceForces
    buffer = forces
    chemical_potential = mu
    grad_phi = grad_phi
  []
  # Hydrodynamics
  [density]
    type = ParsedCompute
    buffer = rho
    extra_symbols = true
    expression = 'phi*(rho_l - rho_g) + rho_g'
    constant_names = 'rho_l rho_g'
    constant_expressions = '${rho_l} ${rho_g}'
    inputs = phi
  []
  [velocity]
    type = LBMComputeVelocity
    buffer = velocity
    f = f
    rho = rho
    enable_forces = true
    forces = forces
  []
  [h_eq]
    type = LBMPhaseEquilibrium
    buffer = h_eq
    phi = phi
    velocity = velocity
  []
  [phase_collision]
    type = LBMBGKCollision
    buffer = h_post_collision
    f = h
    feq = h_eq
    tau0 = tau_h
  []
  [apply_forces_phase]
    type = LBMAllenCahnSource
    buffer = h_post_collision
    phi = phi
    velocity = velocity
    grad_phi = grad_phi
    tau = tau_h
    thickness = D
  []
  [relaxation_tensor]
    type = ParsedCompute
    buffer = relaxation_tensor
    extra_symbols = true
    expression = '(phi*(nu_l - nu_g) + nu_g)/cs2+0.5'
    constant_names = 'nu_l nu_g cs2'
    constant_expressions = '${nu_l} ${nu_g} 0.3333'
    inputs = phi
  []
  [pressure]
    type = LBMPhaseFieldPressure
    buffer = pressure
    f = f
    velocity = velocity
    grad_phi = grad_phi
    rho = rho
    rho_l = '${rho_l}'
    rho_g = '${rho_g}'
  []
  [f_eq]
    type = LBMPressureCorrectedEquilibrium
    buffer = f_eq
    rho = rho
    velocity = velocity
    pressure = pressure
  []
  [collision]
    type = LBMBGKCollision
    buffer = f_post_collision
    f = f
    feq = f_eq
    tau0 = 1.0
    is_dynamic_relaxation = true
    tau_tensor = relaxation_tensor
  []
  [apply_forces_hydro]
    type = LBMForceDistribution
    buffer = f_post_collision
    grad_phi = grad_phi
    velocity = velocity
    forces = forces
    tau_tensor = relaxation_tensor
    tau = 1.0
    rho_l = '${rho_l}'
    rho_g = '${rho_g}'
    is_dynamic_relaxation = true
  []
  [residual]
    type = LBMComputeResidual
    buffer = phi
    speed = phi
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = 'h f'
  f_old = 'h_post_collision f_post_collision'
  root_compute = residual
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = rho
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = rho
    value_type = MAX
  []
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 10
  print_debug_output = true
  scalar_constant_names = 'tau_h D sigma'
  scalar_constant_values = '${tau_h} ${D} ${sigma}'
[]

[Executioner]
  type = Transient
  num_steps = 50
[]

[TensorOutputs]
  [xdmf]
    type = XDMFTensorOutput
    buffer = 'phi rho velocity forces'
    output_mode = 'Cell Cell Cell Cell'
    enable_hdf5 = true
    # transpose = false
  []
[]

(test/tests/lbm/vertical_density_bcs.i)

[Domain]
  dim = 2
  nx = 10
  ny = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_bounce_back]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type=LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type=LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0001 0.0005'
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_f_bb]
      type = LBMEquilibrium
      buffer = f_bounce_back
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
    []
  []
  [Boundary]
    [left]
      type = LBMBounceBack
      buffer = f
      f_old = f_bounce_back
      boundary = left
    []
    [right]
      type = LBMBounceBack
      buffer = f
      f_old = f_bounce_back
      boundary = right
    []
    [top]
      type = LBMFixedZerothOrderBC
      buffer = f
      f = f
      value = 1.1
      boundary = top
    []
    [bottom]
      type = LBMFixedZerothOrderBC
      buffer = f
      f = f
      value = 1.00000
      boundary = bottom
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = f_bounce_back
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 2
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [densty_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = vertical_density_bcs
  csv = true
[]

(examples/lbm/natural_convection_3D/convection.i)

[Domain]
  dim = 3
  nx = 270
  ny = 270
  nz = 405
  xmax = 270
  ymax = 270
  zmax = 405
  device_names='cpu'
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
  floating_precision=SINGLE
[]

[Stencil]
  [d3q19]
    type = LBMD3Q19
  []
[]

[TensorBuffers]
  # Simulation binary media
  [binary_media]
    type = LBMTensorBuffer
    file = binary_media.h5
    buffer_type = ms
    is_integer = true
  []

  # Density distribution functions
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  # Temperature distribution functions
  [g]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [geq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [gpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  # Fluid macroscopic variables: density and velocity
  [density]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [speed]
    type=LBMTensorBuffer
    buffer_type = ms
  []
  # Temperature macroscpic variables: temperature and 'velocity'
  [T]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  # Forces
  [F]
    type = LBMTensorBuffer
    buffer_type = mv
  []
[]

[TensorComputes]
  [Initialize]
    [density]
      type = LBMConstantTensor
      buffer = density
      constants = 'rho0'
    []
    [velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0 0.0'
    []
    [temperature]
      type = LBMConstantTensor
      buffer = T
      constants = T_C
    []
    [equilibrium_fluid]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [equilibrium_fluid_total]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [equilibrium_fluid_pc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
    [equilibrium_temperature]
      type = LBMEquilibrium
      buffer = geq
      bulk = T
      velocity = velocity
    []
    [equilibrium_temperature_total]
      type = LBMEquilibrium
      buffer = g
      bulk = T
      velocity = velocity
    []
    [equilibrium_temperature_pc]
      type = LBMEquilibrium
      buffer = gpc
      bulk = T
      velocity = velocity
    []
  []

  #### Compute ####
  [Solve]

    # For temperature
    [Temperature]
      type = LBMComputeDensity
      buffer = T
      f = g
    []

    # For fluid
    [Fluid_density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []

    [Fluid_velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
      forces = F
      enable_forces = true
    []

    # For temperature
    [Equilibrium_temperature]
      type = LBMEquilibrium
      buffer = geq
      bulk = T
      velocity = velocity
    []

    [Collision_temperature]
      type = LBMBGKCollision
      buffer = gpc
      f = g
      feq = geq
      tau0 = tau_T
    []

    # For fluid
    [Compute_forces]
      type = LBMComputeForces
      buffer = F
      rho0 = 'rho0'
      temperature = T
      T0 = 1.00
      enable_buoyancy = true
      gravity = g
      gravity_direction = 2
    []

    [Equilibrium_fluid]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []

    [Collision_fluid]
      type = LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = tau_f
    []

    [Apply_forces]
      type = LBMApplyForces
      buffer = fpc
      velocity = velocity
      rho = density
      forces = F
      tau0 = tau_f
    []

    [Residual]
      type = LBMComputeResidual
      speed = T
      # TODO this buffer is redundant, but avoids 'missing parameter' error
      buffer = T
    []
  []

  #### Boundary ####
  [Boundary]
    ##### for fluid
    [wall]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = wall
    []
    ##### for temperature
    [heat_source]
      type = LBMNeumannBC
      buffer = g
      f_old = gpc
      feq=geq
      velocity = velocity
      rho = T
      gradient = 0.001
      region_id = 3
      boundary = regional
    []
    [walls]
      type = LBMDirichletBC
      buffer = g
      f_old = gpc
      feq=geq
      velocity = velocity
      rho = T
      value = 1.0
      region_id = 2
      boundary = regional
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = 'f g'
  f_old = 'fpc gpc'
  root_compute = Residual
[]

[Problem]
  type = LatticeBoltzmannProblem
  binary_media = binary_media
  scalar_constant_names = 'rho0 T_C T_H tau_f tau_T g'
  scalar_constant_values = '1.0 1.0 1.05 0.55 0.55 0.01'
  substeps = 100
  print_debug_output = true
[]

[Executioner]
  type = Transient
  num_steps = 1000
[]

[TensorOutputs]
  [xdmf2]
    type = XDMFTensorOutput
    buffer = 'T density velocity'
    output_mode = 'Cell Cell Cell'
    enable_hdf5 = true
  []
[]

(test/tests/lbm/convective_outflow_2d_right.i)

[Domain]
  dim = 2
  nx = 10
  ny = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type = LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0'
    []
    [initial_feq]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_fpc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type = LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
    []
  []
  [Boundary]
    [left]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.005
      uy = 0.0
      order = first
      boundary = left
    []
    [right]
      type = LBMConvectiveOutflow
      buffer = f
      f_old = f
      convection_velocity = auto
      boundary = right
    []
    [top]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = top
    []
    [bottom]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = bottom
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = convective_outflow_2d_right
  csv = true
[]

(test/tests/lbm/nee_3d_box_first.i)

[Domain]
  dim = 3
  nx = 8
  ny = 8
  nz = 8
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [d3q19]
    type = LBMD3Q19
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type = LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0 0.0'
    []
    [initial_feq]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_fpc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type = LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
      add_body_force = true
      body_force_x = 0.0001
    []
  []
  [Boundary]
    [left]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = 0.0
      uz = 0.0
      order = first
      boundary = left
    []
    [right]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = 0.0
      uz = 0.0
      order = first
      boundary = right
    []
    [top]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = 0.0
      uz = 0.0
      order = first
      boundary = top
    []
    [bottom]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = 0.0
      uz = 0.0
      order = first
      boundary = bottom
    []
    [front]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = 0.0
      uz = 0.0
      order = first
      boundary = front
    []
    [back]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = 0.0
      uz = 0.0
      order = first
      boundary = back
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = nee_3d_box_first
  csv = true
[]

(test/tests/lbm/convective_outflow_2d_left.i)

[Domain]
  dim = 2
  nx = 10
  ny = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type = LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0'
    []
    [initial_feq]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_fpc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type = LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
    []
  []
  [Boundary]
    [right]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = -0.005
      uy = 0.0
      order = first
      boundary = right
    []
    [left]
      type = LBMConvectiveOutflow
      buffer = f
      f_old = f
      convection_velocity = auto
      boundary = left
    []
    [top]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = top
    []
    [bottom]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = bottom
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = convective_outflow_2d_left
  csv = true
[]

(test/tests/lbm/mixed_bcs_d3q19_reverse.i)

[Domain]
  dim = 3
  nx = 10
  ny = 10
  nz = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [d3q19]
    type = LBMD3Q19
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_bounce_back]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type=LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type=LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0001 0.0005 0.0'
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_f_bb]
      type = LBMEquilibrium
      buffer = f_bounce_back
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
      add_body_force = true
      body_force_x = 0.0001
    []
  []
  [Boundary]
    [right]
      type = LBMFixedZerothOrderBC
      buffer = f
      f = f
      value = 1.11
      boundary = right
    []
    [left]
      type = LBMFixedFirstOrderBC
      buffer = f
      f = f
      value = 0.0001
      boundary = left
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = f_bounce_back
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [densty_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = mixed_bcs_d3q19_reverse
  csv = true
[]

(examples/lbm/rarefied_gas/channel.i)

[Domain]
  dim = 2
  nx = 102
  ny = 102
  xmax = 102
  ymax = 102
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type=LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type=LBMTensorBuffer
    buffer_type = ms
  []
  [speed]
    type=LBMTensorBuffer
    buffer_type = ms
  []
  [domain]
    type=LBMTensorBuffer
    file = binary_media.h5
    buffer_type = ms
    is_integer = true
  []
  [local_pore]
    type=LBMTensorBuffer
    file = local_pore_size.h5
    buffer_type = ms
    is_integer = false
  []
  [Kn]
    type=LBMTensorBuffer
    file = Kn.h5
    buffer_type = ms
    is_integer = false
  []
  [relaxation_matrix]
    type = LBMTensorBuffer
    buffer_type = df
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 0.2355545440759889
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0'
    []
    [initial_equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_distribution]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_distribution_pc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
    [relaxation_matrix_init]
      type = LBMComputeEffectiveRelaxation
      buffer = relaxation_matrix
      local_pore_size = local_pore
      local_Knudsen_number = Kn
      mfp = 7.904614716131531e-10
      dx = 0.122e-9
      A2 = 0.82
    []
  []
  [Solve]
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
      add_body_force = true
      body_force_x = 1.0e-8
    []
    [equilibrium]
      type=LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type = LBMMRTCollision
      buffer = fpc
      f = f
      feq = feq
      is_dynamic_relaxation = true
      local_relaxation_matrix = relaxation_matrix
      projection = true
    []
    [speed]
      type = LBMComputeVelocityMagnitude
      buffer = speed
      velocity = velocity
    []
    [residual]
      type = LBMComputeResidual
      buffer = speed
      speed = speed
    []
  []
  [Boundary]
    [wall]
      type = LBMSpecularReflectionBoundary
      buffer = f
      f_old = fpc
      local_Knudsen_number = Kn
      boundary = wall
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
  root_compute = residual
[]

[Postprocessors]
  [max_u]
    type = TensorExtremeValuePostprocessor
    buffer = speed
    value_type = MAX
    execute_on = 'INITIAL TIMESTEP_END'
  []
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 2000
  print_debug_output = true
  binary_media = domain
[]

[Executioner]
  type = Transient
  num_steps = 2
[]

[TensorOutputs]
  [xdmf2]
    type = XDMFTensorOutput
    buffer = 'density velocity speed'
    output_mode = 'Cell Cell Cell'
    enable_hdf5 = true
    transpose=false
  []
[]

(test/tests/lbm/convective_outflow_3d_left.i)

[Domain]
  dim = 3
  nx = 8
  ny = 8
  nz = 8
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [d3q19]
    type = LBMD3Q19
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type = LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0 0.0'
    []
    [initial_feq]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_fpc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type = LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
    []
  []
  [Boundary]
    [right]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = -0.005
      uy = 0.0
      uz = 0.0
      order = first
      boundary = right
    []
    [left]
      type = LBMConvectiveOutflow
      buffer = f
      f_old = f
      convection_velocity = Uc
      boundary = left
    []
    [top]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = top
    []
    [bottom]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = bottom
    []
    [front]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = front
    []
    [back]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = back
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  scalar_constant_names = 'Uc'
  scalar_constant_values = '0.005'
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = convective_outflow_3d_left
  csv = true
[]

(test/tests/lbm/convective_outflow_3d_front.i)

[Domain]
  dim = 3
  nx = 8
  ny = 8
  nz = 8
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [d3q19]
    type = LBMD3Q19
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type = LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0 0.0'
    []
    [initial_feq]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_fpc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type = LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
    []
  []
  [Boundary]
    [back]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = 0.0
      uz = -0.005
      order = first
      boundary = back
    []
    [front]
      type = LBMConvectiveOutflow
      buffer = f
      f_old = f
      convection_velocity = Uc
      boundary = front
    []
    [left]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = left
    []
    [right]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = right
    []
    [top]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = top
    []
    [bottom]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = bottom
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  scalar_constant_names = 'Uc'
  scalar_constant_values = '0.005'
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = convective_outflow_3d_front
  csv = true
[]

(test/tests/lbm/channel3D.i)

[Domain]
  dim = 3
  nx = 10
  ny = 10
  nz = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [d3q19]
    type = LBMD3Q19
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type=LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type=LBMTensorBuffer
    buffer_type = ms
  []
  [speed]
    type=LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0 0.0'
    []
    [initial_equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_distribution]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_distribution_pc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type=LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type=LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
      add_body_force = true
      body_force_x = 0.0001
    []
    [speed]
      type = LBMComputeVelocityMagnitude
      buffer = speed
      velocity = velocity
    []
    [residual]
      type = LBMComputeResidual
      buffer = speed
      speed = speed
    []
  []
  [Boundary]
    [top]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = top
    []
    [bottom]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = bottom
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 10
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [speed_min]
    type = TensorExtremeValuePostprocessor
    buffer = speed
    value_type = MIN
  []
  [speed_max]
    type = TensorExtremeValuePostprocessor
    buffer = speed
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [densty_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 11
[]

[Outputs]
  file_base = channel3D
  csv = true
[]

(test/tests/lbm/convective_outflow_3d_bottom.i)

[Domain]
  dim = 3
  nx = 8
  ny = 8
  nz = 8
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [d3q19]
    type = LBMD3Q19
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type = LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0 0.0'
    []
    [initial_feq]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_fpc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type = LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
    []
  []
  [Boundary]
    [top]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = -0.005
      uz = 0.0
      order = first
      boundary = top
    []
    [bottom]
      type = LBMConvectiveOutflow
      buffer = f
      f_old = f
      convection_velocity = Uc
      boundary = bottom
    []
    [left]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = left
    []
    [right]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = right
    []
    [front]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = front
    []
    [back]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = back
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  scalar_constant_names = 'Uc'
  scalar_constant_values = '0.005'
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = convective_outflow_3d_bottom
  csv = true
[]

(test/tests/lbm/nee_3d_box_second.i)

[Domain]
  dim = 3
  nx = 10
  ny = 10
  nz = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [d3q19]
    type = LBMD3Q19
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type = LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0 0.0'
    []
    [initial_feq]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_fpc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type = LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
    []
  []
  [Boundary]
    [left]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.005
      uy = 0.0
      uz = 0.0
      order = second
      boundary = left
    []
    [right]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = density
      prescribed_rho = 1.0
      order = second
      boundary = right
    []
    [top]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = 0.0
      uz = 0.0
      order = second
      boundary = top
    []
    [bottom]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = 0.0
      uz = 0.0
      order = second
      boundary = bottom
    []
    [front]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = density
      prescribed_rho = 1.0
      order = second
      boundary = front
    []
    [back]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = density
      prescribed_rho = 1.0
      order = second
      boundary = back
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = nee_3d_box_second
  csv = true
[]

(test/tests/lbm/mixed_bcs_d3q19.i)

[Domain]
  dim = 3
  nx = 10
  ny = 10
  nz = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [d3q19]
    type = LBMD3Q19
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_bounce_back]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type=LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type=LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0001 0.0005 0.0'
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_f_bb]
      type = LBMEquilibrium
      buffer = f_bounce_back
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
      add_body_force = true
      body_force_x = 0.0001
    []
  []
  [Boundary]
    [left]
      type = LBMFixedZerothOrderBC
      buffer = f
      f = f
      value = 1.11
      boundary = left
    []
    [right]
      type = LBMFixedFirstOrderBC
      buffer = f
      f = f
      value = 0.0001
      boundary = right
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = f_bounce_back
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [densty_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 1
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = mixed_bcs_d3q19
  csv = true
[]

(test/tests/lbm/horizontal_mixed_bcs_d2q9_reverse.i)

[Domain]
  dim = 2
  nx = 10
  ny = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_bounce_back]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type=LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type=LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0001 0.0005'
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_f_bb]
      type = LBMEquilibrium
      buffer = f_bounce_back
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
      add_body_force = true
      body_force_x = 0.0001
    []
  []
  [Boundary]
    [right]
      type = LBMFixedZerothOrderBC
      buffer = f
      f = f
      value = 1.1
      boundary = right
    []
    [left]
      type = LBMMicroscopicZeroGradientBC
      buffer = f
      boundary = left
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = f_bounce_back
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [densty_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = horizontal_mixed_bcs_d2q9_reverse
  csv = true
[]

(examples/lbm/Rayleigh-Benard/solve.i)

#### Compute ####
[Solve]

  # For temperature
  [Temperature]
    type = LBMComputeDensity
    buffer = T
    f = g
  []

  # For fluid
  [Fluid_density]
    type = LBMComputeDensity
    buffer = density
    f = f
  []

  [Fluid_velocity]
    type = LBMComputeVelocity
    buffer = velocity
    f = f
    rho = density
    forces = F
    enable_forces = true
  []

  # For temperature
  [Equilibrium_temperature]
    type = LBMEquilibrium
    buffer = geq
    bulk = T
    velocity = velocity
  []

  [Collision_temperature]
    type = LBMBGKCollision
    buffer = gpc
    f = g
    feq = geq
    tau0 = tau_T
  []

  # For fluid
  [Compute_forces]
    type = LBMComputeForces
    buffer = F
    rho0 = 'rho0'
    temperature = T
    T0 = T_0
    enable_buoyancy = true
    gravity = g
  []

  [Equilibrium_fluid]
    type = LBMEquilibrium
    buffer = feq
    bulk = density
    velocity = velocity
  []

  [Collision_fluid]
    type = LBMBGKCollision
    buffer = fpc
    f = f
    feq = feq
    tau0 = tau_f
  []

  [Apply_forces]
    type = LBMApplyForces
    buffer = fpc
    velocity = velocity
    rho = density
    forces = F
    tau0 = tau_f
  []

  [Residual]
    type = LBMComputeResidual
    speed = T
    # TODO this buffer is redundant, but avoids 'missing parameter' error
    buffer = T
  []
[]

(test/tests/lbm/horizontal_mixed_bcs_d2q9.i)

[Domain]
  dim = 2
  nx = 10
  ny = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_bounce_back]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type=LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type=LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0001 0.0005'
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_f_bb]
      type = LBMEquilibrium
      buffer = f_bounce_back
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
      add_body_force = true
      body_force_x = 0.0001
    []
  []
  [Boundary]
    [left]
      type = LBMFixedZerothOrderBC
      buffer = f
      f = f
      value = 1.1
      boundary = left
    []
    [right]
      type = LBMFixedFirstOrderBC
      buffer = f
      f = f
      value = 0.0001
      boundary = right
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = f_bounce_back
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [densty_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = horizontal_mixed_bcs_d2q9
  csv = true
[]

(test/tests/lbm/nee_2d_channel.i)

[Domain]
  dim = 2
  nx = 10
  ny = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type = LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0'
    []
    [initial_feq]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_fpc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type = LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
    []
  []
  [Boundary]
    [left]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.005
      uy = 0.0
      order = first
      boundary = left
    []
    [right]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = density
      prescribed_rho = 1.0
      order = first
      boundary = right
    []
    [top]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = top
    []
    [bottom]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = bottom
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = nee_2d_channel
  csv = true
[]

(test/tests/lbm/convective_outflow_2d_bottom.i)

[Domain]
  dim = 2
  nx = 10
  ny = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type = LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0'
    []
    [initial_feq]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_fpc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type = LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
    []
  []
  [Boundary]
    [top]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = -0.005
      order = first
      boundary = top
    []
    [bottom]
      type = LBMConvectiveOutflow
      buffer = f
      f_old = f
      convection_velocity = auto
      boundary = bottom
    []
    [left]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = left
    []
    [right]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = right
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = convective_outflow_2d_bottom
  csv = true
[]

(test/tests/lbm/convective_outflow_2d_top.i)

[Domain]
  dim = 2
  nx = 10
  ny = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type = LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0'
    []
    [initial_feq]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_fpc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type = LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
    []
  []
  [Boundary]
    [bottom]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = 0.005
      order = first
      boundary = bottom
    []
    [top]
      type = LBMConvectiveOutflow
      buffer = f
      f_old = f
      convection_velocity = auto
      boundary = top
    []
    [left]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = left
    []
    [right]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = right
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = convective_outflow_2d_top
  csv = true
[]

(examples/lbm/Phase-field/droplet_impact.i)

#
# Droplet Impact on a Thin Liquid Film
# PHYSICAL REVIEW E 97, 033309 (2018) - Section III.D
#

# Domain
Nx = 1500
Ny = 500

# Fluid properties (Re=500, We=8000, 1000:1)
rho_l = 1000.0
rho_g = 1.0
# mu_l = 20.0
# mu_g = 0.2
nu_l = 0.02
nu_g = 0.2
sigma = 0.0625

# Phase field parameters
tau_h = 0.8
D = 5

[Domain]
  dim = 2
  nx = '${Nx}'
  ny = '${Ny}'
  xmax = '${Nx}'
  ymax = '${Ny}'
  device_names = 'cuda'
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [phi]
    type = LBMTensorBuffer
    buffer_type = ms
    file = phi.h5
  []
  [ux]
    type = LBMTensorBuffer
    buffer_type = ms
    file = ux.h5
  []
  [uy]
    type = LBMTensorBuffer
    buffer_type = ms
    file = uy.h5
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []

  [grad_phi]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [laplacian_phi]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [mu]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [forces]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [speed]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [pressure]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [rho]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [h]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [h_post_collision]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [h_eq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [relaxation_tensor]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_post_collision]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_eq]
    type = LBMTensorBuffer
    buffer_type = df
  []
[]

[TensorComputes/Initialize]
  [u_stack]
    type = LBMStackTensors
    buffer = velocity
    inputs = 'ux uy'
  []

  [grad_phi_init]
    type = LBMIsotropicGradient
    buffer = grad_phi
    scalar_field = phi
  []
  [rho_init]
    type = ParsedCompute
    buffer = rho
    expression = 'phi * (rho_l - rho_g) + rho_g'
    constant_names = 'rho_l rho_g'
    constant_expressions = '${rho_l} ${rho_g}'
    inputs = phi
  []
  [pressure_init]
    type = LBMConstantTensor
    buffer = pressure
    constants = 0.3
  []

  # Equilibrium
  [h_eq_init]
    type = LBMPhaseEquilibrium
    buffer = h_eq
    phi = phi
    velocity = velocity
  []
  [h_init]
    type = ParsedCompute
    buffer = h
    expression = 'h_eq'
    inputs = h_eq
  []
  [h_post_collision_init]
    type = ParsedCompute
    buffer = h_post_collision
    expression = 'h_eq'
    inputs = h_eq
  []
  [f_eq_init]
    type = LBMPressureCorrectedEquilibrium
    buffer = f_eq
    rho = rho
    velocity = velocity
    pressure = pressure
  []
  [f_init]
    type = ParsedCompute
    buffer = f
    expression = 'f_eq'
    inputs = f_eq
  []
  [f_post_collision_init]
    type = ParsedCompute
    buffer = f_post_collision
    expression = 'f_eq'
    inputs = f_eq
  []
[]

[TensorComputes/Solve]
  [compute_phi]
    type = LBMComputeDensity
    buffer = phi
    f = h
  []
  [grad_phi]
    type = LBMIsotropicGradient
    buffer = grad_phi
    scalar_field = phi
  []
  [laplacian_phi]
    type = LBMIsotropicLaplacian
    buffer = laplacian_phi
    scalar_field = phi
  []
  [potential]
    type = LBMComputeChemicalPotential
    buffer = mu
    phi = phi
    laplacian_phi = laplacian_phi
    thickness = D
    sigma = sigma
  []
  [compute_forces]
    type = LBMComputeSurfaceForces
    buffer = forces
    chemical_potential = mu
    grad_phi = grad_phi
  []
  [density]
    type = ParsedCompute
    buffer = rho
    expression = 'phi * (rho_l - rho_g) + rho_g'
    constant_names = 'rho_l rho_g'
    constant_expressions = '${rho_l} ${rho_g}'
    inputs = phi
  []
  [velocity]
    type = LBMComputeVelocity
    buffer = velocity
    f = f
    rho = rho
    enable_forces = true
    forces = forces
  []
  [h_eq]
    type = LBMPhaseEquilibrium
    buffer = h_eq
    phi = phi
    velocity = velocity
  []
  [phase_collision]
    type = LBMBGKCollision
    buffer = h_post_collision
    f = h
    feq = h_eq
    tau0 = tau_h
  []
  [apply_forces_phase]
    type = LBMAllenCahnSource
    buffer = h_post_collision
    phi = phi
    velocity = velocity
    grad_phi = grad_phi
    tau = tau_h
    thickness = D
  []
  [relaxation_tensor]
    type = ParsedCompute
    buffer = relaxation_tensor
    expression = '(phi * (nu_l - nu_g) + nu_g) / cs2 + 0.5'
    constant_names = 'nu_l nu_g cs2'
    constant_expressions = '${nu_l} ${nu_g} 0.333333333333'
    inputs = 'phi'
  []
  [pressure]
    type = LBMPhaseFieldPressure
    buffer = pressure
    f = f
    velocity = velocity
    grad_phi = grad_phi
    rho = rho
    rho_l = '${rho_l}'
    rho_g = '${rho_g}'
  []
  [f_eq]
    type = LBMPressureCorrectedEquilibrium
    buffer = f_eq
    rho = rho
    velocity = velocity
    pressure = pressure
  []
  [collision]
    type = LBMBGKCollision
    buffer = f_post_collision
    f = f
    feq = f_eq
    tau0 = 1.0
    is_dynamic_relaxation = true
    tau_tensor = relaxation_tensor
  []
  [apply_forces_hydro]
    type = LBMForceDistribution
    buffer = f_post_collision
    grad_phi = grad_phi
    velocity = velocity
    forces = forces
    tau_tensor = relaxation_tensor
    tau = 1.0
    rho_l = '${rho_l}'
    rho_g = '${rho_g}'
    is_dynamic_relaxation = true
  []
  [residual]
    type = LBMComputeResidual
    buffer = phi
    speed = phi
  []
[]

[TensorComputes/Boundary]
  # Bounce back on both top and bottom
  [top_fluid]
    type = LBMBounceBack
    buffer = f
    f_old = f_post_collision
    boundary = top
  []
  [bottom_fluid]
    type = LBMBounceBack
    buffer = f
    f_old = f_post_collision
    boundary = bottom
  []
  [top_phase]
    type = LBMBounceBack
    buffer = h
    f_old = h_post_collision
    boundary = top
  []
  [bottom_phase]
    type = LBMBounceBack
    buffer = h
    f_old = h_post_collision
    boundary = bottom
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = 'h f'
  f_old = 'h_post_collision f_post_collision'
  root_compute = residual
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 200
  print_debug_output = true
  scalar_constant_names = 'tau_h D sigma'
  scalar_constant_values = '${tau_h} ${D} ${sigma}'
[]

[Executioner]
  type = Transient
  num_steps = 50
[]

[TensorOutputs]
  [xdmf]
    type = XDMFTensorOutput
    buffer = 'phi velocity rho'
    output_mode = 'Cell Cell Cell'
    enable_hdf5 = true
    transpose = false
  []
[]

(examples/lbm/Phase-field/spinodal_decomposition.i)

#
# Spinodal decomposition
# PHYSICAL REVIEW E 97, 033309 (2018) - Section III.C
#

# Domain
Nx = 200
Ny = 200

# Fluid properties
rho_l = 1000.0
rho_g = 1.0
nu_l = 0.1
nu_g = 1.0
sigma = 0.2

# Phase field parameters
tau_h = 0.67
D = 4

[Domain]
  dim = 2
  nx = '${Nx}'
  ny = '${Ny}'
  xmax = '${Nx}'
  ymax = '${Ny}'
  device_names='cpu'
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  # Macroscopic phase field variables
  [phi]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [grad_phi]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [laplacian_phi]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [mu]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [forces]
    type = LBMTensorBuffer
    buffer_type = mv
  []

  # Macroscopic hydrodynamic variables
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [pressure]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [rho]
    type = LBMTensorBuffer
    buffer_type = ms
  []

  # LBM phase field variabels
  [h]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [h_post_collision]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [h_eq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [relaxation_tensor]
    type = LBMTensorBuffer
    buffer_type = ms
  []

  # LBM hydrodynamic variables
  [fdummy]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_post_collision]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_eq]
    type = LBMTensorBuffer
    buffer_type = df
  []
[]

[TensorComputes/Initialize]
  [phi_init]
    type = RandomTensor
    buffer = phi
    min = 0.3233
    max = 0.3433
  []
  [grad_phi_init]
    type = LBMIsotropicGradient
    buffer = grad_phi
    scalar_field = phi
  []
  [rho_init]
    type = ParsedCompute
    buffer = rho
    extra_symbols = true
    expression = 'phi*(rho_l - rho_g) + rho_g'
    constant_names = 'rho_l rho_g'
    constant_expressions = '${rho_l} ${rho_g}'
    inputs = phi
  []
  [pressure_init]
    type = LBMConstantTensor
    buffer = pressure
    constants = 0.3
  []
  # Phase field equilibrium distribution initialization
  [h_eq_init]
    type = LBMPhaseEquilibrium
    buffer = h_eq
    phi = phi
    velocity = velocity
  []
  [h_post_collision_init]
    type = LBMPhaseEquilibrium
    buffer = h_post_collision
    phi = phi
    velocity = velocity
  []
  [h_init]
    type = LBMPhaseEquilibrium
    buffer = h
    phi = phi
    velocity = velocity
  []
  # Hydrodynamic equilibrium distribution initialization
  [f_eq_init]
    type = LBMPressureCorrectedEquilibrium
    buffer = f_eq
    rho = rho
    velocity = velocity
    pressure = pressure
  []
  [f_post_collision_init]
    type = ParsedCompute
    buffer = f_post_collision
    expression = 'f_eq'
    inputs = f_eq
  []
  [f_init]
    type = ParsedCompute
    buffer = f
    expression = 'f_eq'
    inputs = f_eq
  []
[]

[TensorComputes/Solve]
  # Phase Field
  [compute_phi]
    type = LBMComputeDensity
    buffer = phi
    f = h
  []
  [grad_phi]
    type = LBMIsotropicGradient
    buffer = grad_phi
    scalar_field = phi
  []
  [laplacian_phi]
    type = LBMIsotropicLaplacian
    buffer = laplacian_phi
    scalar_field = phi
  []
  [potential]
    type = LBMComputeChemicalPotential
    buffer = mu
    phi = phi
    laplacian_phi = laplacian_phi
    thickness = D
    sigma = sigma
  []
  [forces]
    type = LBMComputeSurfaceForces
    buffer = forces
    chemical_potential = mu
    grad_phi = grad_phi
  []
  # Hydrodynamics
  [density]
    type = ParsedCompute
    buffer = rho
    extra_symbols = true
    expression = 'phi*(rho_l - rho_g) + rho_g'
    constant_names = 'rho_l rho_g'
    constant_expressions = '${rho_l} ${rho_g}'
    inputs = phi
  []
  [velocity]
    type = LBMComputeVelocity
    buffer = velocity
    f = f
    rho = rho
    enable_forces = true
    forces = forces
  []
  # Phase-field
  [h_eq]
    type = LBMPhaseEquilibrium
    buffer = h_eq
    phi = phi
    velocity = velocity
  []
  [phase_collision]
    type = LBMBGKCollision
    buffer = h_post_collision
    f = h
    feq = h_eq
    tau0 = tau_h
  []
  [apply_forces_phase]
    type = LBMAllenCahnSource
    buffer = h_post_collision
    phi = phi
    velocity = velocity
    grad_phi = grad_phi
    tau = tau_h
    thickness = D
  []
  # Hydrodynamics
  [relaxation_tensor]
    type = ParsedCompute
    buffer = relaxation_tensor
    extra_symbols = true
    expression = '(phi*(nu_l - nu_g) + nu_g)/cs2+0.5'
    constant_names = 'nu_l nu_g cs2'
    constant_expressions = '${nu_l} ${nu_g} 0.3333'
    inputs = phi
  []
  [pressure]
    type = LBMPhaseFieldPressure
    buffer = pressure
    f = f
    velocity = velocity
    grad_phi = grad_phi
    rho = rho
    rho_l = '${rho_l}'
    rho_g = '${rho_g}'
  []
  [f_eq]
    type = LBMPressureCorrectedEquilibrium
    buffer = f_eq
    rho = rho
    velocity = velocity
    pressure = pressure
  []
  [collision]
    type = LBMBGKCollision
    buffer = f_post_collision
    f = f
    feq = f_eq
    tau0 = 1.0
    is_dynamic_relaxation = true
    tau_tensor = relaxation_tensor
   []
  [apply_forces_hydro]
    type = LBMForceDistribution
    buffer = f_post_collision
    grad_phi = grad_phi
    velocity = velocity
    forces = forces
    tau_tensor = relaxation_tensor
    tau = 1.0
    rho_l = '${rho_l}'
    rho_g = '${rho_g}'
    is_dynamic_relaxation = true
  []
  [residual]
    type = LBMComputeResidual
    buffer = phi
    speed = phi
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = 'h f'
  f_old = 'h_post_collision f_post_collision'
  root_compute = residual
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 500
  print_debug_output = true
  scalar_constant_names = 'tau_h D sigma'
  scalar_constant_values = '${tau_h} ${D} ${sigma}'
  log_interval = 100
[]

[Executioner]
  type = Transient
  num_steps = 100
[]

[TensorOutputs]
  [xdmf]
    type = XDMFTensorOutput
    buffer = 'phi rho velocity'
    output_mode = 'Cell Cell Cell'
    enable_hdf5 = true
    # transpose = false
  []
[]

(test/tests/lbm/channel2D.i)

[Domain]
  dim = 2
  nx = 10
  ny = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type=LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type=LBMTensorBuffer
    buffer_type = ms
  []
  [speed]
    type=LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0'
    []
    [initial_equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_distribution]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_distribution_pc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type=LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type=LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
      add_body_force = true
      body_force_x = 0.0001
    []
    [speed]
      type = LBMComputeVelocityMagnitude
      buffer = speed
      velocity = velocity
    []
    [residual]
      type = LBMComputeResidual
      buffer = speed
      speed = speed
    []
  []
  [Boundary]
    [top]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = top
    []
    [bottom]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = bottom
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 10
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [speed_min]
    type = TensorExtremeValuePostprocessor
    buffer = speed
    value_type = MIN
  []
  [speed_max]
    type = TensorExtremeValuePostprocessor
    buffer = speed
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [densty_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 11
[]

[Outputs]
  file_base = channel2D
  csv = true
[]

(examples/lbm/Karman-vortex/cylinder.i)

[Domain]
  dim = 2
  nx = 800
  ny = 200
  xmax = 4
  ymax = 1
  device_names='cpu'
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [rho]
    type=LBMTensorBuffer
    buffer_type = ms
  []
  [u]
    type=LBMTensorBuffer
    buffer_type = mv
  []
  [speed]
    type=LBMTensorBuffer
    buffer_type = ms
  []
  [f]
    type=LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type=LBMTensorBuffer
    buffer_type = df
  []
  [f_post_collision]
    type=LBMTensorBuffer
    buffer_type = df
  []
  [binary_media]
    type = LBMTensorBuffer
    file = binary_media.h5
    is_integer = true
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [rho]
      type=LBMConstantTensor
      buffer=rho
      constants = rho0
    []
    [u]
      type=LBMConstantTensor
      buffer=u
      constants = 'Ux Uy'
    []
    [speed]
      type=LBMComputeVelocityMagnitude
      buffer=speed
      velocity=u
    []
    [feq]
      type=LBMEquilibrium
      buffer=feq
      bulk=rho
      velocity=u
    []
    [f]
      type=LBMEquilibrium
      buffer=f
      bulk=rho
      velocity=u
    []
    [f_post_coll]
      type=LBMEquilibrium
      buffer=f_post_collision
      bulk=rho
      velocity=u
    []
  []
  [Solve]
    [Density]
      type = LBMComputeDensity
      buffer=rho
      f = f
    []
    [Velocity]
      type = LBMComputeVelocity
      buffer=u
      f = f
      rho = rho
    []
    [Equilibrium]
      type = LBMEquilibrium
      buffer=feq
      bulk=rho
      velocity=u
    []
    [Collision]
      type = LBMSmagorinskyCollision
      buffer = f_post_collision
      f = f
      feq = feq
      tau0 = tau
      projection=true
    []
    [Speed]
      type = LBMComputeVelocityMagnitude
      buffer = speed
      velocity = u
    []
    [Residual]
      type = LBMComputeResidual
      buffer = speed
      speed = speed
    []
  []
  [Boundary]
    [wall]
      type = LBMBounceBack
      buffer = f
      f_old = f_post_collision
      boundary = wall
    []
    [left]
      type = LBMFixedFirstOrderBC
      buffer=f
      f=f
      value='Ux'
      perturb=true
      boundary=left
    []
    [right]
      type = LBMMicroscopicZeroGradientBC
      buffer=f
      boundary=right
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = f_post_collision
  root_compute=residual
[]

[Postprocessors]
  [rho_avg]
    type = TensorAveragePostprocessor
    buffer = rho
    execute_on = 'TIMESTEP_END'
  []
  [speed_avg]
    type = TensorAveragePostprocessor
    buffer = speed
    execute_on = 'TIMESTEP_END'
  []
  [reynolds]
    type = ComputeReynoldsNumber
    buffer = speed
    tau = tau
    diameter = D
  []
[]

[Problem]
  type = LatticeBoltzmannProblem
  scalar_constant_names = 'rho0 Ux    Uy  tau    dx    D    Cs'
  scalar_constant_values = '1.0 0.01 0.0 0.506 0.001  40  0.15'
  substeps = 100
  print_debug_output=true
  binary_media = binary_media
[]

[Executioner]
  type = Transient
  num_steps = 10000
[]

[TensorOutputs]
  [xdmf2]
    type = XDMFTensorOutput
    buffer = 'rho u speed binary_media'
    output_mode = 'Cell Cell Cell Cell'
    enable_hdf5 = true
  []
[]

(test/tests/lbm/convective_outflow_3d_right.i)

[Domain]
  dim = 3
  nx = 8
  ny = 8
  nz = 8
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [d3q19]
    type = LBMD3Q19
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type = LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0 0.0'
    []
    [initial_feq]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_fpc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type = LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
    []
  []
  [Boundary]
    [left]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.005
      uy = 0.0
      uz = 0.0
      order = first
      boundary = left
    []
    [right]
      type = LBMConvectiveOutflow
      buffer = f
      f_old = f
      convection_velocity = Uc
      boundary = right
    []
    [top]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = top
    []
    [bottom]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = bottom
    []
    [front]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = front
    []
    [back]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = back
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  scalar_constant_names = 'Uc'
  scalar_constant_values = '0.005'
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = convective_outflow_3d_right
  csv = true
[]

(test/tests/lbm/mixed_bcs_d3q27.i)

[Domain]
  dim = 3
  nx = 10
  ny = 10
  nz = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [d3q27]
    type = LBMD3Q27
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_bounce_back]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type=LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type=LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0001 0.0005 0.0'
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_f_bb]
      type = LBMEquilibrium
      buffer = f_bounce_back
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
      add_body_force = true
      body_force_x = 0.0001
    []
  []
  [Boundary]
    [left]
      type = LBMFixedZerothOrderBC
      buffer = f
      f = f
      value = 1.1
      boundary = left
    []
    [right]
      type = LBMFixedFirstOrderBC
      buffer = f
      f = f
      value = 0.0001
      boundary = right
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = f_bounce_back
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [densty_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = mixed_bcs_d3q27
  csv = true
[]

(test/tests/lbm/rarefied_gas.i)

[Domain]
  dim = 2
  nx = 48
  ny = 48
  xmax = 48
  ymax = 48
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type=LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type=LBMTensorBuffer
    buffer_type = ms
  []
  [speed]
    type=LBMTensorBuffer
    buffer_type = ms
  []
  [domain]
    type=LBMTensorBuffer
    file = nanochannel.h5
    buffer_type = ms
    is_integer = true
  []
  [local_pore]
    type=LBMTensorBuffer
    file = local_pore_size.h5
    buffer_type = ms
    is_integer = false
  []
  [Kn]
    type=LBMTensorBuffer
    file = Kn.h5
    buffer_type = ms
    is_integer = false
  []
  [relaxation_matrix]
    type = LBMTensorBuffer
    buffer_type = df
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 0.2355545440759889
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0'
    []
    [initial_equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_distribution]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_distribution_pc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
    [relaxation_matrix_init]
      type = LBMComputeEffectiveRelaxation
      buffer = relaxation_matrix
      local_pore_size = local_pore
      local_Knudsen_number = Kn
      mfp = 7.904614716131531e-10
      dx = 1.0e-9
      A2 = 0.82
    []
  []
  [Solve]
    [equilibrium]
      type=LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type = LBMMRTCollision
      buffer = fpc
      f = f
      feq = feq
      is_dynamic_relaxation = true
      local_relaxation_matrix = relaxation_matrix
      projection = true
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
      add_body_force = true
      body_force_x = 1.0e-9
    []
    [speed]
      type = LBMComputeVelocityMagnitude
      buffer = speed
      velocity = velocity
    []
    [residual]
      type = LBMComputeResidual
      buffer = speed
      speed = speed
    []
  []
  [Boundary]
    [wall]
      type = LBMSpecularReflectionBoundary
      buffer = f
      f_old = fpc
      local_Knudsen_number = Kn
      boundary = wall
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [speed_min]
    type = TensorExtremeValuePostprocessor
    buffer = speed
    value_type = MIN
  []
  [speed_max]
    type = TensorExtremeValuePostprocessor
    buffer = speed
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [densty_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 10
  print_debug_output = true
  binary_media = domain
[]

[Executioner]
  type = Transient
  num_steps = 10
[]

[Outputs]
  file_base = rarefied_gas
  csv = true
[]

(test/tests/lbm/obstacle.i)

[Domain]
  dim = 2
  nx = 11
  ny = 11
  mesh_mode=DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [rho]
    type=LBMTensorBuffer
    buffer_type = ms
  []
  [u]
    type=LBMTensorBuffer
    buffer_type = mv
  []
  [speed]
    type=LBMTensorBuffer
    buffer_type = ms
  []
  [f]
    type=LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type=LBMTensorBuffer
    buffer_type = df
  []
  [f_post_collision]
    type=LBMTensorBuffer
    buffer_type = df
  []
  [binary_media]
    type = LBMTensorBuffer
    file = binary_media.h5
    is_integer = true
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [rho]
      type=LBMConstantTensor
      buffer=rho
      constants = rho0
    []
    [u]
      type=LBMConstantTensor
      buffer=u
      constants = 'Ux Uy'
    []
    [speed]
      type=LBMComputeVelocityMagnitude
      buffer=speed
      velocity=u
    []
    [feq]
      type=LBMEquilibrium
      buffer=feq
      bulk=rho
      velocity=u
    []
    [f]
      type=LBMEquilibrium
      buffer=f
      bulk=rho
      velocity=u
    []
    [f_post_coll]
      type=LBMEquilibrium
      buffer=f_post_collision
      bulk=rho
      velocity=u
    []
  []
  [Solve]
    [Density]
      type = LBMComputeDensity
      buffer=rho
      f = f
    []
    [Velocity]
      type = LBMComputeVelocity
      buffer=u
      f = f
      rho = rho
    []
    [Equilibrium]
      type = LBMEquilibrium
      buffer=feq
      bulk=rho
      velocity=u
    []
    [Collision]
      type = LBMSmagorinskyCollision
      buffer = f_post_collision
      f = f
      feq = feq
      tau0 = tau
      projection=true
    []
    [Speed]
      type = LBMComputeVelocityMagnitude
      buffer = speed
      velocity = u
    []
    [Residual]
      type = LBMComputeResidual
      buffer = speed
      speed = speed
    []
  []
  [Boundary]
    [wall]
      type = LBMBounceBack
      buffer = f
      f_old = f_post_collision
      boundary = wall
    []
    [left]
      type = LBMFixedFirstOrderBC
      buffer=f
      f=f
      value='Ux'
      perturb=true
      boundary=left
    []
    [right]
      type = LBMMicroscopicZeroGradientBC
      buffer=f
      boundary=right
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = f_post_collision
[]

[Problem]
  type = LatticeBoltzmannProblem
  scalar_constant_names = 'rho0 Ux    Uy tau   dx    D    Cs'
  scalar_constant_values = '1.0 0.01 0.0 0.501 0.001 50 0.15'
  substeps = 5
  print_debug_output=true
  binary_media = binary_media
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = u
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = u
    value_type = MAX
  []
  [speed_min]
    type = TensorExtremeValuePostprocessor
    buffer = speed
    value_type = MIN
  []
  [speed_max]
    type = TensorExtremeValuePostprocessor
    buffer = speed
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = rho
    value_type = MIN
  []
  [densty_max]
    type = TensorExtremeValuePostprocessor
    buffer = rho
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = obstacle
  csv = true
[]

(test/tests/lbm/nee_2d_all_walls.i)

[Domain]
  dim = 2
  nx = 10
  ny = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type = LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0'
    []
    [initial_feq]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_fpc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type = LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
    []
  []
  [Boundary]
    [top]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.005
      uy = 0.0
      order = second
      boundary = top
    []
    [bottom]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = 0.0
      order = second
      boundary = bottom
    []
    [left]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = 0.0
      order = second
      boundary = left
    []
    [right]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = 0.0
      order = second
      boundary = right
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = nee_2d_all_walls
  csv = true
[]

(test/tests/lbm/isotropic_stencil_mrt.i)

[Domain]
  dim = 3
  nx = 10
  ny = 10
  nz = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [d3q27]
    type = LBMD3Q27
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type=LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type=LBMTensorBuffer
    buffer_type = ms
  []
  [speed]
    type=LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0 0.0'
    []
    [initial_equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_distribution]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_distribution_pc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type=LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type=LBMMRTCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
      add_body_force = true
      body_force_x = 0.0001
    []
    [speed]
      type = LBMComputeVelocityMagnitude
      buffer = speed
      velocity = velocity
    []
    [residual]
      type = LBMComputeResidual
      buffer = speed
      speed = speed
    []
  []
  [Boundary]
    [top]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = top
    []
    [bottom]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = bottom
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 11
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [speed_min]
    type = TensorExtremeValuePostprocessor
    buffer = speed
    value_type = MIN
  []
  [speed_max]
    type = TensorExtremeValuePostprocessor
    buffer = speed
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [densty_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = isotropic_stencil_mrt
  csv = true
[]

(test/tests/lbm/vertical_velocity_bcs.i)

[Domain]
  dim = 2
  nx = 10
  ny = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_bounce_back]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type=LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type=LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0001 0.0005'
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_f_bb]
      type = LBMEquilibrium
      buffer = f_bounce_back
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
      add_body_force = true
      body_force_x = 0.0001
    []
  []
  [Boundary]
    [left]
      type = LBMBounceBack
      buffer = f
      f_old = f_bounce_back
      boundary = left
    []
    [right]
      type = LBMBounceBack
      buffer = f
      f_old = f_bounce_back
      boundary = right
    []
    [top]
      type = LBMFixedFirstOrderBC
      buffer = f
      f = f
      value = 0.0001
      boundary = top
    []
    [bottom]
      type = LBMFixedFirstOrderBC
      buffer = f
      f = f
      value = 0.00011
      boundary = bottom
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = f_bounce_back
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 2
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [densty_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = vertical_velocity_bcs
  csv = true
[]

(examples/lbm/Phase-field/static_bubble.i)

#
# Static bubble test case
# PHYSICAL REVIEW E 97, 033309 (2018) - Section III.A
# Intended outcome is to have (spurious) velocity magnitude around 1.0e-10
#

# Domain
Nx = 200
Ny = 200

# Bubble parameters
Cx = '${Nx}/2.0'
Cy = '${Ny}/2.0'
R = 50

# Fluid properties
rho_l = 1000.0
rho_g = 1.0
nu_l = 0.1
nu_g = 0.1
sigma = 0.001

# Phase field parameters
tau_h = 0.7
D = 5

[Domain]
  dim = 2
  nx = '${Nx}'
  ny = '${Ny}'
  xmax = '${Nx}'
  ymax = '${Ny}'
  device_names='cpu'
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  # Macroscopic phase field variables
  [phi]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [grad_phi]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [laplacian_phi]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [mu]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [forces]
    type = LBMTensorBuffer
    buffer_type = mv
  []

  # Macroscopic hydrodynamic variables
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [pressure]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [rho]
    type = LBMTensorBuffer
    buffer_type = ms
  []

  # LBM phase field variabels
  [h]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [h_post_collision]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [h_eq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [relaxation_tensor]
    type = LBMTensorBuffer
    buffer_type = ms
  []

  # LBM hydrodynamic variables
  [fdummy]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_post_collision]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_eq]
    type = LBMTensorBuffer
    buffer_type = df
  []
[]

[TensorComputes/Initialize]
  [phi_init]
    type = ParsedCompute
    buffer = phi
    expression = '0.5 + 0.5 * tanh(2*(R - sqrt((x - Cx)^2 + (y - Cy)^2)) / D)'
    constant_names = 'Cx Cy R D'
    constant_expressions = '${Cx} ${Cy} ${R} ${D}'
    extra_symbols = true
  []
  [grad_phi_init]
    type = LBMIsotropicGradient
    buffer = grad_phi
    scalar_field = phi
  []
  [rho_init]
    type = ParsedCompute
    buffer = rho
    extra_symbols = true
    expression = 'phi*(rho_l - rho_g) + rho_g'
    constant_names = 'rho_l rho_g'
    constant_expressions = '${rho_l} ${rho_g}'
    inputs = phi
  []
  [pressure_init]
    type = LBMConstantTensor
    buffer = pressure
    constants = 0.3
  []
  # Phase field equilibrium distribution initialization
  [h_eq_init]
    type = LBMPhaseEquilibrium
    buffer = h_eq
    phi = phi
    velocity = velocity
  []
  [h_post_collision_init]
    type = LBMPhaseEquilibrium
    buffer = h_post_collision
    phi = phi
    velocity = velocity
  []
  [h_init]
    type = LBMPhaseEquilibrium
    buffer = h
    phi = phi
    velocity = velocity
  []
  # Hydrodynamic equilibrium distribution initialization
  [f_eq_init]
    type = LBMPressureCorrectedEquilibrium
    buffer = f_eq
    rho = rho
    velocity = velocity
    pressure = pressure
  []
  [f_post_collision_init]
    type = ParsedCompute
    buffer = f_post_collision
    expression = 'f_eq'
    inputs = f_eq
  []
  [f_init]
    type = ParsedCompute
    buffer = f
    expression = 'f_eq'
    inputs = f_eq
  []
[]

[TensorComputes/Solve]
  # Phase Field
  [compute_phi]
    type = LBMComputeDensity
    buffer = phi
    f = h
  []
  [grad_phi]
    type = LBMIsotropicGradient
    buffer = grad_phi
    scalar_field = phi
  []
  [laplacian_phi]
    type = LBMIsotropicLaplacian
    buffer = laplacian_phi
    scalar_field = phi
  []
  [potential]
    type = LBMComputeChemicalPotential
    buffer = mu
    phi = phi
    laplacian_phi = laplacian_phi
    thickness = D
    sigma = sigma
  []
  [forces]
    type = LBMComputeSurfaceForces
    buffer = forces
    chemical_potential = mu
    grad_phi = grad_phi
  []
  # Hydrodynamics
  [density]
    type = ParsedCompute
    buffer = rho
    extra_symbols = true
    expression = 'phi*(rho_l - rho_g) + rho_g'
    constant_names = 'rho_l rho_g'
    constant_expressions = '${rho_l} ${rho_g}'
    inputs = phi
  []
  [velocity]
    type = LBMComputeVelocity
    buffer = velocity
    f = f
    rho = rho
    enable_forces = true
    forces = forces
  []
  # Phase-field
  [h_eq]
    type = LBMPhaseEquilibrium
    buffer = h_eq
    phi = phi
    velocity = velocity
  []
  [phase_collision]
    type = LBMBGKCollision
    buffer = h_post_collision
    f = h
    feq = h_eq
    tau0 = tau_h
  []
  [apply_forces_phase]
    type = LBMAllenCahnSource
    buffer = h_post_collision
    phi = phi
    velocity = velocity
    grad_phi = grad_phi
    tau = tau_h
    thickness = D
  []
  # Hydrodynamics
  [relaxation_tensor]
    type = ParsedCompute
    buffer = relaxation_tensor
    extra_symbols = true
    expression = '(phi*(nu_l - nu_g) + nu_g)/cs2+0.5'
    constant_names = 'nu_l nu_g cs2'
    constant_expressions = '${nu_l} ${nu_g} 0.3333'
    inputs = phi
  []
  [pressure]
    type = LBMPhaseFieldPressure
    buffer = pressure
    f = f
    velocity = velocity
    grad_phi = grad_phi
    rho = rho
    rho_l = '${rho_l}'
    rho_g = '${rho_g}'
  []
  [f_eq]
    type = LBMPressureCorrectedEquilibrium
    buffer = f_eq
    rho = rho
    velocity = velocity
    pressure = pressure
  []
  [collision]
    type = LBMBGKCollision
    buffer = f_post_collision
    f = f
    feq = f_eq
    tau0 = 1.0
    is_dynamic_relaxation = true
    tau_tensor = relaxation_tensor
   []
  [apply_forces_hydro]
    type = LBMForceDistribution
    buffer = f_post_collision
    grad_phi = grad_phi
    velocity = velocity
    forces = forces
    tau_tensor = relaxation_tensor
    tau = 1.0
    rho_l = '${rho_l}'
    rho_g = '${rho_g}'
    is_dynamic_relaxation = true
  []
  [residual]
    type = LBMComputeResidual
    buffer = phi
    speed = phi
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = 'h f'
  f_old = 'h_post_collision f_post_collision'
  root_compute = residual
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 50
  print_debug_output = true
  scalar_constant_names = 'tau_h D sigma'
  scalar_constant_values = '${tau_h} ${D} ${sigma}'
[]

[Executioner]
  type = Transient
  num_steps = 100
[]

[TensorOutputs]
  [xdmf]
    type = XDMFTensorOutput
    buffer = 'phi rho velocity forces'
    output_mode = 'Cell Cell Cell Cell'
    enable_hdf5 = true
    # transpose = false
  []
[]

(test/tests/lbm/phase.i)

# Domain
Nx = 20
Ny = 20

# Fluid properties
rho_l = 5.0
rho_g = 1.0
nu_l = 0.1
nu_g = 1.0
sigma = 0.2

# Phase field parameters
tau_h = 1.0
D = 4

[Domain]
  dim = 2
  nx = '${Nx}'
  ny = '${Ny}'
  xmax = '${Nx}'
  ymax = '${Ny}'
  device_names='cpu'
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  # Macroscopic phase field variables
  [phi]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [grad_phi]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [laplacian_phi]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [mu]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [forces]
    type = LBMTensorBuffer
    buffer_type = mv
  []

  # Macroscopic hydrodynamic variables
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [pressure]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [rho]
    type = LBMTensorBuffer
    buffer_type = ms
  []

  # LBM phase field variabels
  [h]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [h_post_collision]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [h_eq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [relaxation_tensor]
    type = LBMTensorBuffer
    buffer_type = ms
  []

  # LBM hydrodynamic variables
  [fdummy]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_post_collision]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_eq]
    type = LBMTensorBuffer
    buffer_type = df
  []
[]

[TensorComputes/Initialize]
  [phi_init]
    type = ParsedCompute
    buffer = phi
    extra_symbols = true
    expression = '0.3333 + 0.01*sin((12.9898*x + 78.233*y)*2*pi)'
  []
  [grad_phi_init]
    type = LBMIsotropicGradient
    buffer = grad_phi
    scalar_field = phi
  []
  [rho_init]
    type = ParsedCompute
    buffer = rho
    extra_symbols = true
    expression = 'phi*(rho_l - rho_g) + rho_g'
    constant_names = 'rho_l rho_g'
    constant_expressions = '${rho_l} ${rho_g}'
    inputs = phi
  []
  [pressure_init]
    type = LBMConstantTensor
    buffer = pressure
    constants = 0.3
  []
  # Phase field equilibrium distribution initialization
  [h_eq_init]
    type = LBMPhaseEquilibrium
    buffer = h_eq
    phi = phi
    velocity = velocity
  []
  [h_post_collision_init]
    type = LBMPhaseEquilibrium
    buffer = h_post_collision
    phi = phi
    velocity = velocity
  []
  [h_init]
    type = LBMPhaseEquilibrium
    buffer = h
    phi = phi
    velocity = velocity
  []
  # Hydrodynamic equilibrium distribution initialization
  [f_eq_init]
    type = LBMPressureCorrectedEquilibrium
    buffer = f_eq
    rho = rho
    velocity = velocity
    pressure = pressure
  []
  [f_post_collision_init]
    type = ParsedCompute
    buffer = f_post_collision
    expression = 'f_eq'
    inputs = f_eq
  []
  [f_init]
    type = ParsedCompute
    buffer = f
    expression = 'f_eq'
    inputs = f_eq
  []
[]

[TensorComputes/Solve]
  # Phase Field
  [compute_phi]
    type = LBMComputeDensity
    buffer = phi
    f = h
  []
  [grad_phi]
    type = LBMIsotropicGradient
    buffer = grad_phi
    scalar_field = phi
  []
  [laplacian_phi]
    type = LBMIsotropicLaplacian
    buffer = laplacian_phi
    scalar_field = phi
  []
  [potential]
    type = LBMComputeChemicalPotential
    buffer = mu
    phi = phi
    laplacian_phi = laplacian_phi
    thickness = D
    sigma = sigma
  []
  [forces]
    type = LBMComputeSurfaceForces
    buffer = forces
    chemical_potential = mu
    grad_phi = grad_phi
  []
  # Hydrodynamics
  [density]
    type = ParsedCompute
    buffer = rho
    extra_symbols = true
    expression = 'phi*(rho_l - rho_g) + rho_g'
    constant_names = 'rho_l rho_g'
    constant_expressions = '${rho_l} ${rho_g}'
    inputs = phi
  []
  [velocity]
    type = LBMComputeVelocity
    buffer = velocity
    f = f
    rho = rho
    enable_forces = true
    forces = forces
  []
  # Phase-field
  [h_eq]
    type = LBMPhaseEquilibrium
    buffer = h_eq
    phi = phi
    velocity = velocity
  []
  [phase_collision]
    type = LBMBGKCollision
    buffer = h_post_collision
    f = h
    feq = h_eq
    tau0 = tau_h
  []
  [apply_forces_phase]
    type = LBMAllenCahnSource
    buffer = h_post_collision
    phi = phi
    velocity = velocity
    grad_phi = grad_phi
    tau = tau_h
    thickness = D
  []
  # Hydrodynamics
  [relaxation_tensor]
    type = ParsedCompute
    buffer = relaxation_tensor
    extra_symbols = true
    expression = '(phi*(nu_l - nu_g) + nu_g)/cs2+0.5'
    constant_names = 'nu_l nu_g cs2'
    constant_expressions = '${nu_l} ${nu_g} 0.3333'
    inputs = phi
  []
  [pressure]
    type = LBMPhaseFieldPressure
    buffer = pressure
    f = f
    velocity = velocity
    grad_phi = grad_phi
    rho = rho
    rho_l = '${rho_l}'
    rho_g = '${rho_g}'
  []
  [f_eq]
    type = LBMPressureCorrectedEquilibrium
    buffer = f_eq
    rho = rho
    velocity = velocity
    pressure = pressure
  []
  [collision]
    type = LBMBGKCollision
    buffer = f_post_collision
    f = f
    feq = f_eq
    tau0 = 1.0
    is_dynamic_relaxation = true
    tau_tensor = relaxation_tensor
   []
  [apply_forces_hydro]
    type = LBMForceDistribution
    buffer = f_post_collision
    grad_phi = grad_phi
    velocity = velocity
    forces = forces
    tau_tensor = relaxation_tensor
    tau = 1.0
    rho_l = '${rho_l}'
    rho_g = '${rho_g}'
    is_dynamic_relaxation = true
  []
  [residual]
    type = LBMComputeResidual
    buffer = phi
    speed = phi
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = 'h f'
  f_old = 'h_post_collision f_post_collision'
[]

[Postprocessors]
  [phi_min]
    type = TensorExtremeValuePostprocessor
    buffer = phi
    value_type = MIN
  []
  [phi_max]
    type = TensorExtremeValuePostprocessor
    buffer = phi
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = rho
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = rho
    value_type = MAX
  []
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 5
  print_debug_output = true
  scalar_constant_names = 'tau_h D sigma'
  scalar_constant_values = '${tau_h} ${D} ${sigma}'
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = phase
  csv = true
[]

(test/tests/lbm/convective_outflow_3d_top.i)

[Domain]
  dim = 3
  nx = 8
  ny = 8
  nz = 8
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [d3q19]
    type = LBMD3Q19
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type = LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0 0.0'
    []
    [initial_feq]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_fpc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type = LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
    []
  []
  [Boundary]
    [bottom]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = 0.005
      uz = 0.0
      order = first
      boundary = bottom
    []
    [top]
      type = LBMConvectiveOutflow
      buffer = f
      f_old = f
      convection_velocity = Uc
      boundary = top
    []
    [left]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = left
    []
    [right]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = right
    []
    [front]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = front
    []
    [back]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = back
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  scalar_constant_names = 'Uc'
  scalar_constant_values = '0.005'
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = convective_outflow_3d_top
  csv = true
[]

(test/tests/lbm/smagorinsky_mrt.i)

[Domain]
  dim = 2
  nx = 10
  ny = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type=LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type=LBMTensorBuffer
    buffer_type = ms
  []
  [speed]
    type=LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0'
    []
    [initial_equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_distribution]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_distribution_pc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [root]
      [equilibrium]
        type=LBMEquilibrium
        buffer = feq
        bulk = density
        velocity = velocity
      []
      [collision]
        type=LBMSmagorinskyMRTCollision
        buffer = fpc
        f = f
        feq = feq
        tau0 = 1.0
      []
      [density]
        type = LBMComputeDensity
        buffer = density
        f = f
      []
      [velocity]
        type = LBMComputeVelocity
        buffer = velocity
        f = f
        rho = density
        add_body_force = true
        body_force_x = 0.0001
      []
      [speed]
        type = LBMComputeVelocityMagnitude
        buffer = speed
        velocity = velocity
      []
      [residual]
        type = LBMComputeResidual
        buffer = speed
        speed = speed
      []
    []
  []
  [Boundary]
    [top]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = top
    []
    [bottom]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = bottom
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
  root_compute = root
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 11
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [speed_min]
    type = TensorExtremeValuePostprocessor
    buffer = speed
    value_type = MIN
  []
  [speed_max]
    type = TensorExtremeValuePostprocessor
    buffer = speed
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [densty_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = smagorinsky_mrt
  csv = true
[]

(examples/lbm/Formula1-aerodynamics/f1.i)

[Domain]
  dim = 3
  nx = 781
  ny = 274
  nz = 146
  xmax = 781
  ymax = 274
  zmax = 146
  device_names='cuda'
  floating_precision = 'single'
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [descriptor]
    type = LBMD3Q27
  []
[]

[TensorBuffers]
  [binary_media]
    type = LBMTensorBuffer
    buffer_type = ms
    file = 'binary_media.h5'
    is_integer = true
  []

  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [u]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [speed]
    type=LBMTensorBuffer
    buffer_type = ms
  []
  [rho]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [force]
    type = LBMTensorBuffer
    buffer_type = mv
  []
[]

[TensorComputes]
  [Initialize]
    [density_initial]
      type = LBMConstantTensor
      buffer = rho
      constants = 1.0
    []

    [velocity_initial]
      type = LBMConstantTensor
      buffer = u
      constants = '0 0 0'
    []

    [equilibrium_init]
      type = LBMEquilibrium
      buffer = feq
      bulk = rho
      velocity = u
    []

    [equilibrium_f]
      type = LBMEquilibrium
      buffer = f
      bulk = rho
      velocity = u
    []

    [equilibrium_pc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = rho
      velocity = u
    []
  []

  [Solve]
    [density]
      type = LBMComputeDensity
      buffer = rho
      f = f
    []

    [velocity]
      type = LBMComputeVelocity
      buffer = u
      f = f
      rho = rho
    []

    [equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = rho
      velocity = u
    []

    [collision]
      type = LBMSmagorinskyCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 0.5001
      Cs = 0.15
      projection=true
    []

    [speed]
      type=LBMComputeVelocityMagnitude
      buffer=speed
      velocity=u
    []

    [residual]
      type = LBMComputeResidual
      buffer = speed
      speed = speed
    []
  []

  [Boundary]
    [wall]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = wall
    []

    [left]
      type = LBMFixedFirstOrderBC
      buffer=f
      f=f
      value=0.01
      boundary=left
      # perturb=true
    []

    [right]
      type = LBMMicroscopicZeroGradientBC
      buffer=f
      boundary = right
    []

  []
[]

[TensorSolver]
  type = LBMStream
  root_compute = residual
  buffer = f
  f_old = fpc
[]

[Postprocessors]
  [reynolds]
    type = ComputeReynoldsNumber
    buffer = speed
    tau = 0.5001
    diameter = 781
  []
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 100
  print_debug_output = true
  binary_media = binary_media
[]

[Executioner]
  type = Transient
  num_steps = 2000
[]

[TensorOutputs]
  [xdmf2]
    type = XDMFTensorOutput
    buffer = 'rho u'
    output_mode = 'Cell Cell'
    enable_hdf5 = true
  []
[]

(test/tests/lbm/mixed_bcs_d3q27_reverse.i)

[Domain]
  dim = 3
  nx = 10
  ny = 10
  nz = 10
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [d3q27]
    type = LBMD3Q27
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_bounce_back]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type=LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type=LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0001 0.0005 0.0'
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_f_bb]
      type = LBMEquilibrium
      buffer = f_bounce_back
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
      add_body_force = true
      body_force_x = 0.0001
    []
  []
  [Boundary]
    [right]
      type = LBMFixedZerothOrderBC
      buffer = f
      f = f
      value = 1.1
      boundary = right
    []
    [left]
      type = LBMFixedFirstOrderBC
      buffer = f
      f = f
      value = 0.0001
      boundary = left
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = f_bounce_back
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [densty_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = mixed_bcs_d3q27_reverse
  csv = true
[]

(test/tests/lbm/phase_3D.i)

# Domain
Nx = 20
Ny = 20
Nz = 20

# Fluid properties
rho_l = 5.0
rho_g = 1.0
nu_l = 0.1
nu_g = 1.0
sigma = 0.2

# Phase field parameters
tau_h = 1.0
D = 4

[Domain]
  dim = 3
  nx = '${Nx}'
  ny = '${Ny}'
  nz = '${Nz}'
  xmax = '${Nx}'
  ymax = '${Ny}'
  zmax = '${Nz}'
  device_names='cpu'
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [d3q27]
    type = LBMD3Q27
  []
[]

[TensorBuffers]
  # Macroscopic phase field variables
  [phi]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [grad_phi]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [laplacian_phi]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [mu]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [forces]
    type = LBMTensorBuffer
    buffer_type = mv
  []

  # Macroscopic hydrodynamic variables
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [pressure]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [rho]
    type = LBMTensorBuffer
    buffer_type = ms
  []

  # LBM phase field variabels
  [h]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [h_post_collision]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [h_eq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [relaxation_tensor]
    type = LBMTensorBuffer
    buffer_type = ms
  []

  # LBM hydrodynamic variables
  [fdummy]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_post_collision]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_eq]
    type = LBMTensorBuffer
    buffer_type = df
  []
[]

[TensorComputes/Initialize]
  [phi_init]
    type = ParsedCompute
    buffer = phi
    extra_symbols = true
    expression = '0.3333 + 0.01*sin((12.9898*x + 78.233*y + 43.12*z)*2*pi)'
  []
  [grad_phi_init]
    type = LBMIsotropicGradient
    buffer = grad_phi
    scalar_field = phi
  []
  [rho_init]
    type = ParsedCompute
    buffer = rho
    extra_symbols = true
    expression = 'phi*(rho_l - rho_g) + rho_g'
    constant_names = 'rho_l rho_g'
    constant_expressions = '${rho_l} ${rho_g}'
    inputs = phi
  []
  [pressure_init]
    type = LBMConstantTensor
    buffer = pressure
    constants = 0.3
  []
  # Phase field equilibrium distribution initialization
  [h_eq_init]
    type = LBMPhaseEquilibrium
    buffer = h_eq
    phi = phi
    velocity = velocity
  []
  [h_post_collision_init]
    type = LBMPhaseEquilibrium
    buffer = h_post_collision
    phi = phi
    velocity = velocity
  []
  [h_init]
    type = LBMPhaseEquilibrium
    buffer = h
    phi = phi
    velocity = velocity
  []
  # Hydrodynamic equilibrium distribution initialization
  [f_eq_init]
    type = LBMPressureCorrectedEquilibrium
    buffer = f_eq
    rho = rho
    velocity = velocity
    pressure = pressure
  []
  [f_post_collision_init]
    type = ParsedCompute
    buffer = f_post_collision
    expression = 'f_eq'
    inputs = f_eq
  []
  [f_init]
    type = ParsedCompute
    buffer = f
    expression = 'f_eq'
    inputs = f_eq
  []
[]

[TensorComputes/Solve]
  # Phase Field
  [compute_phi]
    type = LBMComputeDensity
    buffer = phi
    f = h
  []
  [grad_phi]
    type = LBMIsotropicGradient
    buffer = grad_phi
    scalar_field = phi
  []
  [laplacian_phi]
    type = LBMIsotropicLaplacian
    buffer = laplacian_phi
    scalar_field = phi
  []
  [potential]
    type = LBMComputeChemicalPotential
    buffer = mu
    phi = phi
    laplacian_phi = laplacian_phi
    thickness = D
    sigma = sigma
  []
  [forces]
    type = LBMComputeSurfaceForces
    buffer = forces
    chemical_potential = mu
    grad_phi = grad_phi
  []
  # Hydrodynamics
  [density]
    type = ParsedCompute
    buffer = rho
    extra_symbols = true
    expression = 'phi*(rho_l - rho_g) + rho_g'
    constant_names = 'rho_l rho_g'
    constant_expressions = '${rho_l} ${rho_g}'
    inputs = phi
  []
  [velocity]
    type = LBMComputeVelocity
    buffer = velocity
    f = f
    rho = rho
    enable_forces = true
    forces = forces
  []
  # Phase-field
  [h_eq]
    type = LBMPhaseEquilibrium
    buffer = h_eq
    phi = phi
    velocity = velocity
  []
  [phase_collision]
    type = LBMBGKCollision
    buffer = h_post_collision
    f = h
    feq = h_eq
    tau0 = tau_h
  []
  [apply_forces_phase]
    type = LBMAllenCahnSource
    buffer = h_post_collision
    phi = phi
    velocity = velocity
    grad_phi = grad_phi
    tau = tau_h
    thickness = D
  []
  # Hydrodynamics
  [relaxation_tensor]
    type = ParsedCompute
    buffer = relaxation_tensor
    extra_symbols = true
    expression = '(phi*(nu_l - nu_g) + nu_g)/cs2+0.5'
    constant_names = 'nu_l nu_g cs2'
    constant_expressions = '${nu_l} ${nu_g} 0.3333'
    inputs = phi
  []
  [pressure]
    type = LBMPhaseFieldPressure
    buffer = pressure
    f = f
    velocity = velocity
    grad_phi = grad_phi
    rho = rho
    rho_l = '${rho_l}'
    rho_g = '${rho_g}'
  []
  [f_eq]
    type = LBMPressureCorrectedEquilibrium
    buffer = f_eq
    rho = rho
    velocity = velocity
    pressure = pressure
  []
  [collision]
    type = LBMBGKCollision
    buffer = f_post_collision
    f = f
    feq = f_eq
    tau0 = 1.0
    is_dynamic_relaxation = true
    tau_tensor = relaxation_tensor
   []
  [apply_forces_hydro]
    type = LBMForceDistribution
    buffer = f_post_collision
    grad_phi = grad_phi
    velocity = velocity
    forces = forces
    tau_tensor = relaxation_tensor
    tau = 1.0
    rho_l = '${rho_l}'
    rho_g = '${rho_g}'
    is_dynamic_relaxation = true
  []
  [residual]
    type = LBMComputeResidual
    buffer = phi
    speed = phi
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = 'h f'
  f_old = 'h_post_collision f_post_collision'
[]

[Postprocessors]
  [phi_min]
    type = TensorExtremeValuePostprocessor
    buffer = phi
    value_type = MIN
  []
  [phi_max]
    type = TensorExtremeValuePostprocessor
    buffer = phi
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = rho
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = rho
    value_type = MAX
  []
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 5
  print_debug_output = true
  scalar_constant_names = 'tau_h D sigma'
  scalar_constant_values = '${tau_h} ${D} ${sigma}'
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = phase_3D
  csv = true
[]

(test/tests/lbm/convective_outflow_3d_back.i)

[Domain]
  dim = 3
  nx = 8
  ny = 8
  nz = 8
  mesh_mode = DUMMY
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y Z'
[]

[Stencil]
  [d3q19]
    type = LBMD3Q19
  []
[]

[TensorBuffers]
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [density]
    type = LBMTensorBuffer
    buffer_type = ms
  []
[]

[TensorComputes]
  [Initialize]
    [initial_density]
      type = LBMConstantTensor
      buffer = density
      constants = 1.0
    []
    [initial_velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0 0.0'
    []
    [initial_feq]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [initial_f]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []
    [initial_fpc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []
  []
  [Solve]
    [equilibrium]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []
    [collision]
      type = LBMBGKCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = 1.0
    []
    [density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []
    [velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
    []
  []
  [Boundary]
    [front]
      type = LBMNonEquilibriumExtrapolation
      buffer = f
      prescribe_type = velocity
      ux = 0.0
      uy = 0.0
      uz = 0.005
      order = first
      boundary = front
    []
    [back]
      type = LBMConvectiveOutflow
      buffer = f
      f_old = f
      convection_velocity = Uc
      boundary = back
    []
    [left]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = left
    []
    [right]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = right
    []
    [top]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = top
    []
    [bottom]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = bottom
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = f
  f_old = fpc
[]

[Problem]
  type = LatticeBoltzmannProblem
  scalar_constant_names = 'Uc'
  scalar_constant_values = '0.005'
  substeps = 1
[]

[Postprocessors]
  [velocity_min]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MIN
  []
  [velocity_max]
    type = TensorExtremeValuePostprocessor
    buffer = velocity
    value_type = MAX
  []
  [density_min]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MIN
  []
  [density_max]
    type = TensorExtremeValuePostprocessor
    buffer = density
    value_type = MAX
  []
[]

[Executioner]
  type = Transient
  num_steps = 5
[]

[Outputs]
  file_base = convective_outflow_3d_back
  csv = true
[]

(examples/lbm/Phase-field/layered_poiseuille.i)

#
# Layered Poiseuille Flow
# PHYSICAL REVIEW E 97, 033309 (2018) - Section III.B
#

# Domain
Nx = 10
Ny = 400

# Fluid properties
rho_l = 1000.0
rho_g = 1.0
# nu_l = 0.1
# nu_g = 1.0
sigma = 0.001
mu_l = 100.0 # rho_l * nu_l
mu_g = 1.0   # rho_g * nu_g

# Phase field parameters
# M = 0.1
# cs2 = 0.333333333333
tau_h = 0.8 # 0.5 + '${M}' / '${cs2}'
D = 5

# Driving force: Gx = uc * (mu_l + mu_g) / h^2
# uc = 1e-4
# h = # '${Ny} / 2'
Gx = 2.53e-07 #  '${uc} * (${mu_l}  + ${mu_g}) / (${h}^2)'

[Domain]
  dim = 2
  nx = '${Nx}'
  ny = '${Ny}'
  xmax = '${Nx}'
  ymax = '${Ny}'
  device_names = 'cpu'
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d2q9]
    type = LBMD2Q9
  []
[]

[TensorBuffers]
  # Macroscopic phase field variables
  [phi]
    type = LBMTensorBuffer
    buffer_type = ms
    file = phi_init.h5
  []
  [grad_phi]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [laplacian_phi]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [mu]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [forces]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [body_force]
    type = LBMTensorBuffer
    buffer_type = mv
  []

  # Macroscopic hydrodynamic variables
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [speed]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [pressure]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [rho]
    type = LBMTensorBuffer
    buffer_type = ms
  []

  # LBM phase field variables
  [h]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [h_post_collision]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [h_eq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [relaxation_tensor]
    type = LBMTensorBuffer
    buffer_type = ms
  []

  # LBM hydrodynamic variables
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_post_collision]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [f_eq]
    type = LBMTensorBuffer
    buffer_type = df
  []
[]

[TensorComputes/Initialize]
  [grad_phi_init]
    type = LBMIsotropicGradient
    buffer = grad_phi
    scalar_field = phi
  []
  [rho_init]
    type = ParsedCompute
    buffer = rho
    expression = 'phi * (rho_l - rho_g) + rho_g'
    constant_names = 'rho_l rho_g'
    constant_expressions = '${rho_l} ${rho_g}'
    inputs = phi
  []
  [pressure_init]
    type = LBMConstantTensor
    buffer = pressure
    constants = 0.3
  []
  [body_force_init]
    type = LBMConstantTensor
    buffer = body_force
    constants = '${Gx} 0.00'
  []
  # Phase field equilibrium distribution initialization
  [h_eq_init]
    type = LBMPhaseEquilibrium
    buffer = h_eq
    phi = phi
    velocity = velocity
  []
  [h_post_collision_init]
    type = LBMPhaseEquilibrium
    buffer = h_post_collision
    phi = phi
    velocity = velocity
  []
  [h_init]
    type = LBMPhaseEquilibrium
    buffer = h
    phi = phi
    velocity = velocity
  []
  # Hydrodynamic equilibrium distribution initialization
  [f_eq_init]
    type = LBMPressureCorrectedEquilibrium
    buffer = f_eq
    rho = rho
    velocity = velocity
    pressure = pressure
  []
  [f_post_collision_init]
    type = ParsedCompute
    buffer = f_post_collision
    expression = 'f_eq'
    inputs = f_eq
  []
  [f_init]
    type = ParsedCompute
    buffer = f
    expression = 'f_eq'
    inputs = f_eq
  []
[]

[TensorComputes/Solve]
  # Phase Field
  [compute_phi]
    type = LBMComputeDensity
    buffer = phi
    f = h
  []
  [grad_phi]
    type = LBMIsotropicGradient
    buffer = grad_phi
    scalar_field = phi
  []
  [laplacian_phi]
    type = LBMIsotropicLaplacian
    buffer = laplacian_phi
    scalar_field = phi
  []
  [potential]
    type = LBMComputeChemicalPotential
    buffer = mu
    phi = phi
    laplacian_phi = laplacian_phi
    thickness = D
    sigma = sigma
  []
  [compute_forces]
    type = LBMComputeSurfaceForces
    buffer = forces
    chemical_potential = mu
    grad_phi = grad_phi
  []
  [add_body_force]
    type = ParsedCompute
    buffer = forces
    expression = 'forces + body_force'
    inputs = 'forces body_force'
  []
  # Hydrodynamics
  [density]
    type = ParsedCompute
    buffer = rho
    expression = 'phi * (rho_l - rho_g) + rho_g'
    constant_names = 'rho_l rho_g'
    constant_expressions = '${rho_l} ${rho_g}'
    inputs = phi
  []
  [velocity]
    type = LBMComputeVelocity
    buffer = velocity
    f = f
    rho = rho
    enable_forces = true
    forces = forces
    # add_body_force = true
    # body_force_x = '${Gx}'
  []
  # Phase-field
  [h_eq]
    type = LBMPhaseEquilibrium
    buffer = h_eq
    phi = phi
    velocity = velocity
  []
  [phase_collision]
    type = LBMBGKCollision
    buffer = h_post_collision
    f = h
    feq = h_eq
    tau0 = tau_h
  []
  [apply_forces_phase]
    type = LBMAllenCahnSource
    buffer = h_post_collision
    phi = phi
    velocity = velocity
    grad_phi = grad_phi
    tau = tau_h
    thickness = D
  []
  # Hydrodynamics
  [relaxation_tensor]
    type = ParsedCompute
    buffer = relaxation_tensor
    # Implements Eq 26: Sharp step-function for dynamic viscosity
    expression = '(if(phi >= 0.5, mu_l, mu_g) / rho) / cs2 + 0.5'
    constant_names = 'mu_l mu_g cs2'
    constant_expressions = '${mu_l} ${mu_g} 0.333333333333'
    inputs = 'phi rho'
  []
  [pressure]
    type = LBMPhaseFieldPressure
    buffer = pressure
    f = f
    velocity = velocity
    grad_phi = grad_phi
    rho = rho
    rho_l = '${rho_l}'
    rho_g = '${rho_g}'
  []
  [f_eq]
    type = LBMPressureCorrectedEquilibrium
    buffer = f_eq
    rho = rho
    velocity = velocity
    pressure = pressure
  []
  [collision]
    type = LBMBGKCollision
    buffer = f_post_collision
    f = f
    feq = f_eq
    tau0 = 1.0
    is_dynamic_relaxation = true
    tau_tensor = relaxation_tensor
  []
  [apply_forces_hydro]
    type = LBMForceDistribution
    buffer = f_post_collision
    grad_phi = grad_phi
    velocity = velocity
    forces = forces
    tau_tensor = relaxation_tensor
    tau = 1.0
    rho_l = '${rho_l}'
    rho_g = '${rho_g}'
    is_dynamic_relaxation = true
  []
  [speed]
    type = LBMComputeVelocityMagnitude
    buffer = speed
    velocity = velocity
  []
  [residual]
    type = LBMComputeResidual
    buffer = speed
    speed = speed
  []
[]

[TensorComputes/Boundary]
  [top_fluid]
    type = LBMBounceBack
    buffer = f
    f_old = f_post_collision
    boundary = top
  []
  [bottom_fluid]
    type = LBMBounceBack
    buffer = f
    f_old = f_post_collision
    boundary = bottom
  []
  [top_phase]
    type = LBMBounceBack
    buffer = h
    f_old = h_post_collision
    boundary = top
  []
  [bottom_phase]
    type = LBMBounceBack
    buffer = h
    f_old = h_post_collision
    boundary = bottom
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = 'h f'
  f_old = 'h_post_collision f_post_collision'
  root_compute = residual
[]

[Problem]
  type = LatticeBoltzmannProblem
  # Keep this low for this setup: high substeps with top/bottom bounce-back can blow up to NaN.
  substeps = 100000
  print_debug_output = true
  scalar_constant_names = 'tau_h D sigma'
  scalar_constant_values = '${tau_h} ${D} ${sigma}'
[]

[Executioner]
  type = Transient
  num_steps = 2
[]

[TensorOutputs]
  [xdmf]
    type = XDMFTensorOutput
    buffer = 'phi velocity rho'
    output_mode = 'Cell Cell Cell'
    enable_hdf5 = true
    transpose = false
  []
[]

(examples/lbm/Pebbles/pebbles.i)

[Domain]
  dim = 3
  nx = 240
  ny = 50
  nz = 50
  xmax = 240
  ymax = 50
  zmax = 50
  device_names='cpu'
  # floating_precision = 'single'
  parallel_mode = REAL_SPACE
  periodic_directions = 'X Y'
[]

[Stencil]
  [d3q19]
    type = LBMD3Q19
  []
[]

[TensorBuffers]
  # Density distribution functions
  [f]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [feq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [fpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  # Temperature distribution functions
  [g]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [geq]
    type = LBMTensorBuffer
    buffer_type = df
  []
  [gpc]
    type = LBMTensorBuffer
    buffer_type = df
  []
  # Fluid macroscopic variables: density and velocity
  [density]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  [velocity]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  [speed]
    type=LBMTensorBuffer
    buffer_type = ms
  []
  # Temperature macroscpic variables: temperature and 'velocity'
  [T]
    type = LBMTensorBuffer
    buffer_type = ms
  []
  # Forces
  [F]
    type = LBMTensorBuffer
    buffer_type = mv
  []
  # Simulation domain
  [binary_media]
    type = LBMTensorBuffer
    buffer_type = ms
    file = 'binary_media.h5'
    is_integer = true
  []
[]

[TensorComputes]

  #### Initialzie ####
  [Initialize]
    [density]
      type = LBMConstantTensor
      buffer = density
      constants = 'rho0'
    []

    [velocity]
      type = LBMConstantTensor
      buffer = velocity
      constants = '0.0 0.0'
    []

    [temperature]
      type = LBMConstantTensor
      buffer = T
      constants = T_C
    []

    [equilibrium_fluid]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []

    [equilibrium_fluid_total]
      type = LBMEquilibrium
      buffer = f
      bulk = density
      velocity = velocity
    []

    [equilibrium_fluid_pc]
      type = LBMEquilibrium
      buffer = fpc
      bulk = density
      velocity = velocity
    []

    [equilibrium_temperature]
      type = LBMEquilibrium
      buffer = geq
      bulk = T
      velocity = velocity
    []

    [equilibrium_temperature_total]
      type = LBMEquilibrium
      buffer = g
      bulk = T
      velocity = velocity
    []

    [equilibrium_temperature_pc]
      type = LBMEquilibrium
      buffer = gpc
      bulk = T
      velocity = velocity
    []

    # [t_wall_hot_init]
    #   type = LBMDirichletWallBC
    #   buffer = g
    #   f_old = gpc
    #   value = T_H
    #   velocity = velocity
    #   boundary = wall
    # []
  []

  #### Compute ####
  [Solve]

    # For temperature
    [Temperature]
      type = LBMComputeDensity
      buffer = T
      f = g
    []

    # For fluid
    [Fluid_density]
      type = LBMComputeDensity
      buffer = density
      f = f
    []

    [Fluid_velocity]
      type = LBMComputeVelocity
      buffer = velocity
      f = f
      rho = density
      forces = F
      enable_forces = true
    []

    # For temperature
    [Equilibrium_temperature]
      type = LBMEquilibrium
      buffer = geq
      bulk = T
      velocity = velocity
    []

    [Collision_temperature]
      type = LBMMRTCollision
      buffer = gpc
      f = g
      feq = geq
      tau0 = tau_T
    []

    # For fluid
    [Compute_forces]
      type = LBMComputeForces
      buffer = F
      rho0 = 'rho0'
      temperature = T
      T0 = T_H
      enable_buoyancy = true
      gravity = g
      gravity_direction=0
    []

    [Equilibrium_fluid]
      type = LBMEquilibrium
      buffer = feq
      bulk = density
      velocity = velocity
    []

    [Collision_fluid]
      type = LBMMRTCollision
      buffer = fpc
      f = f
      feq = feq
      tau0 = tau_f
    []

    [Apply_forces]
      type = LBMApplyForces
      buffer = fpc
      velocity = velocity
      rho = density
      forces = F
      tau0 = tau_f
    []

    [speed]
      type=LBMComputeVelocityMagnitude
      buffer=speed
      velocity=velocity
    []

    [residual]
      type = LBMComputeResidual
      buffer = speed
      speed = speed
    []
  []

  #### Boundary ####
  [Boundary]
    ##### for fluid
    [inlet]
      type = LBMFixedFirstOrderBC
      buffer = f
      f = f
      value = u0
      boundary = left
    []
    [outlet]
      type = LBMMicroscopicZeroGradientBC
      buffer = f
      boundary = right
    []
    [top]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = top
      # exclude_corners_x = true
    []
    [bottom]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = bottom
      # exclude_corners_x = true
    []
    [front]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = front
      # exclude_corners_x = true
    []
    [back]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = back
      # exclude_corners_x = true
    []
    [wall]
      type = LBMBounceBack
      buffer = f
      f_old = fpc
      boundary = wall
    []

    ##### for temperature
    [t_inlet]
      type = LBMFixedZerothOrderBC
      buffer = g
      f = g
      value = T_C
      boundary = left
    []
    [t_outlet]
      type = LBMMicroscopicZeroGradientBC
      buffer = g
      boundary = right
    []
    [t_top]
      type = LBMBounceBack
      buffer = g
      f_old = gpc
      boundary = top
    []
    [t_bottom]
      type = LBMBounceBack
      buffer = g
      f_old = gpc
      boundary = bottom
    []
    [t_front]
      type = LBMBounceBack
      buffer = g
      f_old = gpc
      boundary = front
    []
    [t_back]
      type = LBMBounceBack
      buffer = g
      f_old = gpc
      boundary = back
    []
    [t_wall_hot]
      type = LBMDirichletBC
      buffer = g
      f_old = gpc
      feq = geq
      velocity = velocity
      rho = T
      value = 1.2
      boundary = wall
    []
  []
[]

[TensorSolver]
  type = LBMStream
  buffer = 'f g'
  f_old = 'fpc gpc'
  root_compute = residual
[]

[Problem]
  type = LatticeBoltzmannProblem
  substeps = 20
  scalar_constant_names = 'rho0 T_C T_H tau_f tau_T g u0'
  scalar_constant_values = '1.0 1.0 1.2 0.9 0.9 0.001 0.001'

  is_binary_media = true
  binary_media = binary_media
[]

[Postprocessors]
  [reynolds]
    type = ComputeReynoldsNumber
    buffer = speed
    tau = tau_f
    diameter = 15
  []
[]

[Executioner]
  type = Transient
  num_steps = 2000
[]

# [TensorOutputs]
#   [xdmf2]
#     type = XDMFTensorOutput
#     buffer = 'binary_media'
#     output_mode = 'Cell'
#     enable_hdf5 = true
#   []
# []

[TensorOutputs]
  [xdmf2]
    type = XDMFTensorOutput
    buffer = 'T velocity density binary_media'
    output_mode = 'Cell Cell Cell Cell'
    enable_hdf5 = true
  []
[]

Overview
Example Input File Syntax
Input Parameters
Input Files