- bufferThe buffers to output
C++ Type:std::vector<std::string>
Controllable:No
Description:The buffers to output
XDMFTensorOutput
Output a tensor in XDMF format.
This TensorOutput object writes tensors to an XDMF data file set. The parameter "enable_hdf5" enables writing of the raw tensor data to an HDF5 file. The resulting data items will be bamed /var.n, where var is the name of the exported tensor and n is the simulation step. Output as cell or node centered data is supported and can be selected using the "output_mode" parameter. Cell centered output results in a value per simulation grid cell (e.g. N[0] * N[1] * N[2] entries), while for node centered output the cell edge nodes are periodically replicated, resulting in (N[0]+1) * (N[1]+1) * (N[2]+1) exported entries.
When running with more than one MPI rank, this output currently supports only CELL mode. Each rank writes its local data to rank-suffixed files and rank 0 emits the global XDMF collection that references those blocks.
Output is written in a separate thread while the simulation continues. When "enable_hdf5" is set and the HDF5 library is not thread-safe, the output is written synchronously instead, as its HDF5 calls would otherwise race with other HDF5 users in the process (such as the Exodus output through netCDF).
Overview
Writes one or more buffers selected by "buffer" to an XDMF/HDF5 pair suitable for Paraview or Visit. Enable HDF5 with "enable_hdf5" and select cell or node data with "output_mode".
Example Input File Syntax
[TensorOutputs<<<{"href": "../../syntax/TensorOutputs/index.html"}>>>]
[xdmf2]
type = XDMFTensorOutput<<<{"description": "Output a tensor in XDMF format.", "href": "XDMFTensorOutput.html"}>>>
buffer<<<{"description": "The buffers to output"}>>> = 'rho u speed binary_media'
output_mode<<<{"description": "Output as cell or node data"}>>> = 'Cell Cell Cell Cell'
enable_hdf5<<<{"description": "Use HDF5 for binary data storage."}>>> = true
[]
[](examples/lbm/Karman-vortex/cylinder.i)Input Parameters
- enable_hdf5FalseUse HDF5 for binary data storage.
Default:False
C++ Type:bool
Controllable:No
Description:Use HDF5 for binary data storage.
- execute_onINITIAL TIMESTEP_ENDThe list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html.
Default:INITIAL TIMESTEP_END
C++ Type:ExecFlagEnum
Options:INITIAL, TIMESTEP_END
Controllable:No
Description:The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html.
- file_baseThe desired solution output name without an extension. If not provided, MOOSE sets it with Outputs/file_base when available. Otherwise, MOOSE uses input file name and this object name for a master input or uses master file_base, the subapp name and this object name for a subapp input to set it.
C++ Type:std::string
Controllable:No
Description:The desired solution output name without an extension. If not provided, MOOSE sets it with Outputs/file_base when available. Otherwise, MOOSE uses input file name and this object name for a master input or uses master file_base, the subapp name and this object name for a subapp input to set it.
- hdf5_compressionTrueEnable DEFLATE compression for HDF5 datasets.
Default:True
C++ Type:bool
Controllable:No
Description:Enable DEFLATE compression for HDF5 datasets.
- output_modeOutput as cell or node data
C++ Type:MultiMooseEnum
Options:CELL, NODE, OVERSIZED_NODAL
Controllable:No
Description:Output as cell or node data
- transposeTrueThe Paraview XDMF reader swaps x-y (x-z in 3d), so we transpose the tensors before we output to make the data look right in Paraview.
Default:True
C++ Type:bool
Controllable:No
Description:The Paraview XDMF reader swaps x-y (x-z in 3d), so we transpose the tensors before we output to make the data look right in Paraview.
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/Karman-vortex/cylinder.i)
- (examples/cahn_hilliard/cahnhilliard3.i)
- (test/tests/mechanics/mech.i)
- (test/tests/tensor_compute/rotating_grain_secant.i)
- (examples/lbm/Rayleigh-Benard/rayleigh-benard.i)
- (examples/swift_hohenberg/swifthohenberg.i)
- (test/tests/anisotropic_grain_growth/analytic_hull.i)
- (examples/anisotropic_grain_growth/five_grains_2D.i)
- (test/tests/kks/KKS_no_flux_bc.i)
- (examples/lbm/Phase-field/static_bubble_3d.i)
- (examples/anisotropic_grain_growth/spherical_grain_Ni.i)
- (test/tests/anisotropic_grain_growth/sigma3_circular_grain.i)
- (test/tests/tensor_compute/parallel_roundtrip_3d.i)
- (examples/phase_field_crystal/pfc_fcc_atomic_structure_3d.i)
- (examples/degeus_mechanics/mech.i)
- (test/tests/real_space/diffusion.i)
- (examples/lbm/rarefied_gas/channel.i)
- (examples/lbm/Phase-field/layered_poiseuille.i)
- (test/tests/kks/KKS_libtorch.i)
- (test/tests/tensor_compute/coupled_pf_mech.i)
- (examples/lbm/Pebbles/pebbles.i)
- (examples/lbm/Phase-field/spinodal_decomposition.i)
- (test/tests/real_space/laplacian.i)
- (test/tests/tensor_compute/coupled_pf_mech_secant.i)
- (examples/libtorch_kks/KKS_libtorch.i)
- (examples/phase_field_crystal/pfc_fcc_atomic_structure.i)
- (test/tests/typed_tensors/gradient.i)
- (test/tests/mechanics/mech3d.i)
- (examples/cahn_hilliard/cahnhilliard2.i)
- (benchmarks/02_oswald_ripening/2a.i)
- (benchmarks/02_oswald_ripening/2a_secant.i)
- (examples/lbm/Formula1-aerodynamics/f1.i)
- (examples/lbm/Phase-field/static_bubble.i)
- (test/tests/tensor_compute/group.i)
- (test/tests/tensor_compute/parallel.i)
- (test/tests/cahnhilliard/cahnhilliard.i)
- (examples/lbm/natural_convection_3D/convection.i)
- (examples/cahn_hilliard/cahnhilliard4.i)
- (test/tests/typed_tensors/gradient_vector.i)
- (examples/anisotropic_grain_growth/isotropic_five_grain.i)
- (benchmarks/01_spinodal_decomposition/1a_implicit_fd.i)
- (test/tests/real_space/laplacian3d.i)
- (test/tests/tensor_compute/smooth_rectangle.i)
- (examples/lbm/Phase-field/droplet_impact.i)
- (test/tests/tensor_compute/pfc_fcc.i)
enable_hdf5
Default:False
C++ Type:bool
Controllable:No
Description:Use HDF5 for binary data storage.
output_mode
C++ Type:MultiMooseEnum
Options:CELL, NODE, OVERSIZED_NODAL
Controllable:No
Description:Output as cell or node data
enable_hdf5
Default:False
C++ Type:bool
Controllable:No
Description:Use HDF5 for binary data storage.
buffer
C++ Type:std::vector<std::string>
Controllable:No
Description:The buffers to output
enable_hdf5
Default:False
C++ Type:bool
Controllable:No
Description:Use HDF5 for binary data storage.
output_mode
C++ Type:MultiMooseEnum
Options:CELL, NODE, OVERSIZED_NODAL
Controllable:No
Description:Output as cell or node data
(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
[]
[]
[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
[]
[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
residual_tensor = speed
[]
[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
[]
[]
(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
[]
[]
[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
[]
[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
residual_tensor = speed
[]
[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
[]
[]
(examples/cahn_hilliard/cahnhilliard3.i)
[Domain]
dim = 3
nx = 100
ny = 100
nz = 100
xmax = ${fparse pi*4}
ymax = ${fparse pi*4}
zmax = ${fparse pi*4}
device_names = 'cuda'
mesh_mode = DUMMY
[]
[TensorBuffers]
[c]
[]
[cbar]
[]
[mu]
[]
[mubar]
[]
[Mbarmubar]
[]
# constant tensors
[Mbar]
[]
[kappabarbar]
[]
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'c mu'
output_mode = 'Node Cell'
enable_hdf5 = true
[]
[]
[TensorComputes]
[Initialize]
[c]
type = RandomTensor
buffer = c
min = 0.44
max = 0.56
[]
[Mbar]
type = ReciprocalLaplacianFactor
factor = 0.2 # Mobility
buffer = Mbar
[]
[kappabarbar]
type = ReciprocalLaplacianSquareFactor
factor = -0.001 # kappa
buffer = kappabarbar
[]
[]
[Solve]
[mu]
type = ParsedCompute
buffer = mu
expression = '0.1*c^2*(c-1)^2'
derivatives = c
# expression = "0.4*c^3-0.6*c^2+0.2*c"
inputs = c
[]
[mubar]
type = ForwardFFT
buffer = mubar
input = mu
[]
[Mbarmubar]
type = ParsedCompute
buffer = Mbarmubar
expression = 'Mbar*mubar'
inputs = 'Mbar mubar'
[]
[cbar]
type = ForwardFFT
buffer = cbar
input = c
[]
[]
[]
[TensorTimeIntegrators]
[c]
type = FFTSemiImplicit
buffer = c
reciprocal_buffer = cbar
linear_reciprocal = kappabarbar
nonlinear_reciprocal = Mbarmubar
[]
[]
[Postprocessors]
[min_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[C]
type = TensorIntegralPostprocessor
buffer = c
execute_on = 'TIMESTEP_END'
[]
[cavg]
type = TensorAveragePostprocessor
buffer = c
execute_on = 'TIMESTEP_END'
[]
[]
[Problem]
type = TensorProblem
spectral_solve_substeps = 1000
[]
[Executioner]
type = Transient
num_steps = 20
[TimeStepper]
type = IterationAdaptiveDT
growth_factor = 1.8
dt = 0.1
[]
dtmax = 1000
[]
[Outputs]
csv = true
perf_graph = true
execute_on = 'TIMESTEP_END'
[]
(test/tests/mechanics/mech.i)
[Domain]
dim = 2
nx = 32
ny = 32
xmax = ${fparse 2*pi}
ymax = ${fparse 2*pi}
zmax = ${fparse 2*pi}
mesh_mode = DUMMY
parallel_mode = FFT_AUTO
[]
[TensorComputes]
[Initialize]
[phase]
type = ParsedCompute
expression = '(cos(x)/2+0.5)^1*(cos(y)/2+0.5)^1*(cos(z)/2+0.5)^1'
extra_symbols = true
buffer = phase
[]
[K]
type = ParsedCompute
buffer = K
expression = '(1-phase)*Ka + phase*Kb'
inputs = phase
constant_names = 'Ka Kb'
constant_expressions = '1 10'
[]
[mu]
type = ParsedCompute
buffer = mu
expression = '(1-phase)*mua + phase*mub'
inputs = phase
constant_names = 'mua mub'
constant_expressions = '0.5 5'
[]
[Finit]
type = RankTwoIdentity
buffer = F
[]
[]
[Solve]
[hyper_elasticity]
type = HyperElasticIsotropic
buffer = stress
F = Fnew
K = K
mu = mu
[]
[root]
[applied_strain]
type = MacroscopicShearTensor
buffer = applied_strain
[]
[mech]
type = FFTMechanics
buffer = Fnew
F = F
K = K
mu = mu
l_max_its = 40
l_tol = 1e-5
nl_rel_tol = 2e-4
nl_abs_tol = 2e-3
constitutive_model = hyper_elasticity
stress = stress
applied_macroscopic_strain = applied_strain
[]
[]
[]
[Postprocess]
[displacements]
type = ComputeDisplacements
buffer = disp
F = F
[]
[vonmises]
type = ComputeVonMisesStress
buffer = sV
[]
[]
[]
[TensorSolver]
# no variables are integrated by this solver (FFTMechanics performs a steady state mechanics solve)
type = ForwardEulerSolver
root_compute = root
# deformation tensor is just forwarded Fnew -> F
forward_buffer = F
forward_buffer_new = Fnew
substeps = 3
[]
[TensorOutputs]
active = 'parallel' # serial
[serial]
type = XDMFTensorOutput
file_base = mech_serial
buffer = 'disp sV F phase'
output_mode = 'OVERSIZED_NODAL CELL CELL NODE'
enable_hdf5 = true
execute_on = 'TIMESTEP_END'
[]
[parallel]
# NODE and OVERSIZED_NODAL are not yet available in parallel output
type = XDMFTensorOutput
file_base = mech_parallel
buffer = 'sV F phase'
output_mode = 'CELL CELL CELL'
enable_hdf5 = true
execute_on = 'TIMESTEP_END'
[]
[]
[Executioner]
type = Transient
num_steps = 3
dt = 0.02
[]
(test/tests/tensor_compute/rotating_grain_secant.i)
w=6
[Domain]
dim = 2
nx = 40
ny = 40
xmax = ${fparse w*pi*2}
ymax = ${fparse w*pi*2/sin(pi/3)}
mesh_mode = DOMAIN
[]
[AuxVariables]
[phi]
[]
[]
[Outputs]
exodus = false
[]
[TensorBuffers]
[psi]
map_to_aux_variable = phi
[]
[psibar]
[]
[psi3]
[]
[psi3bar]
[]
# constant tensors
[linear]
[]
[]
crystal = '-(sin(sin(a)*y/2+cos(a)*x/2)^2 + sin(sin(a+1/3*pi)*y/2+cos(a+1/3*pi)*x/2)^2 + sin(sin(a-1/3*pi)*y/2+cos(a-1/3*pi)*x/2)^2 - 1.5)*0.25'
[Functions]
[grain1]
type = ParsedFunction
expression = 'a := 0; ${crystal}'
[]
[grain2]
type = ParsedFunction
expression = 'a := 0.95; ${crystal}'
[]
[domain]
type = ParsedFunction
expression = 'r := (x-${w}*pi)^2+(y-${w}*pi)^2; if(r<(${w}*2/3*pi)^2, grain2, grain1)'
symbol_names = 'grain1 grain2'
symbol_values = 'grain1 grain2'
[]
[]
[TensorComputes]
[Initialize]
[psi]
type = MooseFunctionTensor
buffer = psi
function = domain
[]
[linear]
type = SwiftHohenbergLinear
buffer = linear
alpha = 1
r = 0.025
[]
[]
[Solve]
[psi3]
type = ParsedCompute
buffer = psi3
expression = "0.20*psi^2-psi^3"
inputs = psi
[]
[psibar]
type = ForwardFFT
buffer = psibar
input = psi
[]
[psi3bar]
type = ForwardFFT
buffer = psi3bar
input = psi3
[]
[]
[]
[TensorSolver]
type = SecantSolver
buffer = psi
substeps = 3
reciprocal_buffer = psibar
linear_reciprocal = linear
nonlinear_reciprocal = psi3bar
[]
[Problem]
type = TensorProblem
[]
[Executioner]
type = Transient
num_steps = 10
[TimeStepper]
type = TensorSolveIterationAdaptiveDT
dt = 1
max_iterations = 400
min_iterations = 100
growth_factor = 1.4
cutback_factor = 0.9
[]
dtmax = 500
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'psi'
enable_hdf5 = true
# Do not transpose output to avoid regolding the test. In practice the default
# of transpose = true should always be used
transpose = false
[]
[]
(examples/lbm/Rayleigh-Benard/rayleigh-benard.i)
Nx = 400
Ny = 200
TH = 1.1
TC = 1.0
frequency = '${Nx}/10.0'
amplitude = '${Ny}/100.0'
[Domain]
dim = 2
nx = '${Nx}'
ny = '${Ny}'
xmax = '${Nx}'
ymax = '${Ny}'
device_names='cpu'
mesh_mode=DUMMY
parallel_mode = REAL_SPACE
periodic_directions = 'X Y'
[]
[Stencil]
[d2q9]
type = LBMD2Q9
[]
[]
!include buffers.i
[TensorComputes]
#### Initialzie ####
[Initialize]
[density]
type = LBMConstantTensor
buffer = density
constants = 'rho0'
[]
[velocity]
type = LBMConstantTensor
buffer = velocity
constants = '0.0 0.0'
[]
[temperature]
type = ParsedCompute
buffer = T
expression = 'a:=abs(y - sin(x / (${frequency} * pi)) * ${amplitude}) + y - sin(x / (${frequency} * pi)) * ${amplitude};
b:= a / (a + 1.0e-14);
${TC} * b - b *${TH} + ${TH}'
extra_symbols = true
[]
[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
[]
[]
!include solve.i
!include boundary.i
[]
[TensorSolver]
type = LBMStream
buffer = 'f g'
f_old = 'fpc gpc'
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'T velocity density'
output_mode = 'Cell Cell Cell'
enable_hdf5 = true
[]
[]
[Problem]
type = LatticeBoltzmannProblem
scalar_constant_names = 'rho0 T_0 T_C T_H tau_f tau_T g'
scalar_constant_values = '1.0 1.05 1.0 1.1 0.7 0.7 0.0001'
substeps = 100
print_debug_output = true
residual_tensor = speed
[]
[Executioner]
type = Transient
num_steps = 10000
[]
(examples/swift_hohenberg/swifthohenberg.i)
#
# Solve a simple Swift-Hohenberg crystal phase field problem. The initial condition is
# a random melt, which crystalizes out.
#
[Domain]
dim = 2
nx = 400
ny = 400
xmax = ${fparse pi*2*30}
ymax = ${fparse pi*2*30}
device_names = 'cuda'
mesh_mode = DOMAIN
[]
[TensorBuffers]
[psi]
map_to_aux_variable = psi
[]
[psibar]
[]
[psi3]
[]
[psi3bar]
[]
# constant tensors
[linear]
[]
# output
[filter]
map_to_aux_variable = filter
[]
[filterbar]
[]
[]
[AuxVariables]
[psi]
[]
[filter]
[]
[]
[TensorComputes]
[Initialize]
[psi]
type = RandomTensor
buffer = psi
min = 0
max = 0.07
[]
[linear]
type = SwiftHohenbergLinear
buffer = linear
alpha = 1
r = 0.025
[]
[]
[Solve]
[psi3]
type = ParsedCompute
buffer = psi3
expression = "0.20*psi^2-psi^3"
inputs = psi
[]
[psibar]
type = ForwardFFT
buffer = psibar
input = psi
[]
[psi3bar]
type = ForwardFFT
buffer = psi3bar
input = psi3
[]
[]
[Postprocess]
[low_pass]
type = ParsedCompute
buffer = filterbar
extra_symbols = true
expression = 'psibar * exp(-k2*10)'
inputs = psibar
[]
[filter]
type = InverseFFT
buffer = filter
input = filterbar
[]
[]
[]
[TensorTimeIntegrators]
[c]
type = FFTSemiImplicit
buffer = psi
reciprocal_buffer = psibar
linear_reciprocal = linear
nonlinear_reciprocal = psi3bar
[]
[]
[Problem]
type = TensorProblem
spectral_solve_substeps = 1000
[]
[Executioner]
type = Transient
num_steps = 300
[TimeStepper]
type = IterationAdaptiveDT
growth_factor = 1.2
dt = 10
[]
dtmax = 1000
[]
[Postprocessors]
[min_psi]
type = TensorExtremeValuePostprocessor
buffer = psi
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_psi]
type = TensorExtremeValuePostprocessor
buffer = psi
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[Psi]
type = TensorIntegralPostprocessor
buffer = psi
[]
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'psi'
output_mode = 'Node'
enable_hdf5 = true
[]
[]
[Outputs]
exodus = true
perf_graph = true
execute_on = 'TIMESTEP_END'
[]
(test/tests/anisotropic_grain_growth/analytic_hull.i)
interface_width = 0.8
r0 = 4
gbe_max = ${units 0.8 J/m^2}
kappa = ${fparse 0.75 * gbe_max * interface_width}
L = 1
g_eta_expr = '(eta^2*(1-eta^2)^2)'
[Domain]
dim = 2
nx = 32
ny = 32
xmax = 16
ymax = 16
mesh_mode = DUMMY
floating_precision = SINGLE
[]
[TensorComputes]
[Initialize]
[eta]
type = ParsedCompute
buffer = 'eta'
extra_symbols = 'true'
expression = 'radius:=sqrt((x-${Domain/xmax}/2)^2+(y-${Domain/ymax}/2)^2);0.5 - 0.5*tanh(2*(radius-${r0})/${interface_width})'
[]
[kappa_linear_term]
type = ReciprocalLaplacianFactor
buffer = kappa_linear_term
factor = ${fparse ${L} * ${kappa}}
[]
[smooth]
type = DeAliasingTensor
method = HOULI
buffer = smooth
[]
[]
[Solve]
[grad_eta_vector]
type = GradientVector
buffer = grad_eta
input = eta
[]
[gb_energy]
type = AnisotropicGBEnergy
buffer = gb_energy
gb_gradient_buffer = grad_eta
dsigma_dn = dsigma_dn
libtorch_model_file = 'marlin:anisotropic_gb_energy/analytic_sigma.pt'
interface_width = ${interface_width}
chunk_size = 256
[]
[bulk_driving_force]
type = ParsedCompute
buffer = 'bulk_driving_force'
expression = '6 * gb_energy * ${g_eta_expr} / ${interface_width}'
inputs = 'mu eta gb_energy'
derivatives = 'eta'
[]
[dmu_dn]
type = ParsedCompute
buffer = 'dmu_dn'
expression = '6 * dsigma_dn / ${interface_width}'
inputs = 'dsigma_dn'
[]
[g]
type = ParsedCompute
buffer = 'g'
expression = '${g_eta_expr}'
inputs = 'eta'
[]
[anisotropy_torque]
type = AnisotropyTorque
buffer = anisotropy_torque
dmu_dn = dmu_dn
g = g
[]
[total_mu_driving_force]
type = ParsedCompute
buffer = 'total_mu_driving_force'
expression = '-${L}*(bulk_driving_force - anisotropy_torque)'
inputs = 'bulk_driving_force anisotropy_torque'
[]
[nonlinear]
type = ForwardFFT
buffer = 'NL'
input = 'total_mu_driving_force'
[]
[smooth_nonlinear]
type = ParsedCompute
buffer = 'smooth_NL'
expression = 'smooth * NL'
inputs = 'smooth NL'
[]
[etabar]
type = ForwardFFT
buffer = etabar
input = eta
[]
[gradient_energy]
type = GradientEnergyDensity
buffer = gradient_energy
gradient_buffers = 'grad_eta'
kappa = ${kappa}
[]
[interface_energy]
type = ParsedCompute
buffer = interface_energy
expression = 'g * 6 * gb_energy / ${interface_width} + gradient_energy'
inputs = 'g gb_energy gradient_energy'
[]
[]
[]
[TensorSolver]
type = AdamsBashforthMoulton
buffer = 'eta'
reciprocal_buffer = 'etabar'
linear_reciprocal = 'kappa_linear_term'
nonlinear_reciprocal = 'smooth_NL'
substeps = 1e2
predictor_order = 1
corrector_order = 1
corrector_steps = 1
[]
[Postprocessors]
[total_gb_energy]
type = TensorIntegralPostprocessor
buffer = interface_energy
[]
[]
[TensorOutputs]
active = ''
[xdmf]
type = XDMFTensorOutput
buffer = 'eta gb_energy'
output_mode = 'NODE NODE'
enable_hdf5 = true
transpose = false
[]
[]
[Executioner]
type = Transient
dt = 0.1
num_steps = 2
[]
[Outputs]
csv = true
execute_on = 'INITIAL TIMESTEP_END'
[]
(examples/anisotropic_grain_growth/five_grains_2D.i)
interface_width = 1.6
side = 80
gbe_max = '${units 1.24 J/m^2}'
L = '${fparse 1.0 * 1.6 / ${interface_width} }'
# ------------------------------------------------------------------
# Periodic Voronoi-like IC for 5 grains.
#
# Each grain is seeded at a point (xi, yi). We use a smooth *periodic*
# distance metric based on sine functions,
# d_i^2 = (side/pi)^2 * [ sin(pi*(x-xi)/side)^2 + sin(pi*(y-yi)/side)^2 ]
# which reduces to the Euclidean distance near the seed and wraps
# correctly across the periodic boundaries.
#
# The order parameters are then a softmax over the negative distances,
# gr_i = exp(-a*d_i) / sum_j exp(-a*d_j), a = 4/interface_width
# For any two competing grains this reduces exactly to the usual
# 0.5*(1 - tanh(2*(d_i - d_j)/interface_width)) profile, so the
# interfaces come out with the correct width, the partition of unity
# sum_i gr_i = 1 holds identically, and the cell structure is the
# (diffuse) Voronoi tessellation of the seeds under the periodic metric.
# ------------------------------------------------------------------
aexp = '${fparse 4.0 / ${interface_width} }'
lp = '${fparse ${side} / pi }'
kp = '${fparse pi / ${side} }'
# seed points (spread out over the periodic 40x40 domain)
d0 = 'sqrt((${lp}*sin(${kp}*(x-8)))^2 + (${lp}*sin(${kp}*(y-6)))^2)'
d1 = 'sqrt((${lp}*sin(${kp}*(x-28)))^2 + (${lp}*sin(${kp}*(y-10)))^2)'
d2 = 'sqrt((${lp}*sin(${kp}*(x-20)))^2 + (${lp}*sin(${kp}*(y-24)))^2)'
d3 = 'sqrt((${lp}*sin(${kp}*(x-4)))^2 + (${lp}*sin(${kp}*(y-28)))^2)'
d4 = 'sqrt((${lp}*sin(${kp}*(x-34)))^2 + (${lp}*sin(${kp}*(y-33)))^2)'
g_gamma0 = '${fparse sqrt(2) / 3 }' # g(gamma=1.5)
f0_gamma0 = ${fparse 1/8}
kappa = '${fparse gbe_max * interface_width * sqrt(f0_gamma0) / g_gamma0 }'
mu = '${fparse gbe_max / (g_gamma0 * interface_width * sqrt(f0_gamma0) ) }'
g_to_gamma_func = '(1.0/(-3.0944 * g^8 -1.8169*g^6 + 10.323 * g^4 - 8.1819*g^2 + 2.0033))'
dgamma_dsigma = '(g * (24.7552*g^6 + 10.9014*g^4 - 41.292*g^2 + 16.3638) / (-3.0944*g^8 - 1.8169*g^6 + 10.323*g^4 - 8.1819*g^2 + 2.0033)^2 / sqrt(${kappa} * ${mu}) )'
f0 = '((gr0^4/4 - gr0^2/2) + (gr1^4/4 - gr1^2/2) + (gr2^4/4 - gr2^2/2) + (gr3^4/4 - gr3^2/2) + (gr4^4/4 - gr4^2/2) + gamma_01*gr0^2*gr1^2 + gamma_02*gr0^2*gr2^2 + gamma_03*gr0^2*gr3^2 + gamma_04*gr0^2*gr4^2 + gamma_12*gr1^2*gr2^2 + gamma_13*gr1^2*gr3^2 + gamma_14*gr1^2*gr4^2 + gamma_23*gr2^2*gr3^2 + gamma_24*gr2^2*gr4^2 + gamma_34*gr3^2*gr4^2 + 0.25)'
[Domain]
dim = 2
nx = 200
ny = 200
xmax = ${side}
ymax = ${side}
mesh_mode = DUMMY
[]
[TensorComputes]
[Initialize]
# unnormalized softmax weights
[c0]
type = ParsedCompute
buffer = c0
expression = 'exp(-${aexp}*${d0})'
extra_symbols = true
[]
[c1]
type = ParsedCompute
buffer = c1
expression = 'exp(-${aexp}*${d1})'
extra_symbols = true
[]
[c2]
type = ParsedCompute
buffer = c2
expression = 'exp(-${aexp}*${d2})'
extra_symbols = true
[]
[c3]
type = ParsedCompute
buffer = c3
expression = 'exp(-${aexp}*${d3})'
extra_symbols = true
[]
[c4]
type = ParsedCompute
buffer = c4
expression = 'exp(-${aexp}*${d4})'
extra_symbols = true
[]
# normalized order parameters (partition of unity by construction)
[gr0]
type = ParsedCompute
buffer = gr0
expression = 'c0/(c0+c1+c2+c3+c4)'
inputs = 'c0 c1 c2 c3 c4'
[]
[gr1]
type = ParsedCompute
buffer = gr1
expression = 'c1/(c0+c1+c2+c3+c4)'
inputs = 'c0 c1 c2 c3 c4'
[]
[gr2]
type = ParsedCompute
buffer = gr2
expression = 'c2/(c0+c1+c2+c3+c4)'
inputs = 'c0 c1 c2 c3 c4'
[]
[gr3]
type = ParsedCompute
buffer = gr3
expression = 'c3/(c0+c1+c2+c3+c4)'
inputs = 'c0 c1 c2 c3 c4'
[]
[gr4]
type = ParsedCompute
buffer = gr4
expression = '1 - gr0 - gr1 - gr2 - gr3'
inputs = 'gr0 gr1 gr2 gr3'
[]
[L_kappa_laplacian]
type = ReciprocalLaplacianFactor
buffer = 'L_kappa_laplacian'
factor = '${fparse ${L} * ${kappa}}'
[]
[kappa_laplacian]
type = ReciprocalLaplacianFactor
buffer = kappa_laplacian
factor = '${kappa}'
[]
[smooth]
type = DeAliasingTensor
method = HOULI
buffer = smooth
[]
[]
[Solve]
[bnds]
type = ParsedCompute
buffer = 'bnds'
inputs = 'gr0 gr1 gr2 gr3 gr4'
expression = 'gr0^2 + gr1^2 + gr2^2 + gr3^2 + gr4^2'
[]
## FFT
[gr0_bar]
type = ForwardFFT
buffer = gr0_bar
input = gr0
[]
[gr1_bar]
type = ForwardFFT
buffer = gr1_bar
input = gr1
[]
[gr2_bar]
type = ForwardFFT
buffer = gr2_bar
input = gr2
[]
[gr3_bar]
type = ForwardFFT
buffer = gr3_bar
input = gr3
[]
[gr4_bar]
type = ForwardFFT
buffer = gr4_bar
input = gr4
[]
[grad_gr0]
type = GradientVector
buffer = grad_gr0
input = gr0
[]
[grad_gr1]
type = GradientVector
buffer = grad_gr1
input = gr1
[]
[grad_gr2]
type = GradientVector
buffer = grad_gr2
input = gr2
[]
[grad_gr3]
type = GradientVector
buffer = grad_gr3
input = gr3
[]
[grad_gr4]
type = GradientVector
buffer = grad_gr4
input = gr4
[]
## Pairwise anisotropic GB energies (10 pairs)
# NOTE: the per-pair libtorch misorientation models (Ri_Rj.pt) will be
# developed separately; placeholder file names below follow the
# existing convention and location.
[sigma_gr0_gr1]
type = PairwiseAnisotropicGBEnergy
buffer = 'sigma_gr0_gr1'
dsigma_dgrad_grain1 = 'dsigma_gr0_gr1_dgrad_gr0'
dsigma_dgrad_grain2 = 'dsigma_gr0_gr1_dgrad_gr1'
grad_grain1_buffer = 'grad_gr0'
grad_grain2_buffer = 'grad_gr1'
interface_width = '${interface_width}'
libtorch_model_file = 'marlin:anisotropic_gb_energy/Ni_R1_R2.pt'
chunk_size = 1e5
[]
[sigma_gr0_gr2]
type = PairwiseAnisotropicGBEnergy
buffer = 'sigma_gr0_gr2'
dsigma_dgrad_grain1 = 'dsigma_gr0_gr2_dgrad_gr0'
dsigma_dgrad_grain2 = 'dsigma_gr0_gr2_dgrad_gr2'
grad_grain1_buffer = 'grad_gr0'
grad_grain2_buffer = 'grad_gr2'
interface_width = '${interface_width}'
libtorch_model_file = 'marlin:anisotropic_gb_energy/Ni_R1_R3.pt'
chunk_size = 1e5
[]
[sigma_gr0_gr3]
type = PairwiseAnisotropicGBEnergy
buffer = 'sigma_gr0_gr3'
dsigma_dgrad_grain1 = 'dsigma_gr0_gr3_dgrad_gr0'
dsigma_dgrad_grain2 = 'dsigma_gr0_gr3_dgrad_gr3'
grad_grain1_buffer = 'grad_gr0'
grad_grain2_buffer = 'grad_gr3'
interface_width = '${interface_width}'
libtorch_model_file = 'marlin:anisotropic_gb_energy/Ni_R1_R4.pt'
chunk_size = 1e5
[]
[sigma_gr0_gr4]
type = PairwiseAnisotropicGBEnergy
buffer = 'sigma_gr0_gr4'
dsigma_dgrad_grain1 = 'dsigma_gr0_gr4_dgrad_gr0'
dsigma_dgrad_grain2 = 'dsigma_gr0_gr4_dgrad_gr4'
grad_grain1_buffer = 'grad_gr0'
grad_grain2_buffer = 'grad_gr4'
interface_width = '${interface_width}'
libtorch_model_file = 'marlin:anisotropic_gb_energy/Ni_R1_R5.pt'
chunk_size = 1e5
[]
[sigma_gr1_gr2]
type = PairwiseAnisotropicGBEnergy
buffer = 'sigma_gr1_gr2'
dsigma_dgrad_grain1 = 'dsigma_gr1_gr2_dgrad_gr1'
dsigma_dgrad_grain2 = 'dsigma_gr1_gr2_dgrad_gr2'
grad_grain1_buffer = 'grad_gr1'
grad_grain2_buffer = 'grad_gr2'
interface_width = '${interface_width}'
libtorch_model_file = 'marlin:anisotropic_gb_energy/Ni_R2_R3.pt'
chunk_size = 1e5
[]
[sigma_gr1_gr3]
type = PairwiseAnisotropicGBEnergy
buffer = 'sigma_gr1_gr3'
dsigma_dgrad_grain1 = 'dsigma_gr1_gr3_dgrad_gr1'
dsigma_dgrad_grain2 = 'dsigma_gr1_gr3_dgrad_gr3'
grad_grain1_buffer = 'grad_gr1'
grad_grain2_buffer = 'grad_gr3'
interface_width = '${interface_width}'
libtorch_model_file = 'marlin:anisotropic_gb_energy/Ni_R2_R4.pt'
chunk_size = 1e5
[]
[sigma_gr1_gr4]
type = PairwiseAnisotropicGBEnergy
buffer = 'sigma_gr1_gr4'
dsigma_dgrad_grain1 = 'dsigma_gr1_gr4_dgrad_gr1'
dsigma_dgrad_grain2 = 'dsigma_gr1_gr4_dgrad_gr4'
grad_grain1_buffer = 'grad_gr1'
grad_grain2_buffer = 'grad_gr4'
interface_width = '${interface_width}'
libtorch_model_file = 'marlin:anisotropic_gb_energy/Ni_R2_R5.pt'
chunk_size = 1e5
[]
[sigma_gr2_gr3]
type = PairwiseAnisotropicGBEnergy
buffer = 'sigma_gr2_gr3'
dsigma_dgrad_grain1 = 'dsigma_gr2_gr3_dgrad_gr2'
dsigma_dgrad_grain2 = 'dsigma_gr2_gr3_dgrad_gr3'
grad_grain1_buffer = 'grad_gr2'
grad_grain2_buffer = 'grad_gr3'
interface_width = '${interface_width}'
libtorch_model_file = 'marlin:anisotropic_gb_energy/Ni_R3_R4.pt'
chunk_size = 1e5
[]
[sigma_gr2_gr4]
type = PairwiseAnisotropicGBEnergy
buffer = 'sigma_gr2_gr4'
dsigma_dgrad_grain1 = 'dsigma_gr2_gr4_dgrad_gr2'
dsigma_dgrad_grain2 = 'dsigma_gr2_gr4_dgrad_gr4'
grad_grain1_buffer = 'grad_gr2'
grad_grain2_buffer = 'grad_gr4'
interface_width = '${interface_width}'
libtorch_model_file = 'marlin:anisotropic_gb_energy/Ni_R3_R5.pt'
chunk_size = 1e5
[]
[sigma_gr3_gr4]
type = PairwiseAnisotropicGBEnergy
buffer = 'sigma_gr3_gr4'
dsigma_dgrad_grain1 = 'dsigma_gr3_gr4_dgrad_gr3'
dsigma_dgrad_grain2 = 'dsigma_gr3_gr4_dgrad_gr4'
grad_grain1_buffer = 'grad_gr3'
grad_grain2_buffer = 'grad_gr4'
interface_width = '${interface_width}'
libtorch_model_file = 'marlin:anisotropic_gb_energy/Ni_R4_R5.pt'
chunk_size = 1e5
[]
## gamma_ij(sigma_ij)
[gamma_01]
type = ParsedCompute
buffer = 'gamma_01'
expression = 'g:=sigma_gr0_gr1/sqrt(${kappa}*${mu});${g_to_gamma_func}'
inputs = 'sigma_gr0_gr1'
[]
[gamma_02]
type = ParsedCompute
buffer = 'gamma_02'
expression = 'g:=sigma_gr0_gr2/sqrt(${kappa}*${mu});${g_to_gamma_func}'
inputs = 'sigma_gr0_gr2'
[]
[gamma_03]
type = ParsedCompute
buffer = 'gamma_03'
expression = 'g:=sigma_gr0_gr3/sqrt(${kappa}*${mu});${g_to_gamma_func}'
inputs = 'sigma_gr0_gr3'
[]
[gamma_04]
type = ParsedCompute
buffer = 'gamma_04'
expression = 'g:=sigma_gr0_gr4/sqrt(${kappa}*${mu});${g_to_gamma_func}'
inputs = 'sigma_gr0_gr4'
[]
[gamma_12]
type = ParsedCompute
buffer = 'gamma_12'
expression = 'g:=sigma_gr1_gr2/sqrt(${kappa}*${mu});${g_to_gamma_func}'
inputs = 'sigma_gr1_gr2'
[]
[gamma_13]
type = ParsedCompute
buffer = 'gamma_13'
expression = 'g:=sigma_gr1_gr3/sqrt(${kappa}*${mu});${g_to_gamma_func}'
inputs = 'sigma_gr1_gr3'
[]
[gamma_14]
type = ParsedCompute
buffer = 'gamma_14'
expression = 'g:=sigma_gr1_gr4/sqrt(${kappa}*${mu});${g_to_gamma_func}'
inputs = 'sigma_gr1_gr4'
[]
[gamma_23]
type = ParsedCompute
buffer = 'gamma_23'
expression = 'g:=sigma_gr2_gr3/sqrt(${kappa}*${mu});${g_to_gamma_func}'
inputs = 'sigma_gr2_gr3'
[]
[gamma_24]
type = ParsedCompute
buffer = 'gamma_24'
expression = 'g:=sigma_gr2_gr4/sqrt(${kappa}*${mu});${g_to_gamma_func}'
inputs = 'sigma_gr2_gr4'
[]
[gamma_34]
type = ParsedCompute
buffer = 'gamma_34'
expression = 'g:=sigma_gr3_gr4/sqrt(${kappa}*${mu});${g_to_gamma_func}'
inputs = 'sigma_gr3_gr4'
[]
## torque prefactors dgamma/dsigma * gr_i^2 * gr_j^2
[coeff_01]
type = ParsedCompute
buffer = 'coeff_01'
expression = 'g:=sigma_gr0_gr1/sqrt(${kappa}*${mu});${dgamma_dsigma} * gr0^2 *gr1^2'
inputs = 'sigma_gr0_gr1 gr0 gr1'
[]
[coeff_02]
type = ParsedCompute
buffer = 'coeff_02'
expression = 'g:=sigma_gr0_gr2/sqrt(${kappa}*${mu});${dgamma_dsigma} * gr0^2 *gr2^2'
inputs = 'sigma_gr0_gr2 gr0 gr2'
[]
[coeff_03]
type = ParsedCompute
buffer = 'coeff_03'
expression = 'g:=sigma_gr0_gr3/sqrt(${kappa}*${mu});${dgamma_dsigma} * gr0^2 *gr3^2'
inputs = 'sigma_gr0_gr3 gr0 gr3'
[]
[coeff_04]
type = ParsedCompute
buffer = 'coeff_04'
expression = 'g:=sigma_gr0_gr4/sqrt(${kappa}*${mu});${dgamma_dsigma} * gr0^2 *gr4^2'
inputs = 'sigma_gr0_gr4 gr0 gr4'
[]
[coeff_12]
type = ParsedCompute
buffer = 'coeff_12'
expression = 'g:=sigma_gr1_gr2/sqrt(${kappa}*${mu});${dgamma_dsigma} * gr1^2 *gr2^2'
inputs = 'sigma_gr1_gr2 gr1 gr2'
[]
[coeff_13]
type = ParsedCompute
buffer = 'coeff_13'
expression = 'g:=sigma_gr1_gr3/sqrt(${kappa}*${mu});${dgamma_dsigma} * gr1^2 *gr3^2'
inputs = 'sigma_gr1_gr3 gr1 gr3'
[]
[coeff_14]
type = ParsedCompute
buffer = 'coeff_14'
expression = 'g:=sigma_gr1_gr4/sqrt(${kappa}*${mu});${dgamma_dsigma} * gr1^2 *gr4^2'
inputs = 'sigma_gr1_gr4 gr1 gr4'
[]
[coeff_23]
type = ParsedCompute
buffer = 'coeff_23'
expression = 'g:=sigma_gr2_gr3/sqrt(${kappa}*${mu});${dgamma_dsigma} * gr2^2 *gr3^2'
inputs = 'sigma_gr2_gr3 gr2 gr3'
[]
[coeff_24]
type = ParsedCompute
buffer = 'coeff_24'
expression = 'g:=sigma_gr2_gr4/sqrt(${kappa}*${mu});${dgamma_dsigma} * gr2^2 *gr4^2'
inputs = 'sigma_gr2_gr4 gr2 gr4'
[]
[coeff_34]
type = ParsedCompute
buffer = 'coeff_34'
expression = 'g:=sigma_gr3_gr4/sqrt(${kappa}*${mu});${dgamma_dsigma} * gr3^2 *gr4^2'
inputs = 'sigma_gr3_gr4 gr3 gr4'
[]
## anisotropy torques (2 per pair)
[torque_01_gr0]
type = AnisotropyTorque
buffer = torque_01_gr0
dmu_dn = dsigma_gr0_gr1_dgrad_gr0
g = 'coeff_01'
[]
[torque_01_gr1]
type = AnisotropyTorque
buffer = torque_01_gr1
dmu_dn = dsigma_gr0_gr1_dgrad_gr1
g = 'coeff_01'
[]
[torque_02_gr0]
type = AnisotropyTorque
buffer = torque_02_gr0
dmu_dn = dsigma_gr0_gr2_dgrad_gr0
g = 'coeff_02'
[]
[torque_02_gr2]
type = AnisotropyTorque
buffer = torque_02_gr2
dmu_dn = dsigma_gr0_gr2_dgrad_gr2
g = 'coeff_02'
[]
[torque_03_gr0]
type = AnisotropyTorque
buffer = torque_03_gr0
dmu_dn = dsigma_gr0_gr3_dgrad_gr0
g = 'coeff_03'
[]
[torque_03_gr3]
type = AnisotropyTorque
buffer = torque_03_gr3
dmu_dn = dsigma_gr0_gr3_dgrad_gr3
g = 'coeff_03'
[]
[torque_04_gr0]
type = AnisotropyTorque
buffer = torque_04_gr0
dmu_dn = dsigma_gr0_gr4_dgrad_gr0
g = 'coeff_04'
[]
[torque_04_gr4]
type = AnisotropyTorque
buffer = torque_04_gr4
dmu_dn = dsigma_gr0_gr4_dgrad_gr4
g = 'coeff_04'
[]
[torque_12_gr1]
type = AnisotropyTorque
buffer = torque_12_gr1
dmu_dn = dsigma_gr1_gr2_dgrad_gr1
g = 'coeff_12'
[]
[torque_12_gr2]
type = AnisotropyTorque
buffer = torque_12_gr2
dmu_dn = dsigma_gr1_gr2_dgrad_gr2
g = 'coeff_12'
[]
[torque_13_gr1]
type = AnisotropyTorque
buffer = torque_13_gr1
dmu_dn = dsigma_gr1_gr3_dgrad_gr1
g = 'coeff_13'
[]
[torque_13_gr3]
type = AnisotropyTorque
buffer = torque_13_gr3
dmu_dn = dsigma_gr1_gr3_dgrad_gr3
g = 'coeff_13'
[]
[torque_14_gr1]
type = AnisotropyTorque
buffer = torque_14_gr1
dmu_dn = dsigma_gr1_gr4_dgrad_gr1
g = 'coeff_14'
[]
[torque_14_gr4]
type = AnisotropyTorque
buffer = torque_14_gr4
dmu_dn = dsigma_gr1_gr4_dgrad_gr4
g = 'coeff_14'
[]
[torque_23_gr2]
type = AnisotropyTorque
buffer = torque_23_gr2
dmu_dn = dsigma_gr2_gr3_dgrad_gr2
g = 'coeff_23'
[]
[torque_23_gr3]
type = AnisotropyTorque
buffer = torque_23_gr3
dmu_dn = dsigma_gr2_gr3_dgrad_gr3
g = 'coeff_23'
[]
[torque_24_gr2]
type = AnisotropyTorque
buffer = torque_24_gr2
dmu_dn = dsigma_gr2_gr4_dgrad_gr2
g = 'coeff_24'
[]
[torque_24_gr4]
type = AnisotropyTorque
buffer = torque_24_gr4
dmu_dn = dsigma_gr2_gr4_dgrad_gr4
g = 'coeff_24'
[]
[torque_34_gr3]
type = AnisotropyTorque
buffer = torque_34_gr3
dmu_dn = dsigma_gr3_gr4_dgrad_gr3
g = 'coeff_34'
[]
[torque_34_gr4]
type = AnisotropyTorque
buffer = torque_34_gr4
dmu_dn = dsigma_gr3_gr4_dgrad_gr4
g = 'coeff_34'
[]
## interpolated local GB energy and nonlinear mobility
[sigma]
type = ParsedCompute
buffer = 'sigma'
expression = '(gr0^2*gr1^2*sigma_gr0_gr1 + gr0^2*gr2^2*sigma_gr0_gr2 + gr0^2*gr3^2*sigma_gr0_gr3 + gr0^2*gr4^2*sigma_gr0_gr4 + gr1^2*gr2^2*sigma_gr1_gr2 + gr1^2*gr3^2*sigma_gr1_gr3 + gr1^2*gr4^2*sigma_gr1_gr4 + gr2^2*gr3^2*sigma_gr2_gr3 + gr2^2*gr4^2*sigma_gr2_gr4 + gr3^2*gr4^2*sigma_gr3_gr4)/(gr0^2*gr1^2 + gr0^2*gr2^2 + gr0^2*gr3^2 + gr0^2*gr4^2 + gr1^2*gr2^2 + gr1^2*gr3^2 + gr1^2*gr4^2 + gr2^2*gr3^2 + gr2^2*gr4^2 + gr3^2*gr4^2 + 1e-4)'
inputs = 'sigma_gr0_gr1 sigma_gr0_gr2 sigma_gr0_gr3 sigma_gr0_gr4 sigma_gr1_gr2 sigma_gr1_gr3 sigma_gr1_gr4 sigma_gr2_gr3 sigma_gr2_gr4 sigma_gr3_gr4 gr0 gr1 gr2 gr3 gr4'
[]
[L_NL]
type = ParsedCompute
buffer = 'L_NL'
expression = '${L} * sigma / ${gbe_max}'
inputs = 'sigma'
[]
## kappa * laplacian(gr_i) in real space (for the mobility-split correction)
[kappa_laplacian_gr0_hat]
type = ParsedCompute
buffer = kappa_laplacian_gr0_hat
inputs = 'gr0_bar kappa_laplacian'
expression = 'gr0_bar * kappa_laplacian'
[]
[kappa_laplacian_gr1_hat]
type = ParsedCompute
buffer = kappa_laplacian_gr1_hat
inputs = 'gr1_bar kappa_laplacian'
expression = 'gr1_bar * kappa_laplacian'
[]
[kappa_laplacian_gr2_hat]
type = ParsedCompute
buffer = kappa_laplacian_gr2_hat
inputs = 'gr2_bar kappa_laplacian'
expression = 'gr2_bar * kappa_laplacian'
[]
[kappa_laplacian_gr3_hat]
type = ParsedCompute
buffer = kappa_laplacian_gr3_hat
inputs = 'gr3_bar kappa_laplacian'
expression = 'gr3_bar * kappa_laplacian'
[]
[kappa_laplacian_gr4_hat]
type = ParsedCompute
buffer = kappa_laplacian_gr4_hat
inputs = 'gr4_bar kappa_laplacian'
expression = 'gr4_bar * kappa_laplacian'
[]
[kappa_laplacian_gr0]
type = InverseFFT
buffer = kappa_laplacian_gr0
input = kappa_laplacian_gr0_hat
[]
[kappa_laplacian_gr1]
type = InverseFFT
buffer = kappa_laplacian_gr1
input = kappa_laplacian_gr1_hat
[]
[kappa_laplacian_gr2]
type = InverseFFT
buffer = kappa_laplacian_gr2
input = kappa_laplacian_gr2_hat
[]
[kappa_laplacian_gr3]
type = InverseFFT
buffer = kappa_laplacian_gr3
input = kappa_laplacian_gr3_hat
[]
[kappa_laplacian_gr4]
type = InverseFFT
buffer = kappa_laplacian_gr4
input = kappa_laplacian_gr4_hat
[]
## bulk driving forces
[gr0_bulk_term]
type = ParsedCompute
buffer = 'gr0_bulk_term'
expression = '-L_NL*${mu} * (gr0^3 - gr0 + 2*gr0*(gamma_01*gr1^2 + gamma_02*gr2^2 + gamma_03*gr3^2 + gamma_04*gr4^2) - torque_01_gr0 - torque_02_gr0 - torque_03_gr0 - torque_04_gr0) + (L_NL - ${L}) * kappa_laplacian_gr0'
inputs = 'gr0 gr1 gr2 gr3 gr4 gamma_01 gamma_02 gamma_03 gamma_04 torque_01_gr0 torque_02_gr0 torque_03_gr0 torque_04_gr0 L_NL kappa_laplacian_gr0'
[]
[gr1_bulk_term]
type = ParsedCompute
buffer = 'gr1_bulk_term'
expression = '-L_NL*${mu} * (gr1^3 - gr1 + 2*gr1*(gamma_01*gr0^2 + gamma_12*gr2^2 + gamma_13*gr3^2 + gamma_14*gr4^2) - torque_01_gr1 - torque_12_gr1 - torque_13_gr1 - torque_14_gr1) + (L_NL - ${L}) * kappa_laplacian_gr1'
inputs = 'gr0 gr1 gr2 gr3 gr4 gamma_01 gamma_12 gamma_13 gamma_14 torque_01_gr1 torque_12_gr1 torque_13_gr1 torque_14_gr1 L_NL kappa_laplacian_gr1'
[]
[gr2_bulk_term]
type = ParsedCompute
buffer = 'gr2_bulk_term'
expression = '-L_NL*${mu} * (gr2^3 - gr2 + 2*gr2*(gamma_02*gr0^2 + gamma_12*gr1^2 + gamma_23*gr3^2 + gamma_24*gr4^2) - torque_02_gr2 - torque_12_gr2 - torque_23_gr2 - torque_24_gr2) + (L_NL - ${L}) * kappa_laplacian_gr2'
inputs = 'gr0 gr1 gr2 gr3 gr4 gamma_02 gamma_12 gamma_23 gamma_24 torque_02_gr2 torque_12_gr2 torque_23_gr2 torque_24_gr2 L_NL kappa_laplacian_gr2'
[]
[gr3_bulk_term]
type = ParsedCompute
buffer = 'gr3_bulk_term'
expression = '-L_NL*${mu} * (gr3^3 - gr3 + 2*gr3*(gamma_03*gr0^2 + gamma_13*gr1^2 + gamma_23*gr2^2 + gamma_34*gr4^2) - torque_03_gr3 - torque_13_gr3 - torque_23_gr3 - torque_34_gr3) + (L_NL - ${L}) * kappa_laplacian_gr3'
inputs = 'gr0 gr1 gr2 gr3 gr4 gamma_03 gamma_13 gamma_23 gamma_34 torque_03_gr3 torque_13_gr3 torque_23_gr3 torque_34_gr3 L_NL kappa_laplacian_gr3'
[]
[gr4_bulk_term]
type = ParsedCompute
buffer = 'gr4_bulk_term'
expression = '-L_NL*${mu} * (gr4^3 - gr4 + 2*gr4*(gamma_04*gr0^2 + gamma_14*gr1^2 + gamma_24*gr2^2 + gamma_34*gr3^2) - torque_04_gr4 - torque_14_gr4 - torque_24_gr4 - torque_34_gr4) + (L_NL - ${L}) * kappa_laplacian_gr4'
inputs = 'gr0 gr1 gr2 gr3 gr4 gamma_04 gamma_14 gamma_24 gamma_34 torque_04_gr4 torque_14_gr4 torque_24_gr4 torque_34_gr4 L_NL kappa_laplacian_gr4'
[]
[NL_gr0]
type = ForwardFFT
buffer = NL_gr0
input = 'gr0_bulk_term'
[]
[NL_gr1]
type = ForwardFFT
buffer = NL_gr1
input = 'gr1_bulk_term'
[]
[NL_gr2]
type = ForwardFFT
buffer = NL_gr2
input = 'gr2_bulk_term'
[]
[NL_gr3]
type = ForwardFFT
buffer = NL_gr3
input = 'gr3_bulk_term'
[]
[NL_gr4]
type = ForwardFFT
buffer = NL_gr4
input = 'gr4_bulk_term'
[]
[NL_gr0_smooth]
type = ParsedCompute
buffer = NL_gr0_smooth
expression = 'smooth*NL_gr0'
inputs = 'smooth NL_gr0'
[]
[NL_gr1_smooth]
type = ParsedCompute
buffer = NL_gr1_smooth
expression = 'smooth*NL_gr1'
inputs = 'smooth NL_gr1'
[]
[NL_gr2_smooth]
type = ParsedCompute
buffer = NL_gr2_smooth
expression = 'smooth*NL_gr2'
inputs = 'smooth NL_gr2'
[]
[NL_gr3_smooth]
type = ParsedCompute
buffer = NL_gr3_smooth
expression = 'smooth*NL_gr3'
inputs = 'smooth NL_gr3'
[]
[NL_gr4_smooth]
type = ParsedCompute
buffer = NL_gr4_smooth
expression = 'smooth*NL_gr4'
inputs = 'smooth NL_gr4'
[]
# Properties for output
[grain_ID]
type = ParsedCompute
buffer = 'grain_ID'
expression = 'if(gr0>=max(gr1,max(gr2,max(gr3,gr4))), 0, if(gr1>=max(gr2,max(gr3,gr4)), 1, if(gr2>=max(gr3,gr4), 2, if(gr3>=gr4, 3, 4))))'
inputs = 'gr0 gr1 gr2 gr3 gr4'
[]
[gradient_energy]
type = GradientEnergyDensity
buffer = gradient_energy
gradient_buffers = 'grad_gr0 grad_gr1 grad_gr2 grad_gr3 grad_gr4'
kappa = ${kappa}
[]
[total_energy]
type = ParsedCompute
buffer = 'total_energy'
expression = '${mu} * (${f0}) + gradient_energy'
inputs = 'gr0 gr1 gr2 gr3 gr4 gamma_01 gamma_02 gamma_03 gamma_04 gamma_12 gamma_13 gamma_14 gamma_23 gamma_24 gamma_34 gradient_energy'
[]
[]
[]
[Postprocessors]
[total_gb_energy]
type = TensorIntegralPostprocessor
buffer = total_energy
[]
[]
[TensorSolver]
type = AdamsBashforthMoulton
buffer = 'gr0 gr1 gr2 gr3 gr4'
reciprocal_buffer = 'gr0_bar gr1_bar gr2_bar gr3_bar gr4_bar'
linear_reciprocal = 'L_kappa_laplacian L_kappa_laplacian L_kappa_laplacian L_kappa_laplacian L_kappa_laplacian'
nonlinear_reciprocal = 'NL_gr0_smooth NL_gr1_smooth NL_gr2_smooth NL_gr3_smooth NL_gr4_smooth'
substeps = 1000
predictor_order = 1
corrector_order = 1
corrector_steps = 1
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'gr0 gr1 gr2 gr3 gr4 bnds sigma grain_ID'
output_mode = 'NODE NODE NODE NODE NODE NODE NODE CELL'
enable_hdf5 = true
transpose = true
[]
[]
[Executioner]
type = Transient
dt = 5
num_steps = 50
[]
[Outputs]
csv = true
execute_on = 'INITIAL TIMESTEP_END'
[]
(test/tests/kks/KKS_no_flux_bc.i)
#
# Kim-Kim-Suzuki with no-flux BC imposed using the smooth boundary method (SBM), solved on a 2D grid.
# Mask tensor 'psi' supplies the mask for the solve region to the system.
# Note: c is not directly conserved here - the masked value (psi > 0.0)*c will however be conserved.
#
# Constants for Initial Conditions
r = 30
l = 4.2
# Initial condition function for order parameter
eta_IC = '0.5*(1-tanh(2*(sqrt(x^2+y^2)-${r})/${l}))'
# Phase-field model parameters
kappa_eta = 5
rho_sq = 2
w = 1
M = 5
L = 5
c0_a = 0.3
c0_b = 0.7
# Expressions for switching function and bulk Gibbs energy
h_eta = 'eta^3*(6*eta^2-15*eta+10)'
F = '${h_eta}*(${rho_sq}*((c - (1-${h_eta})*(${c0_b} - ${c0_a}))-${c0_a})^2) + (1-${h_eta})*(${rho_sq}*((c + (${h_eta})*(${c0_b} - ${c0_a}))-${c0_b})^2 ) + ${w}*(eta^2)*(1-eta)^2'
[Domain]
dim = 2
nx = 20
ny = 20
xmin = -50
xmax = 50
ymin = -50
ymax = 50
# run on a CUDA device (adjust this to `cpu` if not available)
device_names = 'cpu'
# automatically create a matching mesh
mesh_mode = DUMMY
[]
[Functions]
[psi_func]
type = ParsedFunction
expression = 'if(x<x_min-${l},0,if(x>x_min+${l},1,0.5-0.5*cos(pi*(x-(x_min-${l}))/2/${l}) )) * if(x<x_max-${l},1,if(x>x_max+${l},0,0.5+0.5*cos(pi*(x-(x_max-${l}))/2/${l}) ))'
symbol_names = 'x_min x_max y_min y_max'
symbol_values = '-20 20 -50 50'
[]
[]
[TensorComputes]
[Initialize]
[c_IC]
type = ParsedCompute
buffer = c
expression = '0.6 + (${c0_a}-0.6)*${eta_IC}'
extra_symbols = 'true'
[]
[eta_IC]
type = ParsedCompute
buffer = eta
expression = '${eta_IC}'
extra_symbols = 'true'
[]
[psi_init]
type = MooseFunctionTensor
function = psi_func
buffer = psi
[]
[zero]
type = ConstantReciprocalTensor
buffer = zero
[]
[M]
type = ConstantTensor
buffer = M
real = ${M}
[]
[L]
type = ConstantTensor
buffer = L
real = ${L}
[]
[L_kappa]
type = ConstantTensor
buffer = L_kappa
real = ${fparse ${L} * ${kappa_eta} }
[]
[]
[Solve]
[cbar]
type = ForwardFFT
buffer = cbar
input = c
[]
[etabar]
type = ForwardFFT
buffer = etabar
input = eta
[]
[mu]
type = ParsedCompute
buffer = 'mu'
expression = '${F}'
inputs = 'c eta'
derivatives = 'c'
[]
[div_J]
type = ReciprocalMatDiffusion
buffer = 'div_J'
chemical_potential = mu
mobility = M
psi = psi
[]
[domega_chem_deta]
type = ParsedCompute
buffer = 'domega_chem_deta'
expression = '${F} - mu*c'
inputs = 'mu c eta'
derivatives = 'eta'
[]
[AC_bulk]
type = ReciprocalAllenCahn
buffer = AC_bulk
dF_chem_deta = domega_chem_deta
L = L
psi = psi
[]
[kappa_grad_eta]
type = ReciprocalMatDiffusion
buffer = 'kappa_grad_eta'
chemical_potential = 'eta'
mobility = 'L_kappa'
psi = psi
[]
[AC_bar]
type = ParsedCompute
buffer = AC_bar
expression = 'kappa_grad_eta + AC_bulk'
inputs = 'AC_bulk kappa_grad_eta'
[]
[]
[]
[TensorSolver]
type = AdamsBashforthMoulton
buffer = 'c eta'
reciprocal_buffer = 'cbar etabar'
linear_reciprocal = 'zero zero'
nonlinear_reciprocal = 'div_J AC_bar'
substeps = 1e3
predictor_order = 3
[]
[Postprocessors]
[total_C]
type = TensorIntegralPostprocessor
buffer = c
execute_on = 'INITIAL TIMESTEP_END'
[]
[total_eta]
type = TensorIntegralPostprocessor
buffer = eta
execute_on = 'INITIAL TIMESTEP_END'
[]
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'eta c mu psi'
enable_hdf5 = true
transpose = false
[]
[]
[Executioner]
type = Transient
dt = 0.1
num_steps = 10
[]
[Outputs]
csv = true
perf_graph = true
execute_on = 'INITIAL TIMESTEP_END'
[]
(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
[]
[]
[TensorSolver]
type = LBMStream
buffer = 'h f'
f_old = 'h_post_collision f_post_collision'
[]
[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}'
residual_tensor = speed
[]
[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
[]
[]
(examples/anisotropic_grain_growth/spherical_grain_Ni.i)
interface_width = 1.6
r0 = 7.5
# Ni
gbe_max = '${units 0.8797226 J/m^2}'
L = '${fparse 1.0 * 1.6 / ${interface_width}}'
g_gamma0 = '${fparse sqrt(2) / 3}' # g(gamma=1.5)
f0_gamma0 = ${fparse 1/8}
kappa = '${fparse gbe_max * interface_width * sqrt(f0_gamma0) / g_gamma0}'
mu = '${fparse gbe_max / (g_gamma0 * interface_width * sqrt(f0_gamma0))}'
g_to_gamma_func = '(1.0/(-3.0944 * g^8 -1.8169*g^6 + 10.323 * g^4 - 8.1819*g^2 + 2.0033))'
dgamma_dsigma = '(g * (24.7552*g^6 + 10.9014*g^4 - 41.292*g^2 + 16.3638) / (-3.0944*g^8 - 1.8169*g^6 + 10.323*g^4 - 8.1819*g^2 + 2.0033)^2 / sqrt(${kappa} * ${mu}))'
f0 = '((gr0^4/4 - gr0^2/2) + (gr1^4/4 - gr1^2/2) + gamma_01*gr0^2*gr1^2 + 0.25)'
[Domain]
dim = 3
nx = 100
ny = 100
nz = 100
xmax = 20
ymax = 20
zmax = 20
mesh_mode = DUMMY
[]
[TensorComputes]
[Initialize]
[gr0]
type = ParsedCompute
buffer = 'gr0'
extra_symbols = 'true'
expression = 'radius:=sqrt((x-${Domain/xmax}/2)^2+(y-${Domain/ymax}/2)^2+(z-${Domain/zmax}/2)^2);0.5 - 0.5*tanh(2*(radius-${r0})/${interface_width})'
[]
[gr1]
type = ParsedCompute
buffer = 'gr1'
extra_symbols = 'true'
expression = 'radius:=sqrt((x-${Domain/xmax}/2)^2+(y-${Domain/ymax}/2)^2+(z-${Domain/zmax}/2)^2);0.5 + 0.5*tanh(2*(radius-${r0})/${interface_width})'
[]
[L_kappa_laplacian]
type = ReciprocalLaplacianFactor
buffer = L_kappa_laplacian
factor = ${fparse ${L} * ${kappa}}
[]
[kappa_laplacian]
type = ReciprocalLaplacianFactor
buffer = kappa_laplacian
factor = '${kappa}'
[]
[smooth]
type = DeAliasingTensor
method = HOULI
buffer = smooth
[]
[]
[Solve]
## FFT
[gr0_bar]
type = ForwardFFT
buffer = gr0_bar
input = gr0
[]
[gr1_bar]
type = ForwardFFT
buffer = gr1_bar
input = gr1
[]
[gr0_bar_smooth]
type = ParsedCompute
buffer = gr0_bar_smooth
expression = 'gr0_bar * smooth'
inputs = 'gr0_bar smooth'
[]
[gr1_bar_smooth]
type = ParsedCompute
buffer = gr1_bar_smooth
expression = 'gr1_bar * smooth'
inputs = 'gr1_bar smooth'
[]
# Gradients
[grad_gr0]
type = GradientVector
buffer = grad_gr0
input = gr0_bar_smooth
input_is_reciprocal = true
[]
[grad_gr1]
type = GradientVector
buffer = grad_gr1
input = gr1_bar_smooth
input_is_reciprocal = true
[]
# Get interface energy
[sigma_gr0_gr1]
type = PairwiseAnisotropicGBEnergy
buffer = 'sigma_gr0_gr1'
dsigma_dgrad_grain1 = 'dsigma_gr0_gr1_dgrad_gr0'
dsigma_dgrad_grain2 = 'dsigma_gr0_gr1_dgrad_gr1'
grad_grain1_buffer = 'grad_gr0'
grad_grain2_buffer = 'grad_gr1'
interface_width = '${interface_width}'
libtorch_model_file = 'marlin:anisotropic_gb_energy/Ni_R1_R3.pt'
chunk_size = 1e5
[]
[gamma_01]
type = ParsedCompute
buffer = 'gamma_01'
expression = 'g:=sigma_gr0_gr1/sqrt(${kappa}*${mu});${g_to_gamma_func}'
inputs = 'sigma_gr0_gr1'
[]
[coeff_01]
type = ParsedCompute
buffer = 'coeff_01'
expression = 'g:=sigma_gr0_gr1/sqrt(${kappa}*${mu});${dgamma_dsigma} * gr0^2 *gr1^2'
inputs = 'sigma_gr0_gr1 gr0 gr1'
[]
[torque_01_gr0]
type = AnisotropyTorque
buffer = torque_01_gr0
dmu_dn = dsigma_gr0_gr1_dgrad_gr0
g = 'coeff_01'
[]
[torque_01_gr1]
type = AnisotropyTorque
buffer = torque_01_gr1
dmu_dn = dsigma_gr0_gr1_dgrad_gr1
g = 'coeff_01'
[]
[kappa_laplacian_gr0_hat]
type = ParsedCompute
buffer = kappa_laplacian_gr0_hat
inputs = 'gr0_bar kappa_laplacian'
expression = 'gr0_bar * kappa_laplacian'
[]
[kappa_laplacian_gr1_hat]
type = ParsedCompute
buffer = kappa_laplacian_gr1_hat
inputs = 'gr1_bar kappa_laplacian'
expression = 'gr1_bar * kappa_laplacian'
[]
[kappa_laplacian_gr0]
type = InverseFFT
buffer = kappa_laplacian_gr0
input = kappa_laplacian_gr0_hat
[]
[kappa_laplacian_gr1]
type = InverseFFT
buffer = kappa_laplacian_gr1
input = kappa_laplacian_gr1_hat
[]
[gr0_bulk_term]
type = ParsedCompute
buffer = 'gr0_bulk_term'
expression = 'L_NL:=${L}*(sigma_gr0_gr1)/${gbe_max};-L_NL*${mu} * (gr0^3 - gr0 + 2*gr0*(gamma_01*gr1^2) - torque_01_gr0) + (L_NL - ${L})*kappa_laplacian_gr0'
inputs = 'gr0 gr1 gamma_01 torque_01_gr0 sigma_gr0_gr1 kappa_laplacian_gr0'
[]
[gr1_bulk_term]
type = ParsedCompute
buffer = 'gr1_bulk_term'
expression = 'L_NL:=${L}*(sigma_gr0_gr1)/${gbe_max};-L_NL*${mu} * (gr1^3 - gr1 + 2*gr1*(gamma_01*gr0^2) - torque_01_gr1) + (L_NL - ${L})*kappa_laplacian_gr1'
inputs = 'gr0 gr1 gamma_01 torque_01_gr1 sigma_gr0_gr1 kappa_laplacian_gr1'
[]
## build non-linear terms and smooth them
[NL_gr0]
type = ForwardFFT
buffer = NL_gr0
input = 'gr0_bulk_term'
[]
[NL_gr1]
type = ForwardFFT
buffer = NL_gr1
input = 'gr1_bulk_term'
[]
[NL_gr0_smooth]
type = ParsedCompute
buffer = NL_gr0_smooth
expression = 'smooth*NL_gr0'
inputs = 'smooth NL_gr0'
[]
[NL_gr1_smooth]
type = ParsedCompute
buffer = NL_gr1_smooth
expression = 'smooth*NL_gr1'
inputs = 'smooth NL_gr1'
[]
# Properties for output
[sigma]
type = ParsedCompute
buffer = 'sigma'
expression = 'gr0^2*gr1^2*sigma_gr0_gr1'
inputs = 'sigma_gr0_gr1 gr0 gr1'
[]
[gradient_energy]
type = GradientEnergyDensity
buffer = gradient_energy
gradient_buffers = 'grad_gr0 grad_gr1'
kappa = ${kappa}
[]
[total_energy]
type = ParsedCompute
buffer = 'total_energy'
expression = '${mu} * (${f0}) + gradient_energy'
inputs = 'gr0 gr1 gamma_01 gradient_energy'
[]
[]
[]
[Postprocessors]
[total_gb_energy]
type = TensorIntegralPostprocessor
buffer = total_energy
[]
[]
[TensorSolver]
type = AdamsBashforthMoulton
buffer = 'gr0 gr1'
reciprocal_buffer = 'gr0_bar gr1_bar'
linear_reciprocal = 'L_kappa_laplacian L_kappa_laplacian'
nonlinear_reciprocal = 'NL_gr0_smooth NL_gr1_smooth'
substeps = 100
predictor_order = 2
corrector_order = 1
corrector_steps = 2
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'gr0 gr1 sigma_gr0_gr1'
output_mode = 'NODE NODE NODE'
enable_hdf5 = true
transpose = false
[]
[]
[Executioner]
type = Transient
dt = 0.1
num_steps = 200
[]
[Outputs]
csv = true
execute_on = 'INITIAL TIMESTEP_END'
[]
(test/tests/anisotropic_grain_growth/sigma3_circular_grain.i)
interface_width = 2.0
r0 = 3
gbe_max = '${units 0.8797226 J/m^2}'
L = '${fparse 1.0 * 1.6 / ${interface_width}}'
g_gamma0 = '${fparse sqrt(2) / 3}' # g(gamma=1.5)
f0_gamma0 = ${fparse 1/8}
kappa = '${fparse gbe_max * interface_width * sqrt(f0_gamma0) / g_gamma0}'
mu = '${fparse gbe_max / (g_gamma0 * interface_width * sqrt(f0_gamma0))}'
g_to_gamma_func = '(1.0/(-3.0944 * g^8 -1.8169*g^6 + 10.323 * g^4 - 8.1819*g^2 + 2.0033))'
dgamma_dsigma = '(g * (24.7552*g^6 + 10.9014*g^4 - 41.292*g^2 + 16.3638) / (-3.0944*g^8 - 1.8169*g^6 + 10.323*g^4 - 8.1819*g^2 + 2.0033)^2 / sqrt(${kappa} * ${mu}))'
f0 = '((gr0^4/4 - gr0^2/2) + (gr1^4/4 - gr1^2/2) + gamma_01*gr0^2*gr1^2 + 0.25)'
[Domain]
dim = 2
nx = 32
ny = 32
xmax = 8
ymax = 8
mesh_mode = DUMMY
floating_precision = SINGLE
[]
[TensorComputes]
[Initialize]
[gr0]
type = ParsedCompute
buffer = 'gr0'
extra_symbols = 'true'
expression = 'radius:=sqrt((x-${Domain/xmax}/2)^2+(y-${Domain/ymax}/2)^2);0.5 - 0.5*tanh(2*(radius-${r0})/${interface_width})'
[]
[gr1]
type = ParsedCompute
buffer = 'gr1'
extra_symbols = 'true'
expression = 'radius:=sqrt((x-${Domain/xmax}/2)^2+(y-${Domain/ymax}/2)^2);0.5 + 0.5*tanh(2*(radius-${r0})/${interface_width})'
[]
[L_kappa_laplacian]
type = ReciprocalLaplacianFactor
buffer = L_kappa_laplacian
factor = ${fparse ${L} * ${kappa}}
[]
[kappa_laplacian]
type = ReciprocalLaplacianFactor
buffer = kappa_laplacian
factor = '${kappa}'
[]
[smooth]
type = DeAliasingTensor
method = HOULI
buffer = smooth
[]
[]
[Solve]
## FFT
[gr0_bar]
type = ForwardFFT
buffer = gr0_bar
input = gr0
[]
[gr1_bar]
type = ForwardFFT
buffer = gr1_bar
input = gr1
[]
# Gradients
[grad_gr0]
type = GradientVector
buffer = grad_gr0
input = gr0_bar
input_is_reciprocal = true
[]
[grad_gr1]
type = GradientVector
buffer = grad_gr1
input = gr1_bar
input_is_reciprocal = true
[]
# Get interface energy
[sigma_gr0_gr1]
type = PairwiseAnisotropicGBEnergy
buffer = 'sigma_gr0_gr1'
dsigma_dgrad_grain1 = 'dsigma_gr0_gr1_dgrad_gr0'
dsigma_dgrad_grain2 = 'dsigma_gr0_gr1_dgrad_gr1'
grad_grain1_buffer = 'grad_gr0'
grad_grain2_buffer = 'grad_gr1'
interface_width = '${interface_width}'
libtorch_model_file = 'marlin:anisotropic_gb_energy/gb_energy_hull_2d.pt'
chunk_size = 256
[]
[gamma_01]
type = ParsedCompute
buffer = 'gamma_01'
expression = 'g:=sigma_gr0_gr1/sqrt(${kappa}*${mu});${g_to_gamma_func}'
inputs = 'sigma_gr0_gr1'
[]
[coeff_01]
type = ParsedCompute
buffer = 'coeff_01'
expression = 'g:=sigma_gr0_gr1/sqrt(${kappa}*${mu});${dgamma_dsigma} * gr0^2 *gr1^2'
inputs = 'sigma_gr0_gr1 gr0 gr1'
[]
[torque_01_gr0]
type = AnisotropyTorque
buffer = torque_01_gr0
dmu_dn = dsigma_gr0_gr1_dgrad_gr0
g = 'coeff_01'
[]
[torque_01_gr1]
type = AnisotropyTorque
buffer = torque_01_gr1
dmu_dn = dsigma_gr0_gr1_dgrad_gr1
g = 'coeff_01'
[]
[kappa_laplacian_gr0_hat]
type = ParsedCompute
buffer = kappa_laplacian_gr0_hat
inputs = 'gr0_bar kappa_laplacian'
expression = 'gr0_bar * kappa_laplacian'
[]
[kappa_laplacian_gr1_hat]
type = ParsedCompute
buffer = kappa_laplacian_gr1_hat
inputs = 'gr1_bar kappa_laplacian'
expression = 'gr1_bar * kappa_laplacian'
[]
[kappa_laplacian_gr0]
type = InverseFFT
buffer = kappa_laplacian_gr0
input = kappa_laplacian_gr0_hat
[]
[kappa_laplacian_gr1]
type = InverseFFT
buffer = kappa_laplacian_gr1
input = kappa_laplacian_gr1_hat
[]
[gr0_bulk_term]
type = ParsedCompute
buffer = 'gr0_bulk_term'
expression = 'L_NL:=${L}*(sigma_gr0_gr1)/${gbe_max};-L_NL*${mu} * (gr0^3 - gr0 + 2*gr0*(gamma_01*gr1^2) - torque_01_gr0) + (L_NL - ${L})*kappa_laplacian_gr0'
inputs = 'gr0 gr1 gamma_01 torque_01_gr0 sigma_gr0_gr1 kappa_laplacian_gr0'
[]
[gr1_bulk_term]
type = ParsedCompute
buffer = 'gr1_bulk_term'
expression = 'L_NL:=${L}*(sigma_gr0_gr1)/${gbe_max};-L_NL*${mu} * (gr1^3 - gr1 + 2*gr1*(gamma_01*gr0^2) - torque_01_gr1) + (L_NL - ${L})*kappa_laplacian_gr1'
inputs = 'gr0 gr1 gamma_01 torque_01_gr1 sigma_gr0_gr1 kappa_laplacian_gr1'
[]
## build non-linear terms and smooth them
[NL_gr0]
type = ForwardFFT
buffer = NL_gr0
input = 'gr0_bulk_term'
[]
[NL_gr1]
type = ForwardFFT
buffer = NL_gr1
input = 'gr1_bulk_term'
[]
[NL_gr0_smooth]
type = ParsedCompute
buffer = NL_gr0_smooth
expression = 'smooth*NL_gr0'
inputs = 'smooth NL_gr0'
[]
[NL_gr1_smooth]
type = ParsedCompute
buffer = NL_gr1_smooth
expression = 'smooth*NL_gr1'
inputs = 'smooth NL_gr1'
[]
# Properties for output
[sigma]
type = ParsedCompute
buffer = 'sigma'
expression = 'gr0^2*gr1^2*sigma_gr0_gr1'
inputs = 'sigma_gr0_gr1 gr0 gr1'
[]
[gradient_energy]
type = GradientEnergyDensity
buffer = gradient_energy
gradient_buffers = 'grad_gr0 grad_gr1'
kappa = ${kappa}
[]
[total_energy]
type = ParsedCompute
buffer = 'total_energy'
expression = '${mu} * (${f0}) + gradient_energy'
inputs = 'gr0 gr1 gamma_01 gradient_energy'
[]
[]
[]
[Postprocessors]
[total_gb_energy]
type = TensorIntegralPostprocessor
buffer = total_energy
[]
[]
[TensorSolver]
type = AdamsBashforthMoulton
buffer = 'gr0 gr1'
reciprocal_buffer = 'gr0_bar gr1_bar'
linear_reciprocal = 'L_kappa_laplacian L_kappa_laplacian'
nonlinear_reciprocal = 'NL_gr0_smooth NL_gr1_smooth'
substeps = 1e2
predictor_order = 1
corrector_order = 1
corrector_steps = 1
[]
[TensorOutputs]
active = ''
[xdmf]
type = XDMFTensorOutput
buffer = 'gr0 gr1 sigma_gr0_gr1 sigma total_energy'
output_mode = 'NODE NODE NODE NODE NODE'
enable_hdf5 = true
transpose = false
[]
[]
[Executioner]
type = Transient
dt = 1
num_steps = 2
[]
[Outputs]
csv = true
execute_on = 'INITIAL TIMESTEP_END'
[]
(test/tests/tensor_compute/parallel_roundtrip_3d.i)
[Domain]
# Test parallel FFT round-trip with slab decomposition in 3D
device_names = "cpu cpu cpu"
device_weights = "1 1 1"
dim = 3
nx = 64
ny = 64
nz = 64
xmax = ${fparse pi*4}
ymax = ${fparse pi*4}
zmax = ${fparse pi*4}
[]
[TensorBuffers]
[eta_gold]
[]
[eta]
[]
[eta_bar]
[]
[eta_roundtrip]
[]
[diff]
[]
[zero]
[]
[]
[TensorComputes]
[Initialize]
[eta_gold]
type = ParsedCompute
buffer = eta_gold
expression = 'sin(x)+sin(y)+sin(z)+cos(2*x)*sin(3*y)*cos(z)'
extra_symbols = true
[]
[eta]
type = ParsedCompute
buffer = eta
expression = eta_gold
inputs = eta_gold
[]
[zero]
type = ConstantReciprocalTensor
buffer = zero
real = 0
imaginary = 0
[]
[]
[Solve]
# Test: eta -> FFT -> iFFT -> eta_roundtrip
# eta_roundtrip should equal eta (within numerical precision)
[eta_bar]
type = ForwardFFT
buffer = eta_bar
input = eta
[]
[eta_roundtrip]
type = InverseFFT
buffer = eta_roundtrip
input = eta_bar
[]
[]
[Postprocess]
[diff]
type = ParsedCompute
buffer = diff
expression = 'abs(eta - eta_roundtrip) + abs(eta - eta_gold)'
inputs = 'eta eta_roundtrip eta_gold'
[]
[]
[]
[Postprocessors]
[max_error]
type = TensorExtremeValuePostprocessor
buffer = diff
value_type = MAX
[]
[l2_error]
type = TensorIntegralPostprocessor
buffer = diff
[]
[]
[TensorSolver]
type = AdamsBashforthMoulton
buffer = eta
reciprocal_buffer = eta_bar
linear_reciprocal = zero
nonlinear_reciprocal = zero
[]
[TensorOutputs]
[eta]
type = XDMFTensorOutput
buffer = 'eta'
output_mode = 'CELL'
enable_hdf5 = true
[]
[]
[Problem]
type = TensorProblem
[]
[Executioner]
type = Transient
num_steps = 1
[]
[Outputs]
csv = true
execute_on = 'INITIAL TIMESTEP_END'
[]
(examples/phase_field_crystal/pfc_fcc_atomic_structure_3d.i)
# 3D FCC PFC test showing atomic-scale density fluctuations
# Based on PhysRevE.81.061601 two-mode FCC PFC model
# Grid parameters - 3D resolution
N = 64 # Grid points per dimension (64^3 = 262k points)
# PFC model parameters - LARGER epsilon for clearer atomic structure
psi_mean = -0.2 # Mean density
epsilon = 0.15 # Larger undercooling for stronger modulations
R1_param = 0.0 # Two-mode coupling (0 for max FCC stability)
Q1_param = ${fparse sqrt(4.0/3.0)} # FCC wave number ratio
# Domain size: want several unit cells
# The (111) wavelength is 2*pi/q0 = 2*pi*sqrt(3/4) ≈ 5.44 (dimensionless)
# Fit 8 wavelengths in each direction
Lx = ${fparse 8 * 2 * ${pi} / ${Q1_param}}
Ly = ${Lx}
Lz = ${Lx}
[Domain]
dim = 3
nx = ${N}
ny = ${N}
nz = ${N}
xmax = ${Lx}
ymax = ${Ly}
zmax = ${Lz}
[]
[TensorComputes]
[Initialize]
# Random perturbations to seed crystal growth
[psi]
type = RandomTensor
buffer = psi
max = ${fparse ${psi_mean} + 0.01}
min = ${fparse ${psi_mean} - 0.01}
seed = 12345
[]
# Linear operator for FCC
[linear]
type = FCCPFCLinear
buffer = 'linear'
eps = ${epsilon}
R1 = ${R1_param}
Q1 = ${Q1_param}
mobility = 1.0
[]
# Dealiasing for cubic nonlinearity
[smooth_operator]
type = DeAliasingTensor
buffer = smooth_operator
method = HOULI
[]
[]
[Solve]
# Nonlinear term
[nl_div_psi_cubed]
type = FCCPFCNonlinear
buffer = NL
psi = psi
dealiasing = smooth_operator
mobility = 1.0
[]
# FFT for spectral solver
[psi_hat]
type = ForwardFFT
buffer = psi_hat
input = psi
[]
[]
[]
[TensorSolver]
type = AdamsBashforthMoulton
buffer = 'psi'
linear_reciprocal = 'linear'
nonlinear_reciprocal = 'NL'
reciprocal_buffer = 'psi_hat'
corrector_order = 1
corrector_steps = 3
predictor_order = 1
substeps = 10000
[]
[Postprocessors]
[max_psi]
type = TensorExtremeValuePostprocessor
buffer = psi
value_type = MAX
execute_on = 'initial timestep_end'
[]
[min_psi]
type = TensorExtremeValuePostprocessor
buffer = psi
value_type = MIN
execute_on = 'initial timestep_end'
[]
[mean_psi]
type = TensorAveragePostprocessor
buffer = psi
execute_on = 'initial timestep_end'
[]
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'psi'
enable_hdf5 = true
[]
[]
[Executioner]
type = Transient
num_steps = 200 # Reduced for 3D (much more expensive)
dt = 5
[]
[Outputs]
perf_graph = true
csv = true
[]
(examples/degeus_mechanics/mech.i)
[Domain]
dim = 3
nx = 32
ny = 32
nz = 32
xmax = ${fparse 2*pi}
ymax = ${fparse 2*pi}
zmax = ${fparse 2*pi}
mesh_mode = DUMMY
[]
[TensorComputes]
[Initialize]
[Finit]
type = RankTwoIdentity
buffer = F
[]
[phase]
type = PhaseMechanicsTest
buffer = phase
[]
[K]
type = ParsedCompute
buffer = K
expression = '(1-phase)*Ka + phase*Kb'
inputs = phase
constant_names = 'Ka Kb'
constant_expressions = '0.833 8.33'
[]
[mu]
type = ParsedCompute
buffer = mu
expression = '(1-phase)*mua + phase*mub'
inputs = phase
constant_names = 'mua mub'
constant_expressions = '0.386 3.86'
[]
[]
[Solve]
[hyper_elasticity]
type = HyperElasticIsotropic
buffer = stress
F = Fnew
K = K
mu = mu
[]
[root]
[applied_strain]
type = MacroscopicShearTensor
buffer = applied_strain
[]
[mech]
type = FFTMechanics
buffer = Fnew
F = F
K = K
mu = mu
l_tol = 1e-2
nl_rel_tol = 2e-2
nl_abs_tol = 2e-2
constitutive_model = hyper_elasticity
stress = stress
applied_macroscopic_strain = applied_strain
[]
[]
[]
[Postprocess]
[displacements]
type = ComputeDisplacements
buffer = disp
F = F
[]
[vonmises]
type = ComputeVonMisesStress
buffer = sV
[]
[]
[]
[TensorSolver]
# no variables are integrated by this solver (FFTMechanics performs a steady state mechanics solve)
type = ForwardEulerSolver
root_compute = root
# deformation tensor is just forwarded Fnew -> F
forward_buffer = F
forward_buffer_new = Fnew
substeps = 10
[]
[TensorOutputs]
[deformation_tensor]
type = XDMFTensorOutput
buffer = 'disp sV F'
output_mode = 'OVERSIZED_NODAL CELL CELL'
enable_hdf5 = true
[]
[]
[Executioner]
type = Transient
num_steps = 100
dt = 0.01
[]
(test/tests/real_space/diffusion.i)
[Domain]
dim = 2
nx = 80
ny = 60
parallel_mode = REAL_SPACE
periodic_directions = 'X Y'
xmin = -1
ymin = -1
xmax = 1
ymax = 1
[]
[TensorComputes]
[Initialize]
[u_ic]
type = ParsedCompute
buffer = u
expression = 'exp((-x^2-y^2)*100)'
extra_symbols = true
[]
[]
[Solve]
[dt]
type = FiniteDifferenceLaplacian
input = u
buffer = dudt
[]
[]
[]
[TensorSolver]
type = RealSpaceForwardEuler
buffer = u
time_derivative = dudt
substeps = 100
[]
[Postprocessors]
[max_u]
type = TensorExtremeValuePostprocessor
buffer = u
value_type = MAX
execute_on = 'INITIAL TIMESTEP_END'
[]
[U]
type = TensorIntegralPostprocessor
buffer = u
execute_on = 'INITIAL TIMESTEP_END'
[]
[]
[TensorOutputs]
# active = ''
[u]
type = XDMFTensorOutput
buffer = 'u dudt'
enable_hdf5 = true
[]
[]
[Executioner]
type = Transient
dt = 1e-2
num_steps = 100
[]
[Outputs]
[out]
type = CSV
[]
[]
(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
[]
[]
[Boundary]
[wall]
type = LBMSpecularReflectionBoundary
buffer = f
f_old = fpc
local_Knudsen_number = Kn
boundary = wall
[]
[]
[]
[TensorSolver]
type = LBMStream
buffer = f
f_old = fpc
[]
[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
residual_tensor = speed
[]
[Executioner]
type = Transient
num_steps = 2
[]
[TensorOutputs]
[xdmf2]
type = XDMFTensorOutput
buffer = 'density velocity speed'
output_mode = 'Cell Cell Cell'
enable_hdf5 = true
transpose=false
[]
[]
(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
[]
[]
[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'
[]
[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}'
residual_tensor = speed
[]
[Executioner]
type = Transient
num_steps = 2
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'phi velocity rho'
output_mode = 'Cell Cell Cell'
enable_hdf5 = true
transpose = false
[]
[]
(test/tests/kks/KKS_libtorch.i)
#
# Kim-Kim-Suzuki with Gibbs energy supplied by a torch model, solved on a 2D grid.
#
# Constants for Initial Conditions
r = 30
l = 4.2
# Initial condition function for order parameter
eta_IC = '0.5*(1-tanh(2*(sqrt(x^2+y^2)-${r})/${l}))'
# Phase-field model parameters
kappa_eta = 5
w = 1
M = 5
L = 5
# Expressions for switching function and bulk Gibbs energy
h_eta = 'eta^3*(6*eta^2-15*eta+10)'
[Domain]
dim = 2
nx = 50
ny = 50
xmin = -50
xmax = 50
ymin = -50
ymax = 50
# automatically create a matching mesh
mesh_mode = DUMMY
[]
[TensorComputes]
[Initialize]
[c_IC]
type = ParsedCompute
buffer = c
expression = '0.7 + (0.3-0.6)*${eta_IC}'
extra_symbols = 'true'
[]
[eta_IC]
type = ParsedCompute
buffer = eta
expression = '${eta_IC}'
extra_symbols = 'true'
[]
[psi_init]
type = ConstantTensor
buffer = psi
real = 1
[]
[M]
type = ConstantTensor
buffer = M
real = ${M}
[]
[L]
type = ConstantTensor
buffer = L
real = ${L}
[]
[L_kappa]
type = ReciprocalLaplacianFactor
buffer = L_kappa
factor = ${fparse ${L} * ${kappa_eta} }
[]
[h_eta_IC]
type = ParsedCompute
buffer = h_eta
expression = '${h_eta}'
inputs = eta
[]
[G_func_IC]
type = LibtorchGibbsEnergy
buffer = 'G'
phase_fractions = 'h_eta'
concentrations = 'c'
domega_detas = 'dG_dh'
chem_pots = 'mu'
libtorch_model_file = 'marlin:libtorch_gibbs_energy/torch_NN_gibbs_model.pt'
[]
[]
[Solve]
[h_eta]
type = ParsedCompute
buffer = h_eta
expression = '${h_eta}'
inputs = eta
[]
[G_func]
type = LibtorchGibbsEnergy
buffer = 'G'
phase_fractions = 'h_eta'
concentrations = 'c'
domega_detas = 'dG_dh'
chem_pots = 'mu'
libtorch_model_file = 'marlin:libtorch_gibbs_energy/torch_NN_gibbs_model.pt'
[]
[dG_deta]
type = ParsedCompute
buffer = 'dG_deta'
inputs = 'eta dG_dh'
expression = 'dG_dh * ${h_eta} + ${w} * eta^2 * (1-eta^2)^2'
derivatives = 'eta'
[]
[etabar]
type = ForwardFFT
buffer = etabar
input = eta
[]
[AC_bulk]
type = ReciprocalAllenCahn
L = L
buffer = AC_bulk
dF_chem_deta = dG_deta
psi = psi
[]
[NL_eta]
type = ParsedCompute
buffer = NL_eta
expression = 'AC_bulk '
inputs = 'AC_bulk'
[]
[cbar]
type = ForwardFFT
buffer = cbar
input = c
[]
[div_J]
type = ReciprocalMatDiffusion
buffer = 'div_J'
chemical_potential = mu
mobility = M
psi = psi
[]
[]
[]
[TensorSolver]
type = AdamsBashforthMoulton
buffer = 'c eta'
reciprocal_buffer = 'cbar etabar'
linear_reciprocal = '0 L_kappa'
nonlinear_reciprocal = 'div_J NL_eta'
substeps = 1e3
predictor_order = 3
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'eta c mu psi dG_deta dG_dh G'
enable_hdf5 = true
transpose = false
[]
[]
[Executioner]
type = Transient
dt = 0.1
num_steps = 10
[]
[Outputs]
csv = true
perf_graph = true
execute_on = 'INITIAL TIMESTEP_END'
[]
(test/tests/tensor_compute/coupled_pf_mech.i)
[Domain]
dim = 3
nx = 128
ny = 128
nz = 128
xmax = ${fparse pi*4}
ymax = ${fparse pi*4}
zmax = ${fparse pi*4}
mesh_mode = DUMMY
[]
[TensorBuffers]
# phase field
[c]
[]
[cbar]
[]
[mu]
[]
[mubar]
[]
[Mbarmubar]
[]
# mechanics
[disp_x]
[]
[disp_y]
[]
[disp_z]
[]
[mumechbar]
[]
[mumech]
[]
# constant tensors
[Mbar]
[]
[kappabarbar]
[]
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'c disp_x disp_y disp_z mu mumech'
output_mode = 'Node Node Node Node Cell Cell'
enable_hdf5 = true
[]
[]
[TensorComputes]
[Initialize]
[c]
type = RandomTensor
buffer = c
min = 0.44
max = 0.56
[]
[disp_x]
type = RandomTensor
buffer = disp_x
min = 0
max = 0
[]
[disp_y]
type = RandomTensor
buffer = disp_y
min = 0
max = 0
[]
[disp_z]
type = RandomTensor
buffer = disp_z
min = 0
max = 0
[]
[Mbar]
type = ReciprocalLaplacianFactor
factor = 0.2 # Mobility
buffer = Mbar
[]
[kappabarbar]
type = ReciprocalLaplacianSquareFactor
factor = -0.001 # kappa
buffer = kappabarbar
[]
[]
[Solve]
[mu]
# chemical potential (real space)
type = ParsedCompute
buffer = mu
expression = '0.1*c^2*(c-1)^2' # + c*sin(x/2)*0.005'
extra_symbols = true
derivatives = c
inputs = c
[]
[mubar]
# chemical potential (reciprocal space)
type = ForwardFFT
buffer = mubar
input = mu
[]
[mumechbar]
# mechanical chemical potential (reciprocal space)
type = FFTElasticChemicalPotential
buffer = mumechbar
cbar = cbar
displacements = 'disp_x disp_y disp_z'
lambda = 100
mu = 50
e0 = 0.02
[]
[mumech]
# chemical potential (reciprocal space)
type = InverseFFT
buffer = mumech
input = mumechbar
[]
[Mbarmubar]
type = ParsedCompute
buffer = Mbarmubar
expression = 'Mbar*(mubar+mumechbar)'
inputs = 'Mbar mubar mumechbar'
[]
[cbar]
type = ForwardFFT
buffer = cbar
input = c
[]
[qsmech]
type = FFTQuasistaticElasticity
displacements = 'disp_x disp_y disp_z'
cbar = cbar
lambda = 100
mu = 50
e0 = 0.02
[]
[]
[]
[TensorTimeIntegrators]
[c]
type = FFTSemiImplicit
buffer = c
reciprocal_buffer = cbar
linear_reciprocal = kappabarbar
nonlinear_reciprocal = Mbarmubar
[]
[]
[Postprocessors]
[min_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[min_disp_x]
type = TensorExtremeValuePostprocessor
buffer = disp_x
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_disp_x]
type = TensorExtremeValuePostprocessor
buffer = disp_x
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[min_disp_y]
type = TensorExtremeValuePostprocessor
buffer = disp_y
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_disp_y]
type = TensorExtremeValuePostprocessor
buffer = disp_y
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[min_disp_z]
type = TensorExtremeValuePostprocessor
buffer = disp_z
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_disp_z]
type = TensorExtremeValuePostprocessor
buffer = disp_z
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[C]
type = TensorIntegralPostprocessor
buffer = c
execute_on = 'TIMESTEP_END'
[]
[cavg]
type = TensorAveragePostprocessor
buffer = c
execute_on = 'TIMESTEP_END'
[]
[]
[Problem]
type = TensorProblem
spectral_solve_substeps = 1000
print_debug_output = true
[]
[Executioner]
type = Transient
num_steps = 100
[TimeStepper]
type = IterationAdaptiveDT
growth_factor = 1.8
dt = 0.1
[]
dtmax = 1000
[]
[Outputs]
csv = true
perf_graph = true
execute_on = 'TIMESTEP_END'
[]
(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
[]
[]
#### 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'
[]
[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
residual_tensor = speed
[]
[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
[]
[]
(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
[]
[]
[TensorSolver]
type = LBMStream
buffer = 'h f'
f_old = 'h_post_collision f_post_collision'
[]
[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
residual_tensor = speed
[]
[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/real_space/laplacian.i)
[Domain]
dim = 2
nx = 80
ny = 60
parallel_mode = REAL_SPACE
periodic_directions = 'X Y'
xmin = -1
ymin = -1
xmax = 1
ymax = 1
[]
[TensorBuffers]
[u]
type = PlainTensorBuffer
[]
[u_gold]
type = PlainTensorBuffer
[]
[lap]
type = PlainTensorBuffer
[]
[]
[TensorComputes]
[Initialize]
[u_ic]
type = ParsedCompute
buffer = u
expression = 'sin(pi*(x+0.5))*sin(pi*(y+0.5))'
extra_symbols = true
[]
[lap_gold_ic]
type = ParsedCompute
buffer = lap_gold
#expression = '0'
expression = '-2*pi*pi*sin(pi*(x+0.5))*sin(pi*(y+0.5))'
extra_symbols = true
expand = REAL
[]
[]
[Solve]
[lap_compute]
type = FiniteDifferenceLaplacian
input = u
buffer = lap
[]
[diff]
type = ParsedCompute
buffer = diff
expression = 'abs(lap - lap_gold)'
inputs = 'lap lap_gold'
[]
[]
[]
[Postprocessors]
[max_error]
type = TensorExtremeValuePostprocessor
buffer = diff
value_type = MAX
[]
[]
[TensorOutputs]
active = ''
[lap]
type = XDMFTensorOutput
buffer = 'lap lap_gold diff'
enable_hdf5 = true
[]
[]
[Executioner]
type = Transient
num_steps = 1
[]
[Outputs]
[out]
type = CSV
[]
[]
(test/tests/tensor_compute/coupled_pf_mech_secant.i)
[Domain]
dim = 3
nx = 128
ny = 128
nz = 128
xmax = ${fparse pi*4}
ymax = ${fparse pi*4}
zmax = ${fparse pi*4}
mesh_mode = DUMMY
[]
[TensorBuffers]
# phase field
[c]
[]
[cbar]
[]
[mu]
[]
[mubar]
[]
[Mbarmubar]
[]
# mechanics
[disp_x]
[]
[disp_y]
[]
[disp_z]
[]
[mumechbar]
[]
[mumech]
[]
# constant tensors
[Mbar]
[]
[kappabarbar]
[]
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'c disp_x disp_y disp_z mu mumech'
output_mode = 'Node Node Node Node Cell Cell'
enable_hdf5 = true
[]
[]
[TensorComputes]
[Initialize]
[c]
type = RandomTensor
buffer = c
min = 0.44
max = 0.56
[]
[mu_init]
type = ConstantTensor
buffer = mu
[]
[mumech_init]
type = ConstantTensor
buffer = mumech
[]
[disp_x]
type = RandomTensor
buffer = disp_x
min = 0
max = 0
[]
[disp_y]
type = RandomTensor
buffer = disp_y
min = 0
max = 0
[]
[disp_z]
type = RandomTensor
buffer = disp_z
min = 0
max = 0
[]
[Mbar]
type = ReciprocalLaplacianFactor
factor = 0.2 # Mobility
buffer = Mbar
[]
[kappabarbar]
type = ReciprocalLaplacianSquareFactor
factor = -0.001 # kappa
buffer = kappabarbar
[]
[]
[Solve]
[mu]
# chemical potential (real space)
type = ParsedCompute
buffer = mu
expression = '0.1*c^2*(c-1)^2' # + c*sin(x/2)*0.005'
extra_symbols = true
derivatives = c
inputs = c
[]
[mubar]
# chemical potential (reciprocal space)
type = ForwardFFT
buffer = mubar
input = mu
[]
[mumechbar]
# mechanical chemical potential (reciprocal space)
type = FFTElasticChemicalPotential
buffer = mumechbar
cbar = cbar
displacements = 'disp_x disp_y disp_z'
lambda = 100
mu = 50
e0 = 0.02
[]
[mumech]
# chemical potential (reciprocal space)
type = InverseFFT
buffer = mumech
input = mumechbar
[]
[Mbarmubar]
type = ParsedCompute
buffer = Mbarmubar
expression = 'Mbar*(mubar+mumechbar)'
inputs = 'Mbar mubar mumechbar'
[]
[cbar]
type = ForwardFFT
buffer = cbar
input = c
[]
[qsmech]
type = FFTQuasistaticElasticity
displacements = 'disp_x disp_y disp_z'
cbar = cbar
lambda = 100
mu = 50
e0 = 0.02
[]
[]
[]
[TensorSolver]
type = SecantSolver
substeps = 1
max_iterations = 1000
# damping = 0.75
relative_tolerance = 1e-6
absolute_tolerance = 1e-6
buffer = c
dt_epsilon = 1e-7
reciprocal_buffer = cbar
linear_reciprocal = kappabarbar
nonlinear_reciprocal = Mbarmubar
verbose = true
[]
[Postprocessors]
[min_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[min_disp_x]
type = TensorExtremeValuePostprocessor
buffer = disp_x
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_disp_x]
type = TensorExtremeValuePostprocessor
buffer = disp_x
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[min_disp_y]
type = TensorExtremeValuePostprocessor
buffer = disp_y
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_disp_y]
type = TensorExtremeValuePostprocessor
buffer = disp_y
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[min_disp_z]
type = TensorExtremeValuePostprocessor
buffer = disp_z
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_disp_z]
type = TensorExtremeValuePostprocessor
buffer = disp_z
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[C]
type = TensorIntegralPostprocessor
buffer = c
execute_on = 'TIMESTEP_END'
[]
[cavg]
type = TensorAveragePostprocessor
buffer = c
execute_on = 'TIMESTEP_END'
[]
[]
[Problem]
type = TensorProblem
[]
[Executioner]
type = Transient
end_time = 100
[TimeStepper]
type = TensorSolveIterationAdaptiveDT
dt = 0.1
max_iterations = 500
min_iterations = 300
growth_factor = 1.1
cutback_factor = 0.9
[]
[]
[Outputs]
csv = true
perf_graph = true
execute_on = 'TIMESTEP_END'
[]
(examples/libtorch_kks/KKS_libtorch.i)
#
# Kim-Kim-Suzuki with Gibbs energy supplied by a torch model, solved on a 2D grid.
#
# Constants for Initial Conditions
r = 30
l = 4.2
# Initial condition function for order parameter
eta_IC = '0.5*(1-tanh(2*(sqrt(x^2+y^2)-${r})/${l}))'
# Phase-field model parameters
kappa_eta = 5
w = 1
M = 5
L = 5
# Expressions for switching function and bulk Gibbs energy
h_eta = 'eta^3*(6*eta^2-15*eta+10)'
[Domain]
dim = 2
nx = 100
ny = 100
xmin = -50
xmax = 50
ymin = -50
ymax = 50
# automatically create a matching mesh
mesh_mode = DUMMY
[]
[TensorComputes]
[Initialize]
[c_IC]
type = ParsedCompute
buffer = c
expression = '0.6 + (0.3-0.6)*${eta_IC}'
extra_symbols = 'true'
[]
[eta_IC]
type = ParsedCompute
buffer = eta
expression = '${eta_IC}'
extra_symbols = 'true'
[]
[psi_init]
type = ConstantTensor
buffer = psi
real = 1
[]
[M]
type = ConstantTensor
buffer = M
real = ${M}
[]
[L]
type = ConstantTensor
buffer = L
real = ${L}
[]
[L_kappa]
type = ReciprocalLaplacianFactor
buffer = L_kappa
factor = ${fparse ${L} * ${kappa_eta} }
[]
[h_eta_IC]
type = ParsedCompute
buffer = h_eta
expression = '${h_eta}'
inputs = eta
[]
[G_func_IC]
type = LibtorchGibbsEnergy
buffer = 'G'
phase_fractions = 'h_eta'
concentrations = 'c'
domega_detas = 'dG_dh'
chem_pots = 'mu'
libtorch_model_file = 'torch_NN_gibbs_model.pt'
[]
[smooth]
type = DeAliasingTensor
method = HOULI
buffer = smooth
[]
[]
[Solve]
[h_eta]
type = ParsedCompute
buffer = h_eta
expression = '${h_eta}'
inputs = eta
[]
[G_func]
type = LibtorchGibbsEnergy
buffer = 'G'
phase_fractions = 'h_eta'
concentrations = 'c'
domega_detas = 'dG_dh'
chem_pots = 'mu'
libtorch_model_file = 'torch_NN_gibbs_model.pt'
[]
[dG_deta]
type = ParsedCompute
buffer = 'dG_deta'
inputs = 'eta dG_dh'
expression = 'dG_dh * ${h_eta} + ${w} * eta^2 * (1-eta^2)^2'
derivatives = 'eta'
[]
[etabar]
type = ForwardFFT
buffer = etabar
input = eta
[]
[AC_bulk]
type = ReciprocalAllenCahn
L = L
buffer = AC_bulk
dF_chem_deta = dG_deta
psi = psi
[]
[NL_eta]
type = ParsedCompute
buffer = NL_eta
expression = 'AC_bulk '
inputs = 'AC_bulk'
[]
[cbar]
type = ForwardFFT
buffer = cbar
input = c
[]
[div_J]
type = ReciprocalMatDiffusion
buffer = 'div_J'
chemical_potential = mu
mobility = M
psi = psi
[]
[NL_c]
type = ParsedCompute
buffer = 'NL_c'
inputs = 'div_J smooth'
expression = 'smooth * div_J'
[]
[]
[]
[TensorSolver]
type = AdamsBashforthMoulton
buffer = 'c eta'
reciprocal_buffer = 'cbar etabar'
linear_reciprocal = '0 L_kappa'
nonlinear_reciprocal = 'NL_c NL_eta'
substeps = 1e3
predictor_order = 3
corrector_order = 1
corrector_steps = 1
[]
[Postprocessors]
[total_c]
type = TensorIntegralPostprocessor
buffer = c
[]
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'eta c mu psi dG_deta dG_dh G'
enable_hdf5 = true
transpose = false
[]
[]
[Executioner]
type = Transient
num_steps = 100
[TimeStepper]
type = IterationAdaptiveDT
growth_factor = 1.25
dt = 0.1
[]
dtmax = 10
[]
[Outputs]
csv = true
perf_graph = true
execute_on = 'INITIAL TIMESTEP_END'
[]
(examples/phase_field_crystal/pfc_fcc_atomic_structure.i)
# FCC PFC test designed to clearly show atomic-scale density fluctuations
# Using larger epsilon for stronger density modulations
# Based on PhysRevE.81.061601 two-mode FCC PFC model
# Grid parameters - fine resolution to resolve atomic peaks
N = 4 # Grid points per dimension
# PFC model parameters - LARGER epsilon for clearer atomic structure
psi_mean = -0.2 # Mean density
epsilon = 0.15 # Larger undercooling for stronger modulations
R1_param = 0.0 # Two-mode coupling (0 for max FCC stability)
Q1_param = ${fparse sqrt(4.0/3.0)} # FCC wave number ratio
# Domain size: want several unit cells but fine enough resolution
# The (111) wavelength is 2*pi/q0 = 2*pi*sqrt(3/4) ≈ 5.44 (dimensionless)
# Let's fit about 'N' wavelengths in each direction
Lx = ${fparse N * 2 * pi / ${Q1_param}}
Ly = ${Lx}
[Domain]
dim = 2
nx = ${fparse ${N} * 8 }
ny = ${fparse ${N} * 8 }
xmax = ${Lx}
ymax = ${Ly}
[]
[TensorComputes]
[Initialize]
# Random perturbations to seed crystal growth
[psi]
type = RandomTensor
buffer = psi
max = ${fparse ${psi_mean} + 0.01}
min = ${fparse ${psi_mean} - 0.01}
seed = 12345
[]
# Linear operator for FCC
[linear]
type = FCCPFCLinear
buffer = 'linear'
eps = ${epsilon}
R1 = ${R1_param}
Q1 = ${Q1_param}
mobility = 1.0
[]
# Dealiasing for cubic nonlinearity
[smooth_operator]
type = DeAliasingTensor
buffer = smooth_operator
method = HOULI
[]
[]
[Solve]
# Nonlinear term
[nl_div_psi_cubed]
type = FCCPFCNonlinear
buffer = NL
psi = psi
dealiasing = smooth_operator
mobility = 1.0
[]
# FFT for spectral solver
[psi_hat]
type = ForwardFFT
buffer = psi_hat
input = psi
[]
[]
[]
[TensorSolver]
type = AdamsBashforthMoulton
buffer = 'psi'
linear_reciprocal = 'linear'
nonlinear_reciprocal = 'NL'
reciprocal_buffer = 'psi_hat'
corrector_order = 1
corrector_steps = 3
predictor_order = 1
substeps = 10000
[]
[Postprocessors]
[max_psi]
type = TensorExtremeValuePostprocessor
buffer = psi
value_type = MAX
execute_on = 'initial timestep_end'
[]
[min_psi]
type = TensorExtremeValuePostprocessor
buffer = psi
value_type = MIN
execute_on = 'initial timestep_end'
[]
[mean_psi]
type = TensorAveragePostprocessor
buffer = psi
execute_on = 'initial timestep_end'
[]
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'psi'
enable_hdf5 = true
[]
[]
[Executioner]
type = Transient
num_steps = 1000
dt = 20 # Larger timestep for faster evolution
[]
[Outputs]
perf_graph = true
csv = true
[]
(test/tests/typed_tensors/gradient.i)
[Domain]
dim = 3
nx = 20
ny = 10
nz = 5
mesh_mode = DUMMY
device_names = cpu
[]
[TensorComputes]
[Initialize]
[c]
type = ParsedCompute
buffer = c
extra_symbols = true
expression = 'sin(x*8*pi)+cos(y*4*pi)+sin(z*2*pi)'
[]
[grad_c]
type = NEML2GradientVector
buffer = grad_c
input = c
[]
[]
[]
[Problem]
type = TensorProblem
print_debug_output = true
[]
[Executioner]
type = Transient
num_steps = 1
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'c grad_c'
file_base = gradient
output_mode = 'NODE NODE'
enable_hdf5 = true
[]
[]
(test/tests/mechanics/mech3d.i)
[Domain]
dim = 3
nx = 16
ny = 16
nz = 16
xmax = ${fparse 2*pi}
ymax = ${fparse 2*pi}
zmax = ${fparse 2*pi}
mesh_mode = DUMMY
parallel_mode = FFT_AUTO
[]
[TensorComputes]
[Initialize]
[phase]
type = ParsedCompute
expression = '(cos(x)/2+0.5)^1*(cos(y)/2+0.5)^1*(cos(z)/2+0.5)^1'
extra_symbols = true
buffer = phase
[]
[K]
type = ParsedCompute
buffer = K
expression = '(1-phase)*Ka + phase*Kb'
inputs = phase
constant_names = 'Ka Kb'
constant_expressions = '1 10'
[]
[mu]
type = ParsedCompute
buffer = mu
expression = '(1-phase)*mua + phase*mub'
inputs = phase
constant_names = 'mua mub'
constant_expressions = '0.5 5'
[]
[Finit]
type = RankTwoIdentity
buffer = F
[]
[]
[Solve]
[hyper_elasticity]
type = HyperElasticIsotropic
buffer = stress
F = Fnew
K = K
mu = mu
[]
[root]
[applied_strain]
type = MacroscopicShearTensor
buffer = applied_strain
[]
[mech]
type = FFTMechanics
buffer = Fnew
F = F
K = K
mu = mu
l_tol = 1e-2
nl_rel_tol = 2e-2
nl_abs_tol = 2e-2
constitutive_model = hyper_elasticity
stress = stress
applied_macroscopic_strain = applied_strain
[]
[]
[]
[Postprocess]
[displacements]
type = ComputeDisplacements
buffer = disp
F = F
[]
[vonmises]
type = ComputeVonMisesStress
buffer = sV
[]
[]
[]
[TensorSolver]
# no variables are integrated by this solver (FFTMechanics performs a steady state mechanics solve)
type = ForwardEulerSolver
root_compute = root
# deformation tensor is just forwarded Fnew -> F
forward_buffer = F
forward_buffer_new = Fnew
substeps = 10
[]
[TensorOutputs]
active = 'parallel' # serial
[serial]
type = XDMFTensorOutput
file_base = mech3d_serial
buffer = 'disp sV F phase'
output_mode = 'OVERSIZED_NODAL CELL CELL NODE'
enable_hdf5 = true
execute_on = 'TIMESTEP_END'
[]
[parallel]
# NODE and OVERSIZED_NODAL are not yet available in parallel output
type = XDMFTensorOutput
file_base = mech3d_parallel
buffer = 'sV F phase'
output_mode = 'CELL CELL CELL'
enable_hdf5 = true
execute_on = 'TIMESTEP_END'
[]
[]
[Executioner]
type = Transient
num_steps = 3
dt = 0.01
[]
(examples/cahn_hilliard/cahnhilliard2.i)
#
# The same simple Cahn-Hilliard solve as cahnhilliard.i, but on a 3D grid
# and using the faster TensorOutputs system.
#
[Domain]
dim = 3
nx = 200
ny = 200
nz = 200
xmax = ${fparse pi*8}
ymax = ${fparse pi*8}
zmax = ${fparse pi*8}
# run on a CUDA device (adjust this to `cpu` if not available)
device_names = 'cuda'
# create a single element dummy mesh. Output will use the custom XDMF output
# in the `TensorOutputs` system.
mesh_mode = DUMMY
[]
[TensorBuffers]
[c]
[]
[cbar]
[]
[mu]
[]
[mubar]
[]
[Mbarmubar]
[]
# constant tensors
[Mbar]
[]
[kappabarbar]
[]
[]
[TensorOutputs]
# the TensorOutouts system supports asynchronous threaded output.
# for GOU calculations a copy of the solution fields is moved to the CPU,
# and while the output files are written the next time step is already
# starting to compute.
[xdmf]
type = XDMFTensorOutput
buffer = 'c mu'
enable_hdf5 = true
[]
[]
[TensorComputes]
[Initialize]
[c]
type = RandomTensor
buffer = c
min = 0.44
max = 0.56
[]
[Mbar]
type = ReciprocalLaplacianFactor
factor = 0.2 # Mobility
buffer = Mbar
[]
[kappabarbar]
type = ReciprocalLaplacianSquareFactor
factor = -0.001 # kappa
buffer = kappabarbar
[]
[]
[Solve]
[mu]
type = ParsedCompute
buffer = mu
expression = '0.1*c^2*(c-1)^2'
derivatives = c
inputs = c
[]
[mubar]
type = ForwardFFT
buffer = mubar
input = mu
[]
[Mbarmubar]
type = ParsedCompute
buffer = Mbarmubar
expression = 'Mbar*mubar'
inputs = 'Mbar mubar'
[]
[cbar]
type = ForwardFFT
buffer = cbar
input = c
[]
[]
[]
[TensorSolver]
type = AdamsBashforthMoulton
buffer = c
reciprocal_buffer = cbar
linear_reciprocal = kappabarbar
nonlinear_reciprocal = Mbarmubar
substeps = 1000
[]
[Postprocessors]
[min_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[C]
type = TensorIntegralPostprocessor
buffer = c
execute_on = 'TIMESTEP_END'
[]
[cavg]
type = TensorAveragePostprocessor
buffer = c
execute_on = 'TIMESTEP_END'
[]
[]
[Problem]
type = TensorProblem
[]
[Executioner]
type = Transient
num_steps = 100
[TimeStepper]
type = IterationAdaptiveDT
growth_factor = 1.8
dt = 0.1
[]
dtmax = 1000
[]
[Outputs]
csv = true
perf_graph = true
execute_on = 'TIMESTEP_END'
[]
(benchmarks/02_oswald_ripening/2a.i)
[Domain]
dim = 2
nx = 200
ny = 200
xmax = 200
ymax = 200
mesh_mode = DOMAIN
[]
fchem = 'fa:=rho^2*(c-ca)^2;
fb:=rho^2*(cb-c)^2;
h:=n1^3*(6*n1^2-15*n1+10) +
n2^3*(6*n2^2-15*n2+10) +
n3^3*(6*n3^2-15*n3+10) +
n4^3*(6*n4^2-15*n4+10);
g:=n1^2*(1-n1)^2 +
n2^2*(1-n2)^2 +
n3^2*(1-n3)^2 +
n4^2*(1-n4)^2 +
alpha*(
n1^2*n2^2 + n1^2*n3^2 + n1^2*n4^2 +
n2^2*n1^2 + n2^2*n3^2 + n2^2*n4^2 +
n3^2*n1^2 + n3^2*n2^2 + n3^2*n4^2 +
n4^2*n1^2 + n4^2*n2^2 + n4^2*n3^2);
(fa*(1-h) + fb*h + w*g)'
nic = 'epsilon*(cos((0.01*idx)*x-4)*cos((0.007+0.01*idx)*y)
+cos((0.11+0.01*idx)*x)*cos((0.11+0.01*idx)*y)
+psi*(cos((0.046+0.001*idx)*x+(0.0405+0.001*idx)*y)
*cos((0.031+0.001*idx)*x-(0.004+0.001*idx)*y))^2)^2'
cnames = 'rho ca cb alpha w L M'
cvalues = 'sqrt(2) 0.3 0.7 5 1 5 5'
[TensorBuffers]
# variables
[c]
# map_to_aux_variable = c
[]
[n1]
[]
[n2]
[]
[n3]
[]
[n4]
[]
[c_bar]
[]
[n1_bar]
[]
[n2_bar]
[]
[n3_bar]
[]
[n4_bar]
[]
[mu_c]
# map_to_aux_variable = mu
[]
[mu_n1]
[]
[mu_n2]
[]
[mu_n3]
[]
[mu_n4]
[]
[mu_c_bar]
[]
[mu_n1_bar]
[]
[mu_n2_bar]
[]
[mu_n3_bar]
[]
[mu_n4_bar]
[]
[Mbar_mu_c_bar]
[]
# constant tensors
[Lbar] # FFT(M*laplacian)
[]
[MkappaL2bar] # FFT(-M*kappa*laplacian^2)
[]
[kappaLbar] # FFT(L*kappa*laplacian)
[]
# postprocessing
[F]
[]
[Fgrad_c]
[]
[Fgrad_n1]
[]
[Fgrad_n2]
[]
[Fgrad_n3]
[]
[Fgrad_n4]
[]
[bnds]
#map_to_aux_variable = bnds
[]
[]
[TensorComputes]
[Initialize]
[c]
type = ParsedCompute
buffer = c
extra_symbols = true
expression = 'c0+epsilon*(cos(0.105*x)*cos(0.11*y)+(cos(0.13*x)*cos(0.087*y))^2+cos(0.025*x-0.15*y)*cos(0.07*x-0.02*y))'
constant_names = 'c0 epsilon'
constant_expressions = '0.5 0.01'
[]
[Lbar]
type = ReciprocalLaplacianFactor
# Mobility is pulled into the chemical potential below
buffer = Lbar
[]
[MkappaL2bar]
type = ReciprocalLaplacianSquareFactor
factor = -15 # -kappa_c*M
buffer = MkappaL2bar
[]
[kappaLbar]
type = ReciprocalLaplacianFactor
buffer = kappaLbar
factor = 15 # kappa_ni*L
[]
[n1]
type = ParsedCompute
buffer = n1
expression = ${nic}
extra_symbols = true
constant_names = 'idx epsilon psi'
constant_expressions = ' 1 0.1 1.5'
[]
[n2]
type = ParsedCompute
buffer = n2
expression = ${nic}
extra_symbols = true
constant_names = 'idx epsilon psi'
constant_expressions = ' 2 0.1 1.5'
[]
[n3]
type = ParsedCompute
buffer = n3
expression = ${nic}
extra_symbols = true
constant_names = 'idx epsilon psi'
constant_expressions = ' 3 0.1 1.5'
[]
[n4]
type = ParsedCompute
buffer = n4
expression = ${nic}
extra_symbols = true
constant_names = 'idx epsilon psi'
constant_expressions = ' 4 0.1 1.5'
[]
[]
[Solve]
[mu_c]
type = ParsedCompute
buffer = mu_c
expression = '${fchem}*M'
constant_names = ${cnames}
constant_expressions = ${cvalues}
derivatives = c
inputs = 'c n1 n2 n3 n4'
[]
[mu_n1]
type = ParsedCompute
buffer = mu_n1
expression = '${fchem}*(-L)'
constant_names = ${cnames}
constant_expressions = ${cvalues}
derivatives = n1
inputs = 'c n1 n2 n3 n4'
[]
[mu_n2]
type = ParsedCompute
buffer = mu_n2
expression = '${fchem}*(-L)'
constant_names = ${cnames}
constant_expressions = ${cvalues}
derivatives = n2
inputs = 'c n1 n2 n3 n4'
[]
[mu_n3]
type = ParsedCompute
buffer = mu_n3
expression = '${fchem}*(-L)'
constant_names = ${cnames}
constant_expressions = ${cvalues}
derivatives = n3
inputs = 'c n1 n2 n3 n4'
[]
[mu_n4]
type = ParsedCompute
buffer = mu_n4
expression = '${fchem}*(-L)'
constant_names = ${cnames}
constant_expressions = ${cvalues}
derivatives = n4
inputs = 'c n1 n2 n3 n4'
[]
[mu_c_bar]
type = ForwardFFT
buffer = mu_c_bar
input = mu_c
[]
[mu_n1_bar]
type = ForwardFFT
buffer = mu_n1_bar
input = mu_n1
[]
[mu_n2_bar]
type = ForwardFFT
buffer = mu_n2_bar
input = mu_n2
[]
[mu_n3_bar]
type = ForwardFFT
buffer = mu_n3_bar
input = mu_n3
[]
[mu_n4_bar]
type = ForwardFFT
buffer = mu_n4_bar
input = mu_n4
[]
[Mbar_mu_c_bar]
type = ParsedCompute
buffer = Mbar_mu_c_bar
expression = 'Lbar*mu_c_bar'
inputs = 'Lbar mu_c_bar'
[]
[c_bar]
type = ForwardFFT
buffer = c_bar
input = c
[]
[n1_bar]
type = ForwardFFT
buffer = n1_bar
input = n1
[]
[n2_bar]
type = ForwardFFT
buffer = n2_bar
input = n2
[]
[n3_bar]
type = ForwardFFT
buffer = n3_bar
input = n3
[]
[n4_bar]
type = ForwardFFT
buffer = n4_bar
input = n4
[]
[]
[Postprocess]
[Fgrad_c]
type = FFTGradientSquare
buffer = Fgrad_c
input = c
factor = 1.5 # kappa/2
[]
[Fgrad_n1]
type = FFTGradientSquare
buffer = Fgrad_n1
input = n1
factor = 1.5 # kappa/2
[]
[Fgrad_n2]
type = FFTGradientSquare
buffer = Fgrad_n2
input = n2
factor = 1.5 # kappa/2
[]
[Fgrad_n3]
type = FFTGradientSquare
buffer = Fgrad_n3
input = n3
factor = 1.5 # kappa/2
[]
[Fgrad_n4]
type = FFTGradientSquare
buffer = Fgrad_n4
input = n4
factor = 1.5 # kappa/2
[]
[F]
type = ParsedCompute
buffer = F
expression = '${fchem} + Fgrad_c + Fgrad_n1 + Fgrad_n2 + Fgrad_n3 + Fgrad_n4'
constant_names = ${cnames}
constant_expressions = ${cvalues}
inputs = 'c n1 n2 n3 n4 Fgrad_c Fgrad_n1 Fgrad_n2 Fgrad_n3 Fgrad_n4'
[]
[bnds]
type = ParsedCompute
buffer = bnds
expression = 'n1^2 + n2^2 + n3^2 + n4^2'
inputs = 'n1 n2 n3 n4'
[]
[]
[]
[TensorSolver]
type = AdamsBashforthMoulton
buffer = 'c n1 n2 n3 n4'
reciprocal_buffer = 'c_bar n1_bar n2_bar n3_bar n4_bar'
linear_reciprocal = 'MkappaL2bar kappaLbar kappaLbar kappaLbar kappaLbar'
nonlinear_reciprocal = 'Mbar_mu_c_bar mu_n1_bar mu_n2_bar mu_n3_bar mu_n4_bar'
substeps = 2000
predictor_order = 2
corrector_order = 2
corrector_steps = 0
[]
[AuxVariables]
[c]
[]
[bnds]
[]
[]
[Postprocessors]
[min_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[F]
type = TensorIntegralPostprocessor
buffer = F
[]
# [stable_dt]
# type = SemiImplicitCriticalTimeStep
# buffer = MkappaL2bar
# []
[]
[Problem]
type = TensorProblem
[]
[Executioner]
type = Transient
num_steps = 1030
[TimeStepper]
type = IterationAdaptiveDT
growth_factor = 1.1
dt = 0.001
[]
dtmax = 10
[]
[Outputs]
csv = true
perf_graph = true
execute_on = 'TIMESTEP_END'
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'c bnds'
output_mode = 'CELL CELL'
[]
[]
(benchmarks/02_oswald_ripening/2a_secant.i)
[Domain]
dim = 2
nx = 120
ny = 120
xmax = 200
ymax = 200
device_names = 'cuda:1'
mesh_mode = DOMAIN
[]
fchem='fa:=rho^2*(c-ca)^2;
fb:=rho^2*(cb-c)^2;
h:=n1^3*(6*n1^2-15*n1+10) +
n2^3*(6*n2^2-15*n2+10) +
n3^3*(6*n3^2-15*n3+10) +
n4^3*(6*n4^2-15*n4+10);
g:=n1^2*(1-n1)^2 +
n2^2*(1-n2)^2 +
n3^2*(1-n3)^2 +
n4^2*(1-n4)^2 +
alpha*(
n1^2*n2^2 + n1^2*n3^2 + n1^2*n4^2 +
n2^2*n1^2 + n2^2*n3^2 + n2^2*n4^2 +
n3^2*n1^2 + n3^2*n2^2 + n3^2*n4^2 +
n4^2*n1^2 + n4^2*n2^2 + n4^2*n3^2);
(fa*(1-h) + fb*h + w*g)'
nic = 'epsilon*(cos((0.01*idx)*x-4)*cos((0.007+0.01*idx)*y)
+cos((0.11+0.01*idx)*x)*cos((0.11+0.01*idx)*y)
+psi*(cos((0.046+0.001*idx)*x+(0.0405+0.001*idx)*y)
*cos((0.031+0.001*idx)*x-(0.004+0.001*idx)*y))^2)^2'
cnames = 'rho ca cb alpha w L M'
cvalues= 'sqrt(2) 0.3 0.7 5 1 5 5'
[TensorBuffers]
# variables
[c]
[]
[n1]
[]
[n2]
[]
[n3]
[]
[n4]
[]
[c_bar]
[]
[n1_bar]
[]
[n2_bar]
[]
[n3_bar]
[]
[n4_bar]
[]
[mu_c]
[]
[mu_n1]
[]
[mu_n2]
[]
[mu_n3]
[]
[mu_n4]
[]
[mu_c_bar]
[]
[mu_n1_bar]
[]
[mu_n2_bar]
[]
[mu_n3_bar]
[]
[mu_n4_bar]
[]
[Mbar_mu_c_bar]
[]
# constant tensors
[Lbar] # FFT(M*laplacian)
[]
[MkappaL2bar] # FFT(-M*kappa*laplacian^2)
[]
[kappaLbar] # FFT(L*kappa*laplacian)
[]
# postprocessing
[F]
[]
[Fgrad_c]
[]
[Fgrad_n1]
[]
[Fgrad_n2]
[]
[Fgrad_n3]
[]
[Fgrad_n4]
[]
[bnds]
[]
[]
[TensorComputes]
[Initialize]
[c]
type = ParsedCompute
buffer = c
extra_symbols = true
expression = 'c0+epsilon*(cos(0.105*x)*cos(0.11*y)+(cos(0.13*x)*cos(0.087*y))^2+cos(0.025*x-0.15*y)*cos(0.07*x-0.02*y))'
constant_names = 'c0 epsilon'
constant_expressions = '0.5 0.01'
[]
[Lbar]
type = ReciprocalLaplacianFactor
# Mobility is pulled into the chemical potential below
buffer = Lbar
[]
[MkappaL2bar]
type = ReciprocalLaplacianSquareFactor
factor = -15 # -kappa_c*M
buffer = MkappaL2bar
[]
[kappaLbar]
type = ReciprocalLaplacianFactor
buffer = kappaLbar
factor = 15 # kappa_ni*L
[]
[n1]
type = ParsedCompute
buffer = n1
expression = ${nic}
extra_symbols = true
constant_names = 'idx epsilon psi'
constant_expressions = ' 1 0.1 1.5'
[]
[n2]
type = ParsedCompute
buffer = n2
expression = ${nic}
extra_symbols = true
constant_names = 'idx epsilon psi'
constant_expressions = ' 2 0.1 1.5'
[]
[n3]
type = ParsedCompute
buffer = n3
expression = ${nic}
extra_symbols = true
constant_names = 'idx epsilon psi'
constant_expressions = ' 3 0.1 1.5'
[]
[n4]
type = ParsedCompute
buffer = n4
expression = ${nic}
extra_symbols = true
constant_names = 'idx epsilon psi'
constant_expressions = ' 4 0.1 1.5'
[]
[]
[Solve]
[mu_c]
type = ParsedCompute
buffer = mu_c
expression = '${fchem}*M'
constant_names = ${cnames}
constant_expressions = ${cvalues}
derivatives = c
inputs = 'c n1 n2 n3 n4'
[]
[mu_n1]
type = ParsedCompute
buffer = mu_n1
expression = '${fchem}*(-L)'
constant_names = ${cnames}
constant_expressions = ${cvalues}
derivatives = n1
inputs = 'c n1 n2 n3 n4'
[]
[mu_n2]
type = ParsedCompute
buffer = mu_n2
expression = '${fchem}*(-L)'
constant_names = ${cnames}
constant_expressions = ${cvalues}
derivatives = n2
inputs = 'c n1 n2 n3 n4'
[]
[mu_n3]
type = ParsedCompute
buffer = mu_n3
expression = '${fchem}*(-L)'
constant_names = ${cnames}
constant_expressions = ${cvalues}
derivatives = n3
inputs = 'c n1 n2 n3 n4'
[]
[mu_n4]
type = ParsedCompute
buffer = mu_n4
expression = '${fchem}*(-L)'
constant_names = ${cnames}
constant_expressions = ${cvalues}
derivatives = n4
inputs = 'c n1 n2 n3 n4'
[]
[mu_c_bar]
type = ForwardFFT
buffer = mu_c_bar
input = mu_c
[]
[mu_n1_bar]
type = ForwardFFT
buffer = mu_n1_bar
input = mu_n1
[]
[mu_n2_bar]
type = ForwardFFT
buffer = mu_n2_bar
input = mu_n2
[]
[mu_n3_bar]
type = ForwardFFT
buffer = mu_n3_bar
input = mu_n3
[]
[mu_n4_bar]
type = ForwardFFT
buffer = mu_n4_bar
input = mu_n4
[]
[Mbar_mu_c_bar]
type = ParsedCompute
buffer = Mbar_mu_c_bar
expression = 'Lbar*mu_c_bar'
inputs = 'Lbar mu_c_bar'
[]
[c_bar]
type = ForwardFFT
buffer = c_bar
input = c
[]
[n1_bar]
type = ForwardFFT
buffer = n1_bar
input = n1
[]
[n2_bar]
type = ForwardFFT
buffer = n2_bar
input = n2
[]
[n3_bar]
type = ForwardFFT
buffer = n3_bar
input = n3
[]
[n4_bar]
type = ForwardFFT
buffer = n4_bar
input = n4
[]
[]
[Postprocess]
[Fgrad_c]
type = FFTGradientSquare
buffer = Fgrad_c
input = c
factor = 1.5 # kappa/2
[]
[Fgrad_n1]
type = FFTGradientSquare
buffer = Fgrad_n1
input = n1
factor = 1.5 # kappa/2
[]
[Fgrad_n2]
type = FFTGradientSquare
buffer = Fgrad_n2
input = n2
factor = 1.5 # kappa/2
[]
[Fgrad_n3]
type = FFTGradientSquare
buffer = Fgrad_n3
input = n3
factor = 1.5 # kappa/2
[]
[Fgrad_n4]
type = FFTGradientSquare
buffer = Fgrad_n4
input = n4
factor = 1.5 # kappa/2
[]
[F]
type = ParsedCompute
buffer = F
expression = '${fchem} + Fgrad_c + Fgrad_n1 + Fgrad_n2 + Fgrad_n3 + Fgrad_n4'
constant_names = ${cnames}
constant_expressions = ${cvalues}
inputs = 'c n1 n2 n3 n4 Fgrad_c Fgrad_n1 Fgrad_n2 Fgrad_n3 Fgrad_n4'
[]
[bnds]
type = ParsedCompute
buffer = bnds
expression = 'n1^2 + n2^2 + n3^2 + n4^2'
inputs = 'n1 n2 n3 n4'
[]
[]
[]
[TensorSolver]
type = SecantSolver
substeps = 1
max_iterations = 1000
# damping = 0.75
relative_tolerance = 1e-6
absolute_tolerance = 1e-6
buffer = 'c n1 n2 n3 n4'
dt_epsilon = 1e-7
reciprocal_buffer = 'c_bar n1_bar n2_bar n3_bar n4_bar'
linear_reciprocal = 'MkappaL2bar kappaLbar kappaLbar kappaLbar kappaLbar'
nonlinear_reciprocal = 'Mbar_mu_c_bar mu_n1_bar mu_n2_bar mu_n3_bar mu_n4_bar'
# verbose = true
[]
[Postprocessors]
[min_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[F]
type = TensorIntegralPostprocessor
buffer = F
[]
[]
[Problem]
type = TensorProblem
[]
[Executioner]
type = Transient
# num_steps = 2
end_time = 1e+5
[TimeStepper]
type = TensorSolveIterationAdaptiveDT
dt = 0.1
max_iterations = 500
min_iterations = 300
growth_factor = 1.1
cutback_factor = 0.9
[]
[]
[Outputs]
csv = true
perf_graph = true
execute_on = 'TIMESTEP_END'
[]
[TensorOutputs]
[c]
type = XDMFTensorOutput
buffer = 'c bnds'
output_mode = 'NODE NODE'
enable_hdf5 = 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
[]
[]
[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
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
residual_tensor = speed
[]
[Executioner]
type = Transient
num_steps = 2000
[]
[TensorOutputs]
[xdmf2]
type = XDMFTensorOutput
buffer = 'rho u'
output_mode = 'Cell Cell'
enable_hdf5 = 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
[]
[]
[TensorSolver]
type = LBMStream
buffer = 'h f'
f_old = 'h_post_collision f_post_collision'
[]
[Problem]
type = LatticeBoltzmannProblem
substeps = 50
print_debug_output = true
scalar_constant_names = 'tau_h D sigma'
scalar_constant_values = '${tau_h} ${D} ${sigma}'
residual_tensor = speed
[]
[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/tensor_compute/group.i)
[Domain]
dim = 3
nx = 128
ny = 128
nz = 128
xmax = ${fparse pi*4}
ymax = ${fparse pi*4}
zmax = ${fparse pi*4}
mesh_mode = DUMMY
[]
[TensorBuffers]
# phase field
[c]
[]
[cbar]
[]
[mu]
[]
[mubar]
[]
[Mbarmubar]
[]
# mechanics
[disp_x]
[]
[disp_y]
[]
[disp_z]
[]
[mumechbar]
[]
[mumech]
[]
# constant tensors
[Mbar]
[]
[kappabarbar]
[]
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'c disp_x disp_y disp_z mu mumech'
output_mode = 'Node Node Node Node Cell Cell'
enable_hdf5 = true
[]
[]
[TensorComputes]
[Initialize]
[c]
type = RandomTensor
buffer = c
min = 0.44
max = 0.56
[]
[disp_x]
type = RandomTensor
buffer = disp_x
min = 0
max = 0
[]
[disp_y]
type = RandomTensor
buffer = disp_y
min = 0
max = 0
[]
[disp_z]
type = RandomTensor
buffer = disp_z
min = 0
max = 0
[]
[Mbar]
type = ReciprocalLaplacianFactor
factor = 0.2 # Mobility
buffer = Mbar
[]
[kappabarbar]
type = ReciprocalLaplacianSquareFactor
factor = -0.001 # kappa
buffer = kappabarbar
[]
[]
[Solve]
[mu]
# chemical potential (real space)
type = ParsedCompute
buffer = mu
expression = '0.1*c^2*(c-1)^2' # + c*sin(x/2)*0.005'
extra_symbols = true
derivatives = c
inputs = c
[]
[mubar]
# chemical potential (reciprocal space)
type = ForwardFFT
buffer = mubar
input = mu
[]
[mumechbar]
# mechanical chemical potential (reciprocal space)
type = FFTElasticChemicalPotential
buffer = mumechbar
cbar = cbar
displacements = 'disp_x disp_y disp_z'
lambda = 100
mu = 50
e0 = 0.02
[]
[mumech]
# chemical potential (reciprocal space)
type = InverseFFT
buffer = mumech
input = mumechbar
[]
[Mbarmubar]
type = ParsedCompute
buffer = Mbarmubar
expression = 'Mbar*(mubar+mumechbar)'
inputs = 'Mbar mubar mumechbar'
[]
[cbar]
type = ForwardFFT
buffer = cbar
input = c
[]
[qsmech]
type = FFTQuasistaticElasticity
displacements = 'disp_x disp_y disp_z'
cbar = cbar
lambda = 100
mu = 50
e0 = 0.02
[]
[group]
type = ComputeGroup
computes = 'cbar mumech mu mubar'
[]
[]
[]
[TensorSolver]
type = AdamsBashforthMoulton
buffer = c
reciprocal_buffer = cbar
linear_reciprocal = kappabarbar
nonlinear_reciprocal = Mbarmubar
[]
[Postprocessors]
[min_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[min_disp_x]
type = TensorExtremeValuePostprocessor
buffer = disp_x
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_disp_x]
type = TensorExtremeValuePostprocessor
buffer = disp_x
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[min_disp_y]
type = TensorExtremeValuePostprocessor
buffer = disp_y
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_disp_y]
type = TensorExtremeValuePostprocessor
buffer = disp_y
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[min_disp_z]
type = TensorExtremeValuePostprocessor
buffer = disp_z
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_disp_z]
type = TensorExtremeValuePostprocessor
buffer = disp_z
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[C]
type = TensorIntegralPostprocessor
buffer = c
execute_on = 'TIMESTEP_END'
[]
[cavg]
type = TensorAveragePostprocessor
buffer = c
execute_on = 'TIMESTEP_END'
[]
[]
[Problem]
type = TensorProblem
spectral_solve_substeps = 1000
print_debug_output = true
[]
[Executioner]
type = Transient
num_steps = 100
[TimeStepper]
type = IterationAdaptiveDT
growth_factor = 1.8
dt = 0.1
[]
dtmax = 1000
[]
[Outputs]
csv = true
perf_graph = true
execute_on = 'TIMESTEP_END'
[]
(test/tests/tensor_compute/parallel.i)
[Domain]
device_names = "cuda:1 cuda:0 cpu"
device_weights = "10 10 1"
parallel_mode = FFT_SLAB
dim = 2
nx = 400
ny = 400
xmax = ${fparse pi*2*30}
ymax = ${fparse pi*2*30}
[]
# [Mesh]
# type = GeneratedMesh
# dim =1
# []
[TensorBuffers]
[psi]
# map_to_aux_variable = psi
[]
[psibar]
[]
[psi3]
[]
[psi3bar]
[]
# constant tensors
[linear]
[]
# output
[filter]
# map_to_aux_variable = filter
[]
[filterbar]
[]
[]
[AuxVariables]
[psi]
[]
[filter]
[]
[]
[TensorComputes]
[Initialize]
[psi]
type = RandomTensor
buffer = psi
min = 0
max = 0.07
[]
[linear]
type = SwiftHohenbergLinear
buffer = linear
alpha = 1
r = 0.025
[]
[]
[Solve]
[psi3]
type = ParsedCompute
buffer = psi3
expression = "0.20*psi^2-psi^3"
inputs = psi
[]
[psibar]
type = ForwardFFT
buffer = psibar
input = psi
[]
[psi3bar]
type = ForwardFFT
buffer = psi3bar
input = psi3
[]
[]
[Postprocess]
[low_pass]
type = ParsedCompute
buffer = filterbar
extra_symbols = true
expression = 'psibar * exp(-k2*10)'
inputs = psibar
[]
[filter]
type = InverseFFT
buffer = filter
input = filterbar
[]
[]
[]
[TensorTimeIntegrators]
[c]
type = FFTSemiImplicit
buffer = psi
reciprocal_buffer = psibar
linear_reciprocal = linear
nonlinear_reciprocal = psi3bar
[]
[]
[Problem]
type = TensorProblem
spectral_solve_substeps = 1000
[]
[Executioner]
type = Transient
num_steps = 300
[TimeStepper]
type = IterationAdaptiveDT
growth_factor = 1.2
dt = 10
[]
dtmax = 1000
[]
[Postprocessors]
[min_psi]
type = TensorExtremeValuePostprocessor
buffer = psi
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_psi]
type = TensorExtremeValuePostprocessor
buffer = psi
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[Psi]
type = TensorIntegralPostprocessor
buffer = psi
[]
[]
# [TensorOutputs]
# [xdmf]
# type = XDMFTensorOutput
# buffer = 'psi'
# output_mode = 'Node'
# enable_hdf5 = true
# []
# []
[Outputs]
exodus = true
perf_graph = true
execute_on = 'TIMESTEP_END'
[]
(test/tests/cahnhilliard/cahnhilliard.i)
#
# Simple Cahn-Hilliard solve on a 2D grid. We create a matching (conforming)
# MOOSE mesh (with one element per FFT grid cell) and project the solution onto
# the MOOSE mesh to utilize the exodus output object.
#
[Domain]
dim = 2
nx = 20
ny = 20
xmax = 3
ymax = 3
mesh_mode = DOMAIN
[]
# In this input we fully trely on implicit TensorBuffer declaration
[TensorComputes]
[Initialize]
[c]
# Random initial condition around a concentration of 1/2
type = RandomTensor
buffer = c
min = 0.44
max = 0.56
seed = 0
[]
[mu_init]
type = ConstantTensor
buffer = mu
[]
# precompute fixed factors for the solve
[Mbar]
type = ReciprocalLaplacianFactor
factor = 0.2 # Mobility
buffer = Mbar
[]
[kappabarbar]
type = ReciprocalLaplacianSquareFactor
factor = -0.001 # kappa
buffer = kappabarbar
[]
[]
[Solve]
[cahn_hilliard]
[mu]
type = ParsedCompute
buffer = mu
expression = '0.1*c^2*(c-1)^2'
derivatives = c
inputs = c
[]
[mubar]
type = ForwardFFT
buffer = mubar
input = mu
[]
[Mbarmubar]
type = ParsedCompute
buffer = Mbarmubar
expression = 'Mbar*mubar'
inputs = 'Mbar mubar'
[]
[cbar]
type = ForwardFFT
buffer = cbar
input = c
[]
[]
[]
[]
[TensorSolver]
type = AdamsBashforthMoulton
root_compute = cahn_hilliard
buffer = c
reciprocal_buffer = cbar
linear_reciprocal = kappabarbar
nonlinear_reciprocal = Mbarmubar
substeps = 10
[]
[AuxVariables]
[mu]
# the mu tensor is projected onto this elemental variable
family = MONOMIAL
order = CONSTANT
[]
[c]
# the c tensor is projected onto this nodal variable
[]
[]
[AuxKernels]
active = ''
[c]
type = ProjectTensorAux
buffer = c
variable = c
execute_on = 'INITIAL TIMESTEP_END'
[]
[mu]
type = ProjectTensorAux
buffer = mu
variable = mu
execute_on = 'INITIAL TIMESTEP_END'
[]
[]
[Postprocessors]
[min_c]
type = SemiImplicitCriticalTimeStep
buffer = kappabarbar
execute_on = 'INITIAL TIMESTEP_END'
[]
[delta_int_c]
type = TensorIntegralChangePostprocessor
buffer = c
[]
[]
[Problem]
type = TensorProblem
[]
[Executioner]
type = Transient
num_steps = 10
dt = 1e-3
[]
[TensorOutputs]
active = ''
[xdmf]
type = XDMFTensorOutput
buffer = 'c mu'
output_mode = 'Node Cell'
enable_hdf5 = true
# Do not transpose output to avoid regolding the test. In practice the default
# of transpose = true should always be used
transpose = false
[]
[xdmf2]
# second output to trigger the hdf5 thread safety error
type = XDMFTensorOutput
buffer = 'c'
output_mode = 'Cell'
enable_hdf5 = true
# Do not transpose output to avoid regolding the test. In practice the default
# of transpose = true should always be used
transpose = false
[]
[]
[Outputs]
exodus = true
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
[]
[]
#### 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'
[]
[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
residual_tensor = speed
[]
[Executioner]
type = Transient
num_steps = 1000
[]
[TensorOutputs]
[xdmf2]
type = XDMFTensorOutput
buffer = 'T density velocity'
output_mode = 'Cell Cell Cell'
enable_hdf5 = true
[]
[]
(examples/cahn_hilliard/cahnhilliard4.i)
[Domain]
dim = 3
nx = 100
ny = 100
nz = 100
xmax = ${fparse pi*4}
ymax = ${fparse pi*4}
zmax = ${fparse pi*4}
device_names = 'cuda'
mesh_mode = DUMMY
debug = true
[]
[TensorBuffers]
[c]
[]
[cbar]
[]
[mu]
[]
[mubar]
[]
[Mbarmubar]
[]
# constant tensors
[Mbar]
[]
[kappabarbar]
[]
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'c mu'
output_mode = 'Node Cell'
enable_hdf5 = true
[]
[]
[TensorComputes]
[Initialize]
[c]
type = RandomTensor
buffer = c
min = 0.44
max = 0.56
[]
[Mbar]
type = ReciprocalLaplacianFactor
factor = 0.2 # Mobility
buffer = Mbar
[]
[kappabarbar]
type = ReciprocalLaplacianSquareFactor
factor = -0.001 # kappa
buffer = kappabarbar
[]
[]
[Solve]
[mu]
type = ParsedCompute
buffer = mu
expression = '0.1*c^2*(c-1)^2 + c*sin(x/2)*0.005'
extra_symbols = true
derivatives = c
# expression = "0.4*c^3-0.6*c^2+0.2*c"
inputs = c
[]
[mubar]
type = ForwardFFT
buffer = mubar
input = mu
[]
[Mbarmubar]
type = ParsedCompute
buffer = Mbarmubar
expression = 'Mbar*mubar'
inputs = 'Mbar mubar'
[]
[cbar]
type = ForwardFFT
buffer = cbar
input = c
[]
[]
[]
[TensorTimeIntegrators]
[c]
type = FFTSemiImplicit
buffer = c
reciprocal_buffer = cbar
linear_reciprocal = kappabarbar
nonlinear_reciprocal = Mbarmubar
[]
[]
[Postprocessors]
[min_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[C]
type = TensorIntegralPostprocessor
buffer = c
execute_on = 'TIMESTEP_END'
[]
[cavg]
type = TensorAveragePostprocessor
buffer = c
execute_on = 'TIMESTEP_END'
[]
[]
[Problem]
type = TensorProblem
spectral_solve_substeps = 1000
[]
[Executioner]
type = Transient
num_steps = 50
[TimeStepper]
type = IterationAdaptiveDT
growth_factor = 1.8
dt = 0.1
[]
dtmax = 500
[]
[Outputs]
csv = true
perf_graph = true
execute_on = 'TIMESTEP_END'
[]
(test/tests/typed_tensors/gradient_vector.i)
[Domain]
dim = 3
nx = 20
ny = 10
nz = 5
mesh_mode = DUMMY
device_names = cpu
[]
[TensorComputes]
[Initialize]
[c]
type = ParsedCompute
buffer = c
extra_symbols = true
expression = 'sin(x*8*pi)+cos(y*4*pi)+sin(z*2*pi)'
[]
[grad_c]
type = GradientVector
buffer = grad_c
input = c
[]
[]
[]
[Problem]
type = TensorProblem
print_debug_output = true
[]
[Executioner]
type = Transient
num_steps = 1
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'c grad_c'
file_base = gradient
output_mode = 'NODE NODE'
enable_hdf5 = true
[]
[]
(examples/anisotropic_grain_growth/isotropic_five_grain.i)
interface_width = 1.6
side = 80
gbe_max = '${units 0.65 J/m^2}'
L = '${fparse 1.0 * 1.6 / ${interface_width} }'
# ------------------------------------------------------------------
# Isotropic reference case.
#
# The GB energy is fixed at sigma = 0.65 J/m^2 for all boundaries,
# equal to gbe_max. Consequently gamma_ij = 1.5 for all pairs (the
# symmetric-well value), there are no torque terms, the mobility is
# uniform (L_NL == L), and the mobility-split correction
# (L_NL - L)*kappa*laplacian(gr_i) vanishes identically.
# ------------------------------------------------------------------
# ------------------------------------------------------------------
# Periodic Voronoi-like IC for 5 grains.
#
# Each grain is seeded at a point (xi, yi). We use a smooth *periodic*
# distance metric based on sine functions,
# d_i^2 = (side/pi)^2 * [ sin(pi*(x-xi)/side)^2 + sin(pi*(y-yi)/side)^2 ]
# which reduces to the Euclidean distance near the seed and wraps
# correctly across the periodic boundaries.
#
# The order parameters are then a softmax over the negative distances,
# gr_i = exp(-a*d_i) / sum_j exp(-a*d_j), a = 4/interface_width
# For any two competing grains this reduces exactly to the usual
# 0.5*(1 - tanh(2*(d_i - d_j)/interface_width)) profile, so the
# interfaces come out with the correct width, the partition of unity
# sum_i gr_i = 1 holds identically, and the cell structure is the
# (diffuse) Voronoi tessellation of the seeds under the periodic metric.
# ------------------------------------------------------------------
aexp = '${fparse 4.0 / ${interface_width} }'
lp = '${fparse ${side} / pi }'
kp = '${fparse pi / ${side} }'
# seed points (spread out over the periodic domain)
d0 = 'sqrt((${lp}*sin(${kp}*(x-8)))^2 + (${lp}*sin(${kp}*(y-6)))^2)'
d1 = 'sqrt((${lp}*sin(${kp}*(x-28)))^2 + (${lp}*sin(${kp}*(y-10)))^2)'
d2 = 'sqrt((${lp}*sin(${kp}*(x-20)))^2 + (${lp}*sin(${kp}*(y-24)))^2)'
d3 = 'sqrt((${lp}*sin(${kp}*(x-4)))^2 + (${lp}*sin(${kp}*(y-28)))^2)'
d4 = 'sqrt((${lp}*sin(${kp}*(x-34)))^2 + (${lp}*sin(${kp}*(y-33)))^2)'
g_gamma0 = '${fparse sqrt(2) / 3 }' # g(gamma=1.5)
f0_gamma0 = ${fparse 1/8}
kappa = '${fparse gbe_max * interface_width * sqrt(f0_gamma0) / g_gamma0 }'
mu = '${fparse gbe_max / (g_gamma0 * interface_width * sqrt(f0_gamma0) ) }'
gamma = 1.5
f0 = '((gr0^4/4 - gr0^2/2) + (gr1^4/4 - gr1^2/2) + (gr2^4/4 - gr2^2/2) + (gr3^4/4 - gr3^2/2) + (gr4^4/4 - gr4^2/2) + ${gamma}*(gr0^2*gr1^2 + gr0^2*gr2^2 + gr0^2*gr3^2 + gr0^2*gr4^2 + gr1^2*gr2^2 + gr1^2*gr3^2 + gr1^2*gr4^2 + gr2^2*gr3^2 + gr2^2*gr4^2 + gr3^2*gr4^2) + 0.25)'
[Domain]
dim = 2
nx = 200
ny = 200
xmax = ${side}
ymax = ${side}
mesh_mode = DUMMY
[]
[TensorComputes]
[Initialize]
# unnormalized softmax weights
[c0]
type = ParsedCompute
buffer = c0
expression = 'exp(-${aexp}*${d0})'
extra_symbols = true
[]
[c1]
type = ParsedCompute
buffer = c1
expression = 'exp(-${aexp}*${d1})'
extra_symbols = true
[]
[c2]
type = ParsedCompute
buffer = c2
expression = 'exp(-${aexp}*${d2})'
extra_symbols = true
[]
[c3]
type = ParsedCompute
buffer = c3
expression = 'exp(-${aexp}*${d3})'
extra_symbols = true
[]
[c4]
type = ParsedCompute
buffer = c4
expression = 'exp(-${aexp}*${d4})'
extra_symbols = true
[]
# normalized order parameters (partition of unity by construction)
[gr0]
type = ParsedCompute
buffer = gr0
expression = 'c0/(c0+c1+c2+c3+c4)'
inputs = 'c0 c1 c2 c3 c4'
[]
[gr1]
type = ParsedCompute
buffer = gr1
expression = 'c1/(c0+c1+c2+c3+c4)'
inputs = 'c0 c1 c2 c3 c4'
[]
[gr2]
type = ParsedCompute
buffer = gr2
expression = 'c2/(c0+c1+c2+c3+c4)'
inputs = 'c0 c1 c2 c3 c4'
[]
[gr3]
type = ParsedCompute
buffer = gr3
expression = 'c3/(c0+c1+c2+c3+c4)'
inputs = 'c0 c1 c2 c3 c4'
[]
[gr4]
type = ParsedCompute
buffer = gr4
expression = '1 - gr0 - gr1 - gr2 - gr3'
inputs = 'gr0 gr1 gr2 gr3'
[]
[L_kappa_laplacian]
type = ReciprocalLaplacianFactor
buffer = 'L_kappa_laplacian'
factor = '${fparse ${L} * ${kappa}}'
[]
[smooth]
type = DeAliasingTensor
method = HOULI
buffer = smooth
[]
[]
[Solve]
[bnds]
type = ParsedCompute
buffer = 'bnds'
inputs = 'gr0 gr1 gr2 gr3 gr4'
expression = 'gr0^2 + gr1^2 + gr2^2 + gr3^2 + gr4^2'
[]
## FFT
[gr0_bar]
type = ForwardFFT
buffer = gr0_bar
input = gr0
[]
[gr1_bar]
type = ForwardFFT
buffer = gr1_bar
input = gr1
[]
[gr2_bar]
type = ForwardFFT
buffer = gr2_bar
input = gr2
[]
[gr3_bar]
type = ForwardFFT
buffer = gr3_bar
input = gr3
[]
[gr4_bar]
type = ForwardFFT
buffer = gr4_bar
input = gr4
[]
## bulk driving forces (constant mobility L, constant gamma = ${gamma})
[gr0_bulk_term]
type = ParsedCompute
buffer = 'gr0_bulk_term'
expression = '-${L}*${mu} * (gr0^3 - gr0 + 2*${gamma}*gr0*(gr1^2 + gr2^2 + gr3^2 + gr4^2))'
inputs = 'gr0 gr1 gr2 gr3 gr4'
[]
[gr1_bulk_term]
type = ParsedCompute
buffer = 'gr1_bulk_term'
expression = '-${L}*${mu} * (gr1^3 - gr1 + 2*${gamma}*gr1*(gr0^2 + gr2^2 + gr3^2 + gr4^2))'
inputs = 'gr0 gr1 gr2 gr3 gr4'
[]
[gr2_bulk_term]
type = ParsedCompute
buffer = 'gr2_bulk_term'
expression = '-${L}*${mu} * (gr2^3 - gr2 + 2*${gamma}*gr2*(gr0^2 + gr1^2 + gr3^2 + gr4^2))'
inputs = 'gr0 gr1 gr2 gr3 gr4'
[]
[gr3_bulk_term]
type = ParsedCompute
buffer = 'gr3_bulk_term'
expression = '-${L}*${mu} * (gr3^3 - gr3 + 2*${gamma}*gr3*(gr0^2 + gr1^2 + gr2^2 + gr4^2))'
inputs = 'gr0 gr1 gr2 gr3 gr4'
[]
[gr4_bulk_term]
type = ParsedCompute
buffer = 'gr4_bulk_term'
expression = '-${L}*${mu} * (gr4^3 - gr4 + 2*${gamma}*gr4*(gr0^2 + gr1^2 + gr2^2 + gr3^2))'
inputs = 'gr0 gr1 gr2 gr3 gr4'
[]
[NL_gr0]
type = ForwardFFT
buffer = NL_gr0
input = 'gr0_bulk_term'
[]
[NL_gr1]
type = ForwardFFT
buffer = NL_gr1
input = 'gr1_bulk_term'
[]
[NL_gr2]
type = ForwardFFT
buffer = NL_gr2
input = 'gr2_bulk_term'
[]
[NL_gr3]
type = ForwardFFT
buffer = NL_gr3
input = 'gr3_bulk_term'
[]
[NL_gr4]
type = ForwardFFT
buffer = NL_gr4
input = 'gr4_bulk_term'
[]
[NL_gr0_smooth]
type = ParsedCompute
buffer = NL_gr0_smooth
expression = 'smooth*NL_gr0'
inputs = 'smooth NL_gr0'
[]
[NL_gr1_smooth]
type = ParsedCompute
buffer = NL_gr1_smooth
expression = 'smooth*NL_gr1'
inputs = 'smooth NL_gr1'
[]
[NL_gr2_smooth]
type = ParsedCompute
buffer = NL_gr2_smooth
expression = 'smooth*NL_gr2'
inputs = 'smooth NL_gr2'
[]
[NL_gr3_smooth]
type = ParsedCompute
buffer = NL_gr3_smooth
expression = 'smooth*NL_gr3'
inputs = 'smooth NL_gr3'
[]
[NL_gr4_smooth]
type = ParsedCompute
buffer = NL_gr4_smooth
expression = 'smooth*NL_gr4'
inputs = 'smooth NL_gr4'
[]
# Properties for output
[grain_ID]
type = ParsedCompute
buffer = 'grain_ID'
expression = 'if(gr0>=max(gr1,max(gr2,max(gr3,gr4))), 0, if(gr1>=max(gr2,max(gr3,gr4)), 1, if(gr2>=max(gr3,gr4), 2, if(gr3>=gr4, 3, 4))))'
inputs = 'gr0 gr1 gr2 gr3 gr4'
[]
[grad_gr0]
type = GradientVector
buffer = grad_gr0
input = gr0
[]
[grad_gr1]
type = GradientVector
buffer = grad_gr1
input = gr1
[]
[grad_gr2]
type = GradientVector
buffer = grad_gr2
input = gr2
[]
[grad_gr3]
type = GradientVector
buffer = grad_gr3
input = gr3
[]
[grad_gr4]
type = GradientVector
buffer = grad_gr4
input = gr4
[]
[gradient_energy]
type = GradientEnergyDensity
buffer = gradient_energy
gradient_buffers = 'grad_gr0 grad_gr1 grad_gr2 grad_gr3 grad_gr4'
kappa = ${kappa}
[]
[total_energy]
type = ParsedCompute
buffer = 'total_energy'
expression = '${mu} * (${f0}) + gradient_energy'
inputs = 'gr0 gr1 gr2 gr3 gr4 gradient_energy'
[]
[]
[]
[Postprocessors]
[total_gb_energy]
type = TensorIntegralPostprocessor
buffer = total_energy
[]
[]
[TensorSolver]
type = AdamsBashforthMoulton
buffer = 'gr0 gr1 gr2 gr3 gr4'
reciprocal_buffer = 'gr0_bar gr1_bar gr2_bar gr3_bar gr4_bar'
linear_reciprocal = 'L_kappa_laplacian L_kappa_laplacian L_kappa_laplacian L_kappa_laplacian L_kappa_laplacian'
nonlinear_reciprocal = 'NL_gr0_smooth NL_gr1_smooth NL_gr2_smooth NL_gr3_smooth NL_gr4_smooth'
substeps = 1000
predictor_order = 1
corrector_order = 1
corrector_steps = 1
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'gr0 gr1 gr2 gr3 gr4 bnds grain_ID'
output_mode = 'NODE NODE NODE NODE NODE NODE CELL'
enable_hdf5 = true
transpose = true
[]
[]
[Executioner]
type = Transient
dt = 5
num_steps = 50
[]
[Outputs]
csv = true
execute_on = 'INITIAL TIMESTEP_END'
[]
(benchmarks/01_spinodal_decomposition/1a_implicit_fd.i)
#
# PFHub benchmark 1a (spinodal decomposition, square periodic domain) with the
# implicit finite-difference solver. Same problem as 1a.i:
#
# dc/dt = div(M grad(f'(c) - kappa Laplacian(c))), f = rho_s (c - c_alpha)^2 (c_beta - c)^2
#
# Backward Euler, with the biharmonic implicit and the double-well flux
# divergence explicit:
#
# (c - c^n)/dt + M kappa Laplacian(Laplacian(c)) = M Laplacian(f'(c^n))
#
# M = 5, kappa = 2, rho_s = 5, c_alpha = 0.3, c_beta = 0.7.
#
# substeps per time step, as in 1a_solver.i (each substep is one linear solve)
substeps = 1000
[Domain]
dim = 2
nx = 200
ny = 200
xmax = 200
ymax = 200
parallel_mode = REAL_SPACE
periodic_directions = 'X Y'
[]
[TensorBuffers]
[c]
[]
[dcdt]
[]
[bih_c]
[]
[mu]
[]
[rhs]
[]
# postprocessing
[lap_c]
[]
[F]
[]
[]
[TensorComputes]
[Initialize]
[c]
type = ParsedCompute
buffer = c
extra_symbols = true
expression = 'c0+epsilon*(cos(0.105*x)*cos(0.11*y)+(cos(0.13*x)*cos(0.087*y))^2+cos(0.025*x-0.15*y)*cos(0.07*x-0.02*y))'
constant_names = 'c0 epsilon'
constant_expressions = '0.5 0.01'
[]
[]
[Solve]
# implicit operators
[dcdt]
type = TimeDerivativeOperator
buffer = dcdt
input = c
[]
[bih]
type = FiniteDifferenceBiharmonic
buffer = bih_c
input = c
factor = 10 # M * kappa
[]
# explicit source M Laplacian(f'(c^n))
[explicit]
[mu]
type = ParsedCompute
buffer = mu
expression = 'rho_s*(c-c_alpha)^2*(c_beta-c)^2'
constant_names = 'rho_s c_alpha c_beta'
constant_expressions = '5 0.3 0.7'
derivatives = c
inputs = c
[]
[rhs]
type = FiniteDifferenceLaplacian
buffer = rhs
input = mu
factor = 5 # M
[]
[]
[]
[Postprocess]
# On the periodic grid sum |grad c|^2 = -sum c Laplacian(c) for the 3-point
# Laplacian, so the gradient energy (kappa/2)|grad c|^2 integrates to the same
# total as -(kappa/2) c Laplacian(c).
[lap_c]
type = FiniteDifferenceLaplacian
buffer = lap_c
input = c
[]
[F]
type = ParsedCompute
buffer = F
expression = 'rho_s * (c-c_alpha)^2 * (c_beta-c)^2 - kappa / 2 * c * lap_c'
constant_names = 'rho_s c_alpha c_beta kappa'
constant_expressions = '5 0.3 0.7 2'
inputs = 'c lap_c'
[]
[]
[]
[TensorSolver]
type = ImplicitFDSolver
buffer = c
operators = 'dcdt bih'
source = rhs
root_compute = explicit
substeps = ${substeps}
[]
[UserObjects]
[terminator]
type = Terminator
expression = change<1e-4
[]
[]
[Postprocessors]
[min_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MIN
execute_on = 'TIMESTEP_END'
[]
[max_c]
type = TensorExtremeValuePostprocessor
buffer = c
value_type = MAX
execute_on = 'TIMESTEP_END'
[]
[F]
type = TensorIntegralPostprocessor
buffer = F
[]
[change]
type = TensorIntegralChangePostprocessor
buffer = c
[]
# cumulative wall time of the run (s)
[wall_time]
type = PerfGraphData
section_name = Root
data_type = TOTAL
[]
[]
[Problem]
type = TensorProblem
[]
[Executioner]
type = Transient
num_steps = 1000
[TimeStepper]
type = IterationAdaptiveDT
growth_factor = 1.1
dt = 1
[]
dtmax = 300
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'c mu'
enable_hdf5 = true
[]
[]
[Outputs]
csv = true
perf_graph = true
execute_on = 'TIMESTEP_END'
[]
(test/tests/real_space/laplacian3d.i)
[Domain]
dim = 3
nx = 30
ny = 24
nz = 17
parallel_mode = REAL_SPACE
periodic_directions = 'X Y Z'
xmin = -1
ymin = -1
zmin = -1
xmax = 1
ymax = 1
zmax = 1
[]
[TensorBuffers]
[u]
type = PlainTensorBuffer
[]
[u_gold]
type = PlainTensorBuffer
[]
[lap]
type = PlainTensorBuffer
[]
[]
[TensorComputes]
[Initialize]
[u_ic]
type = ParsedCompute
buffer = u
expression = 'sin(pi*(x+0.5))*sin(pi*(y+0.5))*sin(pi*(z+0.5))'
extra_symbols = true
[]
[lap_gold_ic]
type = ParsedCompute
buffer = lap_gold
#expression = '0'
expression = '-3*pi*pi*sin(pi*(x+0.5))*sin(pi*(y+0.5))*sin(pi*(z+0.5))'
extra_symbols = true
expand = REAL
[]
[]
[Solve]
[lap_compute]
type = FiniteDifferenceLaplacian
input = u
buffer = lap
[]
[diff]
type = ParsedCompute
buffer = diff
expression = 'abs(lap - lap_gold)'
inputs = 'lap lap_gold'
[]
[]
[]
[Postprocessors]
[max_error]
type = TensorExtremeValuePostprocessor
buffer = diff
value_type = MAX
[]
[]
[TensorOutputs]
# active = ''
[lap]
type = XDMFTensorOutput
buffer = 'lap lap_gold diff'
enable_hdf5 = true
[]
[]
[Executioner]
type = Transient
num_steps = 1
[]
[Outputs]
[out]
type = CSV
[]
[]
(test/tests/tensor_compute/smooth_rectangle.i)
[Domain]
dim = 2
nx = 100
ny = 100
xmax = 20
ymax = 20
mesh_mode = DUMMY
device_names = cpu
[]
[TensorComputes]
[Initialize]
[rectangle_sharp]
type = SmoothRectangleCompute
buffer = rectangle_sharp
x1 = 5
x2 = 15
y1 = 5
y2 = 15
inside = -1
outside = 3
[]
[rectangle_cos]
type = SmoothRectangleCompute
buffer = rectangle_cos
x1 = 5
x2 = 15
y1 = 5
y2 = 15
inside = -1
outside = 3
profile = COS
int_width = 1
[]
[rectangle_tanh]
type = SmoothRectangleCompute
buffer = rectangle_tanh
x1 = 5
x2 = 15
y1 = 5
y2 = 15
inside = -1
outside = 3
profile = TANH
int_width = 1
[]
[]
[]
[Problem]
type = TensorProblem
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'rectangle_sharp rectangle_cos rectangle_tanh'
enable_hdf5 = true
[]
[]
[Executioner]
type = Transient
num_steps = 0
[]
[Outputs]
perf_graph = 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
[]
[]
[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'
[]
[Problem]
type = LatticeBoltzmannProblem
substeps = 200
print_debug_output = true
scalar_constant_names = 'tau_h D sigma'
scalar_constant_values = '${tau_h} ${D} ${sigma}'
residual_tensor = speed
[]
[Executioner]
type = Transient
num_steps = 50
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'phi velocity rho'
output_mode = 'Cell Cell Cell'
enable_hdf5 = true
transpose = false
[]
[]
(test/tests/tensor_compute/pfc_fcc.i)
# FCC PFC test designed to clearly show atomic-scale density fluctuations
# Using larger epsilon for stronger density modulations
# Based on PhysRevE.81.061601 two-mode FCC PFC model
# Grid parameters - fine resolution to resolve atomic peaks
N = 1 # Grid points per dimension
# PFC model parameters - LARGER epsilon for clearer atomic structure
psi_mean = -0.2 # Mean density
epsilon = 0.15 # Larger undercooling for stronger modulations
R1_param = 0.0 # Two-mode coupling (0 for max FCC stability)
Q1_param = ${fparse sqrt(4.0/3.0)} # FCC wave number ratio
# Domain size: want several unit cells but fine enough resolution
# The (111) wavelength is 2*pi/q0 = 2*pi*sqrt(3/4) ≈ 5.44 (dimensionless)
# Let's fit about 'N' wavelengths in each direction
Lx = ${fparse N * 2 * pi / ${Q1_param}}
Ly = ${Lx}
[Domain]
dim = 2
nx = ${fparse ${N} * 8 }
ny = ${fparse ${N} * 8 }
xmax = ${Lx}
ymax = ${Ly}
[]
[TensorComputes]
[Initialize]
# Random perturbations to seed crystal growth
[psi]
type = RandomTensor
buffer = psi
max = ${fparse ${psi_mean} + 0.01}
min = ${fparse ${psi_mean} - 0.01}
seed = 12345
[]
# Linear operator for FCC
[linear]
type = FCCPFCLinear
buffer = 'linear'
eps = ${epsilon}
R1 = ${R1_param}
Q1 = ${Q1_param}
mobility = 1.0
[]
# Dealiasing for cubic nonlinearity
[smooth_operator]
type = DeAliasingTensor
buffer = smooth_operator
method = HOULI
[]
[]
[Solve]
# Nonlinear term
[nl_div_psi_cubed]
type = FCCPFCNonlinear
buffer = NL
psi = psi
dealiasing = smooth_operator
mobility = 1.0
[]
# FFT for spectral solver
[psi_hat]
type = ForwardFFT
buffer = psi_hat
input = psi
[]
[]
[]
[TensorSolver]
type = AdamsBashforthMoulton
buffer = 'psi'
linear_reciprocal = 'linear'
nonlinear_reciprocal = 'NL'
reciprocal_buffer = 'psi_hat'
corrector_order = 1
corrector_steps = 3
predictor_order = 1
substeps = 1000
[]
[Postprocessors]
[max_psi]
type = TensorExtremeValuePostprocessor
buffer = psi
value_type = MAX
execute_on = 'initial timestep_end'
[]
[min_psi]
type = TensorExtremeValuePostprocessor
buffer = psi
value_type = MIN
execute_on = 'initial timestep_end'
[]
[mean_psi]
type = TensorAveragePostprocessor
buffer = psi
execute_on = 'initial timestep_end'
[]
[]
[TensorOutputs]
[xdmf]
type = XDMFTensorOutput
buffer = 'psi'
enable_hdf5 = true
[]
[]
[Executioner]
type = Transient
num_steps = 100
dt = 0.2 # Larger timestep for faster evolution
[]
[Outputs]
perf_graph = true
csv = true
[]