Reference
XCALibre.Mesh.Boundary
— Typestruct Boundary{S<:Symbol, UR<:UnitRange{<:Integer}}
name::S # Boundary patch name
IDs_range::UR # range to access boundary info (faces and boundary_cellsID)
end
XCALibre.Mesh.Cell
— Typestruct Cell{F<:AbstractFloat, SV3<:SVector{3,F},UR<:UnitRange{<:Integer}}
centre::SV3 # coordinate of cell centroid
volume::F # cell volume
nodes_range::UR # range to access cell nodes in Mesh3.cell_nodes
faces_range::UR # range to access cell faces info (faces, neighbours cells, etc.)
end
XCALibre.Mesh.Face3D
— Typestruct Face3D{
F<:AbstractFloat,
SV2<:SVector{2,<:Integer},
SV3<:SVector{3,F},
UR<:UnitRange{<:Integer}
}
nodes_range::UR # range to access face nodes in Mesh3.face_nodes
ownerCells::SV2 # IDs of face owner cells (always 2)
centre::SV3 # coordinates of face centre
normal::SV3 # face normal unit vector
e::SV3 # unit vector in the direction between owner cells
area::F # face area
delta::F # distance between owner cells centres
weight::F # linear interpolation weight
end
XCALibre.Mesh.Mesh3
— Typestruct Mesh3{VC, VI, VF<:AbstractArray{<:Face3D}, VB, VN, SV3, UR} <: AbstractMesh
cells::VC # vector of cells
cell_nodes::VI # vector of indices to access cell nodes
cell_faces::VI # vector of indices to access cell faces
cell_neighbours::VI # vector of indices to access cell neighbours
cell_nsign::VI # vector of indices to with face normal correction (1 or -1 )
faces::VF # vector of faces
face_nodes::VI # vector of indices to access face nodes
boundaries::VB # vector of boundaries
nodes::VN # vector of nodes
node_cells::VI # vector of indices to access node cells
get_float::SV3 # store mesh float type
get_int::UR # store mesh integer type
boundary_cellsID::VI # vector of indices of boundary cell IDs
end
XCALibre.Mesh.Node
— Typestruct Node{SV3<:SVector{3,<:AbstractFloat}, UR<:UnitRange{<:Integer}}
coords::SV3 # node coordinates
cells_range::UR # range to access neighbour cells in Mesh3.node_cells
end
XCALibre.UNV2.UNV2D_mesh
— MethodUNV2D_mesh(meshFile; scale=1, integer_type=Int64, float_type=Float64)
Read and convert 2D UNV mesh file into XCALibre.jl
Input
meshFile
– path to mesh file.
Optional arguments
scale
– used to scale mesh file e.g. scale=0.001 will convert mesh from mm to metres defaults to 1 i.e. no scalinginteger_type
- select interger type to use in the mesh (Int32 may be useful on GPU runs)float_type
- select interger type to use in the mesh (Float32 may be useful on GPU runs)
XCALibre.UNV3.UNV3D_mesh
— MethodUNV3D_mesh(unv_mesh; scale=1, integer_type=Int64, float_type=Float64)
Read and convert 3D UNV mesh file into XCALibre.jl. Note that a limitation of the .unv mesh format is that it only supports the following 3D cells:
- Tetahedrals
- Prisms
- Hexahedrals
Input
unv_mesh
– path to mesh file.
Optional arguments
scale
– used to scale mesh file e.g. scale=0.001 will convert mesh from mm to metres defaults to 1 i.e. no scalinginteger_type
- select interger type to use in the mesh (Int32 may be useful on GPU runs)float_type
- select interger type to use in the mesh (Float32 may be useful on GPU runs)
XCALibre.FoamMesh.FOAM3D_mesh
— MethodFOAM3D_mesh(mesh_file; scale=1, integer_type=Int64, float_type=Float64)
Read and convert 3D OpenFOAM mesh file into XCALibre.jl. Note that, at present, it is not recommended to run 2D cases using meshes imported using this function.
Input
mesh_file
– path to mesh file.
Optional arguments
scale
– used to scale mesh file e.g. scale=0.001 will convert mesh from mm to metres defaults to 1 i.e. no scalinginteger_type
- select interger type to use in the mesh (Int32 may be useful on GPU runs)float_type
- select interger type to use in the mesh (Float32 may be useful on GPU runs)
XCALibre.Multithread.activate_multithread
— Methodactivate_multithread(backend::CPU; nthreads=1) = BLAS.set_num_threads(nthreads)
Convenience function to set number of BLAS threads.
Input arguments
backend
is the only required input which must beCPU()
fromKernelAbstractions.jl
nthreads
can be used to set the number of BLAS cores (defaultnthreads=1
)
XCALibre.Fields.ScalarField
— Typestruct ScalarField{VF,M<:AbstractMesh,BC} <: AbstractScalarField
values::VF # scalar values at cell centre
mesh::M # reference to mesh
BCs::BC # store user-provided boundary conditions
end
XCALibre.Fields.VectorField
— Typestruct VectorField{S1<:ScalarField,S2,S3,M<:AbstractMesh,BC} <: AbstractVectorField
x::S1 # x-component is itself a `ScalarField`
y::S2 # y-component is itself a `ScalarField`
z::S3 # z-component is itself a `ScalarField`
mesh::M
BCs::BC
end
XCALibre.Fields.initialise!
— Methodfunction initialise!(field, value) # dummy function for documentation
# Assign `value` to field in-place
nothing
end
This function will set the given field
to the value
provided in-place. Useful for initialising fields prior to running a simulation.
Input arguments
field
specifies the field to be initialised. The field must be either aAbractScalarField
orAbstractVectorField
value
defines the value to be set. This should be a scalar or vector (3 components) depending on the field to be modified e.g. for anAbstractVectorField
we can specify asvalue=[10,0,0]
Note: in most cases the fields to be modified are stored within a physics model i.e. a Physics
object. Thus, the argument value
must fully qualify the model. For example, if we have created a Physics
model named mymodel
to set the velocity field, U
, we would set the argument field
to mymodel.momentum.U
. See the example below.
Example
initialise!(mymodel.momentum.U, [2.5, 0, 0])
initialise!(mymodel.momentum.p, 1.25)
XCALibre.Discretise.Dirichlet
— TypeDirichlet <: AbstractDirichlet
Dirichlet boundary condition model.
Inputs
ID
Name of the boundary given as a symbol (e.g. :inlet). Internally it gets replaced with the boundary index IDvalue
Scalar or Vector value for Dirichlet boundary condition
XCALibre.Discretise.DirichletFunction
— TypeDirichletFunction(ID, value) <: AbstractDirichlet
Dirichlet boundary condition defined with user-provided function.
Input
ID
Name of the boundary given as a symbol (e.g. :inlet). Internally it gets replaced with the boundary index IDvalue
Custom function or struct (<:XCALibreUserFunctor) for Dirichlet boundary conditionIDs_range
Range of indices to access boundary patch faces
Function requirements
The function passed to this boundary condition must have the following signature:
f(coords, time, index) = SVector{3}(ux, uy, uz)
Where, coords
is a vector containing the coordinates of a face, time
is the current time in transient simulations (and the iteration number in steady simulations), and index
is the local face index (from 1 to N
, where N
is the number of faces in a given boundary). The function must return an SVector (from StaticArrays.jl) representing the velocity vector.
XCALibre.Discretise.Empty
— TypeEmpty <: AbstractPhysicalConstraint
Empty boundary condition model (currently only configured for zero gradient)
Fields
- 'ID' – Boundary ID
value
– Scalar or Vector value for Empty boundary condition.
XCALibre.Discretise.Extrapolated
— TypeExtrapolated <: AbstractNeumann
This boundary condition extrapolates the face value using the interior cell value. Equivalent to setting a zero gradient boundary condition (semi-implicitly). It applies to both scalar and vector fields.
Example
Extrapolated(:outlet)
XCALibre.Discretise.FixedTemperature
— TypeFixedTemperature(name::Symbol, model::Enthalpy; T::Number)
Fixed temperature boundary condition model, which allows the user to specify wall temperature that can be translated to the energy specific model, such as sensible enthalpy.
Inputs
name
: Name of the boundary given as a symbol (e.g. :inlet). Internally it gets replaced with the boundary index IDT
: keyword argument use to define the boundary temperaturemodel
: defines the underlyingEnergy
model to be used (currently onlyEnthalpy
is available)
Example
FixedTemperature(:inlet, T=300.0, Enthalpy(cp=cp, Tref=288.15))
XCALibre.Discretise.Neumann
— TypeNeumann <: AbstractNeumann
Neumann boundary condition model to set the gradient at the boundary explicitly (currently only configured for zero gradient)
Inputs
ID
Name of the boundary given as a symbol (e.g. :inlet). Internally it gets replaced with the boundary index IDvalue
Scalar providing face normal gradient
Example
Neumann(:outlet, 0)
XCALibre.Discretise.NeumannFunction
— TypeNeumannFunction(ID, value) <: AbstractNeumann
Neumann boundary condition defined with user-provided function.
Input
ID
Name of the boundary given as a symbol (e.g. :inlet). Internally it gets replaced with the boundary index IDvalue
Custom function defining the desired Neumann boundary condition.
Function requirements
The function passed to this boundary condition has not yet been implemented. However, users can pass a custom struct to specialise the internal implementations of many functions. By default, at present, this function will assign a zero gradient boundary condition.
XCALibre.Discretise.Periodic
— Typestruct Periodic{I,V,R<:UnitRange} <: AbstractPhysicalConstraint
ID::I
value::V
end
Periodic boundary condition model.
Fields
ID
is the name of the boundary given as a symbol (e.g. :inlet). Internally it gets replaced with the boundary index IDvalue
tuple containing information needed to apply this boundary
XCALibre.Discretise.Symmetry
— TypeSymmetry <: AbstractBoundary
Symmetry boundary condition vector and scalar fields. Notice that for scalar fields, this boundary condition applies an explicit zero gradient condition. In some rare cases, the use of an Extrapolated
condition for scalars may be beneficial (to assign a semi-implicit zero gradient condition)
Input
ID
Name of the boundary given as a symbol (e.g. :freestream). Internally it gets replaced with the boundary index ID
Example
Symmetry(:freestream)
XCALibre.Discretise.Wall
— TypeWall <: AbstractDirichlet
Wall boundary condition model for no-slip or moving walls (linear motion). It should be applied to the velocity vector, and in most cases, its scalar variant should be applied to scalars.
Inputs
ID
represents the name of the boundary given as a symbol (e.g. :inlet). Internally it gets replaced with the boundary index IDvalue
should be given as a vector for the velocity e.g. [10,0,0]. For scalar fields such as the pressure the value entry can be omitted or set to zero explicitly.
Examples
Wall(:plate, [0, 0, 0]) # no-slip wall condition for velocity
Wall(:plate) # corresponding definition for scalars, e.g. pressure
Wall(:plate, 0) # alternative definition for scalars, e.g. pressure
XCALibre.Discretise.Zerogradient
— TypeZerogradient <: AbstractNeumann
Zerogradient boundary condition model (explicitly applied to the boundary)
Input
ID
Name of the boundary given as a symbol (e.g. :inlet). Internally it gets replaced with the boundary index ID
Example
Zerogradient(:inlet)
XCALibre.Discretise.construct_periodic
— Methodconstruct_periodic(mesh, backend, patch1::Symbol, patch2::Symbol)
Function for construction of periodic boundary conditions.
Input
mesh
– Mesh.backend
– Backend configuraton.patch1
– Primary periodic patch ID.patch2
– Neighbour periodic patch ID.
Output
periodic::Tuple - tuple containing boundary defintions for
patch1
andpatch2
i.e. (periodic1, periodic2). The fields ofperiodic1
andperiodic2
areID
– Index to access boundary information in mesh objectvalue
represents aNamedTuple
with the following keyword arguments:- index – ID used to find boundary geometry information in the mesh object
- distance – perpendicular distance between the patches
- face_map – vector providing indeces to faces of match patch
- ismaster – flat to identify one of the patch pairs as the main patch
Example
- `periodic = construct_periodic(mesh, CPU(), :top, :bottom)` - Example using CPU
backend with periodic boundaries named `top` and `bottom`.
XCALibre.Discretise.@define_boundary
— Macromacro define_boundary(boundary, operator, definition)
quote
@inline (bc::$boundary)(
term::Operator{F,P,I,$operator}, cellID, zcellID, cell, face, fID, i, component, time
) where {F,P,I} = $definition
end |> esc
end
Macro to reduce boilerplate code when defining boundary conditions (implemented as functors) and provides access to key fields needed in the implementation of boundary conditions, such as the boundary cell and face objects (more details below)
Input arguments
boundary
specifies the boundary type being definedoperator
specifies the operator to which the condition applies e.g.Laplacian
definition
provides the implementation details
Available fields
term
reference to operator on which the boundary applies (gives access to the field and mesh)cellID
ID of the corresponding boundary cellzcellID
sparse matrix linear index for the cellcell
gives access to boundary cell object and corresponding informationface
gives access to boundary face object and corresponding informationfID
ID of the boundary face (to indexMesh2.faces
vector)i
local index of the boundary faces within a kernel or loopcomponent
for vectors this specifies the components being evaluated (access ascomponent.value
). For scalarscomponent = nothing
time
provides the current simulation time. This only applies to time dependent boundary implementation defined as functions or neural networks.
Example
Below the use of this macro is illustrated for the implementation of a Dirichlet
boundary condition acting on the Laplacian
using the Linear
scheme:
@define_boundary Dirichlet Laplacian{Linear} begin
J = term.flux[fID] # extract operator flux
(; area, delta) = face # extract boundary face information
flux = J*area/delta # calculate the face flux
ap = term.sign*(-flux) # diagonal (cell) matrix coefficient
ap, ap*bc.value # return `ap` and `an`
end
When called, this functor will return two values ap
and an
, where ap
is the cell contribution for approximating the boundary face value, and an
is the explicit part of the face value approximation i.e. ap
contributes to the diagonal of the sparse matrix (left-hand side) and an
is the explicit contribution assigned to the solution vector b
on the right-hand of the linear system of equations $Ax = b$
XCALibre.Solve.JacobiSmoother
— Typestruct JacobiSmoother{L,F,V} <: AbstractSmoother
loops::L
omega::F
x_temp::V
end
Structure to hold information for using the weighted Jacobi smoother.
Fields
loops
represents the number of smoothing iterations.omega
represents the relaxation weight, 1 corresponds to no weighting. Typically a weight of 2/3 is used to filter high frequencies in the residual field.x_temp
is a vector used internally to store intermediate solutions.
Constructors
JacobiSmoother(mesh::AbstractMesh)
This constructor takes a mesh object as input to create the smoother. The number of cells in the mesh is used to allocate x_temp
. The number of loops and smoother factor are set to 5 and 1, respectively.
JacobiSmoother(; domain, loops, omega=2/3)
Alternative constructor for JacobiSmoother
.
keyword arguments
domain
represents a mesh object of typeAbstractMesh
.loops
is the number of iterations to be usedomega
represents the weighting factor, 1 does not relax the system, 2/3 is found to work well for smoothing high frequencies in the residual field
XCALibre.Solve.Runtime
— MethodRuntime(;
# keyword arguments
iterations::I,
write_interval::I,
time_step::N
) where {I<:Integer,N<:Number} = begin
# returned Runtime struct
Runtime{I<:Integer,F<:AbstractFloat}
(
iterations=iterations,
dt=time_step,
write_interval=write_interval
)
end
This is a convenience function to set the top-level runtime information. The inputs are all keyword arguments and provide basic information to flow solvers just before running a simulation.
Input arguments
iterations::Integer
: specifies the number of iterations in a simulation run.write_interval::Integer
: defines how often simulation results are written to file (on the current working directory). The interval is currently based on number of iterations. Set to-1
to run without writing results to file.time_step::AbstractFloat
: the time step to use in the simulation. Notice that for steady solvers this is simply a counter and it is recommended to simply use1
.
Example
runtime = Runtime(iterations=2000, time_step=1, write_interval=2000)
XCALibre.Solve.Schemes
— TypeSchemes(;
# keyword arguments and their default values
time=SteadyState,
divergence=Linear,
laplacian=Linear,
gradient=Gauss,
limiter=nothing) = begin
# Returns Schemes struct used to configure discretisation
Schemes(
time=time,
divergence=divergence,
laplacian=laplacian,
gradient=gradient,
limiter=limiter
)
end
The Schemes
struct is used at the top-level API to help users define discretisation schemes for every field solved.
Inputs
time
: is used to set the time schemes (default isSteadyState
)divergence
: is used to set the divergence scheme (default isLinear
)laplacian
: is used to set the laplacian scheme (default isLinear
)gradient
: is used to set the gradient scheme (default isGauss
)limiter
: is used to specify if gradient limiters should be used, currently supported limiters includeFaceBased
andMFaceBased
(default isnothing
)
XCALibre.Solve.SolverSetup
— MethodSolverSetup(;
# required keyword arguments
solver::S,
preconditioner::PT,
convergence,
relax,
# optional keyword arguments
float_type=Float64,
smoother=nothing,
limit=nothing,
itmax::Integer=1000,
atol=(eps(_get_float(region)))^0.9,
rtol=_get_float(region)(1e-1)
) where {S,PT<:PreconditionerType} = begin
return SolverSetup(kwargs...)
end
This function is used to provide solver settings that will be used internally in XCALibre.jl. It returns a SolverSetup
object with solver settings that are used internally by the flow solvers.
Input arguments
solver
: solver object from Krylov.jl and it could be one ofBicgstab()
,Cg()
,Gmres()
which are re-exported in XCALibre.jlpreconditioner
: instance of preconditioner to be used e.g. Jacobi()convergence
sets the stopping criteria of this fieldrelax
: specifies the relaxation factor to be used e.g. set to 1 for no relaxationsmoother
: specifies smoothing method to be applied before discretisation.JacobiSmoother
: is currently the only choice (defaults tonothing
)limit
: used in some solvers to bound the solution within these limits e.g. (min, max). It defaults tonothing
itmax
: maximum number of iterations in a single solver pass (defaults to 1000)atol
: absolute tolerance for the solver (default to eps(FloatType)^0.9)rtol
: set relative tolerance for the solver (defaults to 1e-1)float_type
: specifies the floating point type to be used by the solver. It is also used to estimate the absolute tolerance for the solver (defaults toFloat64
)
XCALibre.ModelPhysics.Compressible
— TypeCompressible <: AbstractCompressible
Compressible fluid model containing fluid field parameters for compressible flows with constant parameters - ideal gas with constant viscosity.
Fields
- 'nu' – Fluid kinematic viscosity.
- 'cp' – Fluid specific heat capacity.
gamma
– Ratio of specific heats.Pr
– Fluid Prantl number.
Examples
Fluid{Compressible}(; nu=1E-5, cp=1005.0, gamma=1.4, Pr=0.7)
- Constructur with default values.
XCALibre.ModelPhysics.Fluid
— TypeFluid <: AbstractFluid
Abstract fluid model type for constructing new fluid models.
Fields
- 'args' – Model arguments.
XCALibre.ModelPhysics.Incompressible
— TypeIncompressible <: AbstractIncompressible
Incompressible fluid model containing fluid field parameters for incompressible flows.
Fields
- 'nu' – Fluid kinematic viscosity.
- 'rho' – Fluid density.
Examples
Fluid{Incompressible}(nu=0.001, rho=1.0)
- Constructor with default values.
XCALibre.ModelPhysics.KEquation
— TypeKEquation <: AbstractTurbulenceModel
KEquation LES model containing all Smagorinksy field parameters.
Fields
nut
– Eddy viscosity ScalarField.nutf
– Eddy viscosity FaceScalarField.coeffs
– Model coefficients.
XCALibre.ModelPhysics.KOmega
— TypeKOmega <: AbstractTurbulenceModel
kOmega model containing all kOmega field parameters.
Fields
k
– Turbulent kinetic energy ScalarField.omega
– Specific dissipation rate ScalarField.nut
– Eddy viscosity ScalarField.kf
– Turbulent kinetic energy FaceScalarField.omegaf
– Specific dissipation rate FaceScalarField.nutf
– Eddy viscosity FaceScalarField.coeffs
– Model coefficients.
XCALibre.ModelPhysics.KOmegaLKE
— TypeKOmegaLKE <: AbstractTurbulenceModel
kOmega model containing all kOmega field parameters.
Fields
k
– Turbulent kinetic energy ScalarField.omega
– Specific dissipation rate ScalarField.kl
– ScalarField.nut
– Eddy viscosity ScalarField.kf
– Turbulent kinetic energy FaceScalarField.omegaf
– Specific dissipation rate FaceScalarField.klf
– FaceScalarField.nutf
– Eddy viscosity FaceScalarField.coeffs
– Model coefficients.Tu
– Freestream turbulence intensity for model.y
– Near-wall distance for model.
XCALibre.ModelPhysics.LES
— TypeLES <: AbstractLESModel
Abstract LES model type for constructing LES models.
Fields
- 'args' – Model arguments.
XCALibre.ModelPhysics.Laminar
— TypeLaminar <: AbstractTurbulenceModel
Laminar model definition for physics API.
XCALibre.ModelPhysics.Momentum
— Typestruct Momentum{V,S,SS} <: AbstractMomentumModel
U::V
p::S
sources::SS
end
Momentum model containing key momentum fields.
Fields
- 'U' – Velocity VectorField.
- 'p' – Pressure ScalarField.
- 'sources' – Momentum model sources.
Examples
- `Momentum(mesh::AbstractMesh)
XCALibre.ModelPhysics.Physics
— Typestruct Physics{T,F,M,Tu,E,D,BI}
time::T
fluid::F
momentum::M
turbulence::Tu
energy::E
domain::D
boundary_info::BI
end
XCALibre's parametric Physics type for user-level API. Also used to dispatch flow solvers.
Fields
- 'time::Union{Steady, Transient}' – User-provided time model.
- 'fluid::AbstractFluid' – User-provided
Fluid
` model. - 'momentum' – Momentum model. Currently this is auto-generated by the
Physics
constructor - 'turbulence::AbstractTurbulenceModel' – User-provided
Turbulence
` model. - 'energy:AbstractEnergyModel' – User-provided
Energy
model. - 'domain::AbstractMesh ' – User-provided
Mesh
. Must be adapted to target device before constructing a Physics object. - 'boundaryinfo::boundaryinfo' – Mesh boundary information. Auto-generated by the
Physics
constructor
XCALibre.ModelPhysics.Physics
— MethodPhysics(; time, fluid, turbulence, energy, domain)::Physics{T,F,M,Tu,E,D,BI}
Physics
constructor part of the top-level API. It can be used to define the Physics and models relevant to a simulation. This constructor uses keyword arguments to allow users to fine-tune their simulations, whilst some fields are auto-generated behind the scenes for convenience (Momentum
and boundary_info
). Where:
time
- specified the time model (Steady or Transient)fluid
- specifies the type of fluid (Incompressible, etc.)turbulence
- specified the Turbulence modelenergy
- specifies the Energy treatmentdomain
- provides the mesh to used (must be adapted to the target backend device)
XCALibre.ModelPhysics.RANS
— TypeRANS <: AbstractRANSModel
Abstract RANS model type for consturcting RANS models.
Fields
- 'args' – Model arguments.
XCALibre.ModelPhysics.Smagorinsky
— TypeSmagorinsky <: AbstractTurbulenceModel
Smagorinsky LES model containing all Smagorinksy field parameters.
Fields
nut
– Eddy viscosity ScalarField.nutf
– Eddy viscosity FaceScalarField.coeffs
– Model coefficients.
XCALibre.ModelPhysics.Steady
— TypeSteady
Steady model for Physics model API.
Examples
Steady()
XCALibre.ModelPhysics.Transient
— TypeTransient
Transient model for Physics model API.
Examples
Transient()
XCALibre.ModelPhysics.WeaklyCompressible
— TypeWeaklyCompressible <: AbstractCompressible
Weakly compressible fluid model containing fluid field parameters for weakly compressible flows with constant parameters - ideal gas with constant viscosity.
Fields
- 'nu' – Fluid kinematic viscosity.
- 'cp' – Fluid specific heat capacity.
gamma
– Ratio of specific heats.Pr
– Fluid Prandtl number.
Examples
Fluid{WeaklyCompressible}(; nu=1E-5, cp=1005.0, gamma=1.4, Pr=0.7)
- Constructor with
default values.
XCALibre.ModelPhysics.Ttoh!
— MethodTtoh!(model::Physics{T1,F,M,Tu,E,D,BI}, T::ScalarField, h::ScalarField
) where {T1,F<:AbstractCompressible,M,Tu,E,D,BI}
Function coverts temperature ScalarField to sensible enthalpy ScalarField.
Input
model
– Physics model defined by user.T
– Temperature ScalarField.h
– Sensible enthalpy ScalarField.
XCALibre.ModelPhysics.change
— Methodchange(model::Physics, property, value) => updatedModel::Physics
A convenience function to change properties of an exisitng Physics
model.
Input arguments
model::Physics
aPhysics
model to modifyproperty
is a symbol specifying the property to changevalue
is the new setting for the specifiedproperty
Output
This function return a new Physics
object
Example
To change a model to run a transient simulation e.g. after converging in steady state
modelTransient = change(model, :time, Transient())
XCALibre.ModelPhysics.energy!
— Methodenergy::SensibleEnthalpyModel, model::Physics{T1,F,M,Tu,E,D,BI}, prev, mdotf, rho, mueff, time, config
) where {T1,F,M,Tu,E,D,BI,E1}
Run energy transport equations.
Input
energy
: Energy model.model
: Physics model defined by user.prev
: Previous energy cell values.mdtof
: Face mass flow.rho
: Density ScalarField.mueff
: Effective viscosity FaceScalarField.time
: Simulation runtime.config
: Configuration structure defined by user with solvers, schemes, runtime and hardware structures set.
XCALibre.ModelPhysics.htoT!
— MethodhtoT!(model::Physics{T1,F,M,Tu,E,D,BI}, h::ScalarField, T::ScalarField
) where {T1,F<:AbstractCompressible,M,Tu,E,D,BI}
Function coverts sensible enthalpy ScalarField to temperature ScalarField.
Input
model
– Physics model defined by user.h
– Sensible enthalpy ScalarField.T
– Temperature ScalarField.
XCALibre.ModelPhysics.initialise
— Methodinitialise(energy::SensibleEnthalpy, model::Physics{T1,F,M,Tu,E,D,BI}, mdotf, rho, peqn, config
) where {T1,F,M,Tu,E,D,BI})
Initialisation of energy transport equations.
Input
energy
: Energy model.model
: Physics model defined by user.mdtof
: Face mass flow.rho
: Density ScalarField.peqn
: Pressure equation.config
: Configuration structure defined by user with solvers, schemes, runtime and hardware structures set.
Output
SensibleEnthalpyModel
: Energy model struct containing energy equation.
XCALibre.ModelPhysics.initialise
— Methodinitialise(turbulence::KEquation, model::Physics{T,F,M,Tu,E,D,BI}, mdotf, peqn, config
) where {T,F,M,Tu,E,D,BI}
Initialisation of turbulent transport equations.
Input
turbulence
– turbulence model.model
– Physics model defined by user.mdtof
– Face mass flow.peqn
– Pressure equation.config
– Configuration structure defined by user with solvers, schemes, runtime and hardware structures set.
Output
KEquationModel( turbulence, Δ, magS, ModelState((), false) )
– Turbulence model structure.
XCALibre.ModelPhysics.initialise
— Methodinitialise(turbulence::KOmega, model::Physics{T,F,M,Tu,E,D,BI}, mdotf, peqn, config
) where {T,F,M,Tu,E,D,BI}
Initialisation of turbulent transport equations.
Input
turbulence
– turbulence model.model
– Physics model defined by user.mdtof
– Face mass flow.peqn
– Pressure equation.config
– Configuration structure defined by user with solvers, schemes, runtime and hardware structures set.
Output
KOmegaModel( turbulence, k_eqn, ω_eqn, state )
– Turbulence model structure.
XCALibre.ModelPhysics.initialise
— Methodinitialise(turbulence::KOmegaLKE, model::Physics{T,F,M,Tu,E,D,BI}, mdotf, peqn, config
) where {T,F,M,Tu,E,D,BI}
Initialisation of turbulent transport equations.
Input
turbulence
– turbulence model.model
– Physics model defined by user.mdtof
– Face mass flow.peqn
– Pressure equation.config
– Configuration structure defined by user with solvers, schemes, runtime and hardware structures set.
Output
KOmegaLKEModel( turbulence, k_eqn, ω_eqn, kl_eqn, nueffkLS, nueffkS, nueffωS, nuL, nuts, Ω, γ, fv, normU, Reυ, ∇k, ∇ω, state )
– Turbulence model structure.
XCALibre.ModelPhysics.initialise
— Methodfunction initialise(
turbulence::Laminar, model::Physics{T,F,M,Tu,E,D,BI}, mdotf, peqn, config
) where {T,F,M,Tu,E,D,BI}
return LaminarModel()
end
Initialisation of turbulent transport equations.
Input
turbulence
– turbulence model.model
– Physics model defined by user.mdtof
– Face mass flow.peqn
– Pressure equation.config
– Configuration structure defined by user with solvers, schemes, runtime and hardware structures set.
Output
LaminarModel()
– Turbulence model structure.
XCALibre.ModelPhysics.initialise
— Methodinitialise(turbulence::Smagorinsky, model::Physics{T,F,M,Tu,E,D,BI}, mdotf, peqn, config
) where {T,F,M,Tu,E,D,BI}
Initialisation of turbulent transport equations.
Input
turbulence
: turbulence model.model
: Physics model defined by user.mdtof
: Face mass flow.peqn
: Pressure equation.config
: Configuration structure defined by user with solvers, schemes, runtime and hardware structures set.
Output
Returns a structure holding the fields and data needed for this model
SmagorinskyModel(
turbulence,
Δ,
ModelState((), false)
)
XCALibre.ModelPhysics.thermo_Psi!
— Methodthermo_Psi!(model::Physics{T,F,M,Tu,E,D,BI}, Psif::FaceScalarField)
where {T,F<:AbstractCompressible,M,Tu,E,D,BI}
Function updates the value of Psi.
Input
model
– Physics model defined by user.Psif
– Compressibility factor FaceScalarField.
Algorithm
Weakly compressible currently uses the ideal gas equation for establishing the compressibility factor where $\rho = p * \Psi$. $\Psi$ is calculated from the sensible enthalpy, reference temperature and fluid model specified $C_p$ and $R$ value where $R$ is calculated from $C_p$ and $\gamma$ specified in the fluid model.
XCALibre.ModelPhysics.thermo_Psi!
— Methodthermo_Psi!(model::Physics{T,F,M,Tu,E,D,BI}, Psi::ScalarField)
where {T,F<:AbstractCompressible,M,Tu,E,D,BI}
Model updates the value of Psi.
Input
model
– Physics model defined by user.Psi
– Compressibility factor ScalarField.
Algorithm
Weakly compressible currently uses the ideal gas equation for establishing the compressibility factor where $\rho = p * \Psi$. $\Psi$ is calculated from the sensible enthalpy, reference temperature and fluid model specified $C_p$ and $R$ value where $R$ is calculated from $C_p$ and $\gamma$ specified in the fluid model.
XCALibre.ModelPhysics.turbulence!
— Methodturbulence!(rans::KOmegaLKEModel, model::Physics{T,F,M,Turb,E,D,BI}, S, prev, time, config ) where {T,F,M,Turb<:AbstractTurbulenceModel,E,D,BI}
Run turbulence model transport equations.
Input
rans::KOmegaLKEModel
– KOmega turbulence model.model
– Physics model defined by user.S
– Strain rate tensor.prev
– Previous field.time
–config
– Configuration structure defined by user with solvers, schemes, runtime and hardware structures set.
XCALibre.ModelPhysics.turbulence!
— Methodturbulence!(les::KEquationModel, model::Physics{T,F,M,Tu,E,D,BI}, S, S2, prev, time, config
) where {T,F,M,Tu<:AbstractTurbulenceModel,E,D,BI}
Run turbulence model transport equations.
Input
les::KEquationModel
– KEquation LES turbulence model.model
– Physics model defined by user.S
– Strain rate tensor.S2
– Square of the strain rate magnitude.prev
– Previous field.time
–config
– Configuration structure defined by user with solvers, schemes, runtime and hardware structures set.
XCALibre.ModelPhysics.turbulence!
— Methodturbulence!(rans::KOmegaModel, model::Physics{T,F,M,Tu,E,D,BI}, S, prev, time, config
) where {T,F,M,Tu<:AbstractTurbulenceModel,E,D,BI}
Run turbulence model transport equations.
Input
rans::KOmegaModel
– KOmega turbulence model.model
– Physics model defined by user.S
– Strain rate tensor.prev
– Previous field.time
–config
– Configuration structure defined by user with solvers, schemes, runtime and hardware structures set.
XCALibre.ModelPhysics.turbulence!
— Methodturbulence!(rans::LaminarModel, model::Physics{T,F,M,Tu,E,D,BI}, S, prev, time, config
) where {T,F,M,Tu<:Laminar,E,D,BI}
Run turbulence model transport equations.
Input
rans::LaminarModel
– Laminar turbulence model.model
– Physics model defined by user.S
– Strain rate tensor.prev
– Previous field.time
–config
– Configuration structure defined by user with solvers, schemes, runtime and hardware structures set.
XCALibre.ModelPhysics.turbulence!
— Methodturbulence!(les::SmagorinskyModel, model::Physics{T,F,M,Tu,E,D,BI}, S, S2, prev, time, config
) where {T,F,M,Tu<:AbstractTurbulenceModel,E,D,BI}
Run turbulence model transport equations.
Input
les::SmagorinskyModel
:Smagorinsky
LES turbulence model.model
: Physics model defined by user.S
: Strain rate tensor.S2
: Square of the strain rate magnitude.prev
: Previous field.time
: current simulation timeconfig
: Configuration structure defined by user with solvers, schemes, runtime and hardware structures set.
XCALibre.Simulate.Configuration
— Type@kwdef struct Configuration{SC,SL,RT,HW}
schemes::SC
solvers::SL
runtime::RT
hardware::HW
end
The Configuration
type is passed to all flow solvers and provides all the relevant information to run a simulation.
Inputs
schemes::NamedTuple
this keyword argument is used to pass distretisation scheme information to flow solvers. See Numerical setup for details.solvers::NamedTuple
this keyword argument is used to pass the configurations for the linear solvers for each field information to flow solvers. See Runtime and solvers for details.runtime::NamedTuple
this keyword argument is used to pass runtime information to the flow solvers. See Runtime and solvers for details.hardware::NamedTuple
this keyword argument is used to pass the hardware configuration and backend settings to the flow solvers. See Pre-processing for details.
Example
config = Configuration(
solvers=solvers, schemes=schemes, runtime=runtime, hardware=hardware, boundaries=BCs)
XCALibre.Simulate.Hardware
— Typehardware = Hardware(backend, workgroup)
Struct used to configure the backend.
Inputs
backend
: used to specify the backend e.g.CPU()
,CUDABackend()
or other backends supported byKernelAbstraction.jl
workgroup::Int
this is an integer specifying the number of workers that cooperate in a parallel run. For GPUs this could be set to the size of the device's warp e.g.workgroup = 32
. On CPUs, the default value inKernelAbstractions.jl
is currentlyworkgroup = 1024
.
Output
This function returns a Hardware
object with the fields backend
and workgroup
which are accessed by internally in XCALibre.jl
to execute a given kernel in the target backend
.
XCALibre.Solvers.cpiso!
— Methodcpiso!(model, config;
output=VTK(), pref=nothing, ncorrectors=0, inner_loops=0)
Compressible and transient variant of the PISO algorithm with a sensible enthalpy transport equation for the energy.
Input arguments
model
reference to aPhysics
model defined by the user.config
Configuration structure defined by the user with solvers, schemes, runtime and hardware structures configuration details.output
select the format used for simulation results fromVTK()
orOpenFOAM
(default =VTK()
)pref
Reference pressure value for cases that do not have a pressure defining BC. Incompressible solvers only (default =nothing
)ncorrectors
number of non-orthogonality correction loops (default =0
)inner_loops
number to inner loops used in transient solver based on PISO algorithm (default =0
)
Output
Ux
Vector of x-velocity residuals for each iteration.Uy
Vector of y-velocity residuals for each iteration.Uz
Vector of y-velocity residuals for each iteration.p
Vector of pressure residuals for each iteration.
XCALibre.Solvers.csimple!
— Methodcsimple!(
model_in, config;
output=VTK(), pref=nothing, ncorrectors=0, inner_loops=0
)
Compressible variant of the SIMPLE algorithm with a sensible enthalpy transport equation for the energy.
Input arguments
model
reference to aPhysics
model defined by the user.config
Configuration structure defined by the user with solvers, schemes, runtime and hardware structures configuration details.output
select the format used for simulation results fromVTK()
orOpenFOAM
(default =VTK()
)pref
Reference pressure value for cases that do not have a pressure defining BC. Incompressible solvers only (default =nothing
)ncorrectors
number of non-orthogonality correction loops (default =0
)inner_loops
number to inner loops used in transient solver based on PISO algorithm (default =0
)
Output
Ux
Vector of x-velocity residuals for each iteration.Uy
Vector of y-velocity residuals for each iteration.Uz
Vector of y-velocity residuals for each iteration.p
Vector of pressure residuals for each iteration.e
Vector of energy residuals for each iteration.
XCALibre.Solvers.piso!
— Methodcpiso!(model, config;
output=VTK(), pref=nothing, ncorrectors=0, inner_loops=0)
Incompressible and transient variant of the SIMPLE algorithm to solving coupled momentum and mass conservation equations.
Input arguments
model
reference to aPhysics
model defined by the user.config
Configuration structure defined by the user with solvers, schemes, runtime and hardware structures configuration details.output
select the format used for simulation results fromVTK()
orOpenFOAM
(default =VTK()
)pref
Reference pressure value for cases that do not have a pressure defining BC. Incompressible solvers only (default =nothing
)ncorrectors
number of non-orthogonality correction loops (default =0
)inner_loops
number to inner loops used in transient solver based on PISO algorithm (default =0
)
Output
Ux
Vector of x-velocity residuals for each iteration.Uy
Vector of y-velocity residuals for each iteration.Uz
Vector of y-velocity residuals for each iteration.p
Vector of pressure residuals for each iteration.
XCALibre.Solvers.run!
— Methodfunction run!(
model::Physics, config;
output=VTK(), pref=nothing, ncorrectors=0, inner_loops=0
)
# here an internal function is used for solver dispatch
return residuals
end
This is the top level API function to initiate a simulation. It uses the user-provided model
defined as a Physics
object to dispatch to the appropriate solver.
Dispatched flow solvers
- Steady incompressible (SIMPLE algorithm for coupling)
- Transient incompressible (PISO algorithm for coupling)
- Steady weakly compressible (SIMPLE algorithm for coupling)
- Transient weakly compressible (PISO algorithm for coupling)
Input arguments
model
reference to aPhysics
model defined by the user.config
Configuration structure defined by the user with solvers, schemes, runtime and hardware structures configuration details.output
select the format used for simulation results fromVTK()
orOpenFOAM()
(default =VTK()
)pref
Reference pressure value for cases that do not have a pressure defining BC. Incompressible solvers only (default =nothing
)ncorrectors
number of non-orthogonality correction loops (default =0
)inner_loops
number to inner loops used in transient solver based on PISO algorithm (default =0
)
Output
This function returns a NamedTuple
for accessing the residuals (e.g. residuals.Ux
). The fields available within the returned residuals
tuple depend on the solver used. For example, for an incompressible solver, a x-momentum equation residual can be retrieved accessing the Ux
field i.e. residuals.Ux
. Look at reference guide for each dispatch method to find out which fields are available.
Example
residuals = run!(model, config)
# to access the pressure residual
residuals.p
XCALibre.Solvers.run!
— Methodrun!(
model::Physics{T,F,M,Tu,E,D,BI}, config;
output=VTK(), pref=nothing, ncorrectors=0, inner_loops=0
) where{T<:Steady,F<:Incompressible,M,Tu,E,D,BI} =
begin
residuals = simple!(model, config, pref=pref)
return residuals
end
Calls the incompressible steady solver using the SIMPLE algorithm.
Input
model
represents thePhysics
model defined by user.config
Configuration structure defined by user with solvers, schemes, runtime and hardware structures configuration details.output
select the format used for simulation results fromVTK()
orOpenFOAM()
(default =VTK()
)pref
Reference pressure value for cases that do not have a pressure defining BC. Incompressible solvers only.
Output
This function returns a NamedTuple
for accessing the residuals (e.g. residuals.Ux
) with the following entries:
Ux
- Vector of x-velocity residuals for each iteration.Uy
- Vector of y-velocity residuals for each iteration.Uz
- Vector of y-velocity residuals for each iteration.p
- Vector of pressure residuals for each iteration.
XCALibre.Solvers.run!
— Methodrun!(
model::Physics{T,F,M,Tu,E,D,BI}, config;
output=VTK(), pref=nothing, ncorrectors=0, inner_loops=0
) where{T<:Steady,F<:WeaklyCompressible,M,Tu,E,D,BI} =
begin
residuals = csimple!(model, config, pref=pref); #, pref=0.0)
return residuals
end
Calls the compressible steady solver using the SIMPLE algorithm for weakly compressible fluids.
Input
model
represents thePhysics
model defined by user.config
Configuration structure defined by user with solvers, schemes, runtime and hardware structures configuration details.output
select the format used for simulation results fromVTK()
orOpenFOAM()
(default =VTK()
)pref
Reference pressure value for cases that do not have a pressure defining BC. Incompressible solvers only.
Output
This function returns a NamedTuple
for accessing the residuals (e.g. residuals.Ux
) with the following entries:
Ux
- Vector of x-velocity residuals for each iteration.Uy
- Vector of y-velocity residuals for each iteration.Uz
- Vector of y-velocity residuals for each iteration.p
- Vector of pressure residuals for each iteration.e
- Vector of energy residuals for each iteration.
XCALibre.Solvers.run!
— Methodrun!(
model::Physics{T,F,M,Tu,E,D,BI}, config;
output=VTK(), pref=nothing, ncorrectors=0, inner_loops=0
) where{T<:Transient,F<:Incompressible,M,Tu,E,D,BI} =
begin
residuals = piso!(model, config, pref=pref); #, pref=0.0)
return residuals
end
Calls the incompressible transient solver using the PISO algorithm.
Input
model
represents thePhysics
model defined by user.config
Configuration structure defined by user with solvers, schemes, runtime and hardware structures configuration details.output
select the format used for simulation results fromVTK()
orOpenFOAM()
(default =VTK()
)pref
Reference pressure value for cases that do not have a pressure defining BC. Incompressible solvers only.
Output
This function returns a NamedTuple
for accessing the residuals (e.g. residuals.Ux
) with the following entries:
Ux
- Vector of x-velocity residuals for each iteration.Uy
- Vector of y-velocity residuals for each iteration.Uz
- Vector of y-velocity residuals for each iteration.p
- Vector of pressure residuals for each iteration.
XCALibre.Solvers.run!
— Methodrun!(
model::Physics{T,F,M,Tu,E,D,BI};
output=VTK(), pref=nothing, ncorrectors=0, inner_loops=0
) where{T<:Transient,F<:WeaklyCompressible,M,Tu,E,D,BI} =
begin
residuals = cpiso!(model, config)
return residuals
end
Calls the compressible transient solver using the PISO algorithm for weakly compressible fluids.
Input
model
represents thePhysics
model defined by user.config
Configuration structure defined by user with solvers, schemes, runtime and hardware structures configuration details.output
select the format used for simulation results fromVTK()
orOpenFOAM()
(default =VTK()
)pref
Reference pressure value for cases that do not have a pressure defining BC. Incompressible solvers only.
Output
This function returns a NamedTuple
for accessing the residuals (e.g. residuals.Ux
) with the following entries:
Ux
- Vector of x-velocity residuals for each iteration.Uy
- Vector of y-velocity residuals for each iteration.Uz
- Vector of y-velocity residuals for each iteration.p
- Vector of pressure residuals for each iteration.e
- Vector of energy residuals for each iteration.
XCALibre.Solvers.simple!
— Methodsimple!(model_in, config;
output=VTK(), pref=nothing, ncorrectors=0, inner_loops=0)
Incompressible variant of the SIMPLE algorithm to solving coupled momentum and mass conservation equations.
Input arguments
model
reference to aPhysics
model defined by the user.config
Configuration structure defined by the user with solvers, schemes, runtime and hardware structures configuration details.output
select the format used for simulation results fromVTK()
orOpenFOAM
(default =VTK()
)pref
Reference pressure value for cases that do not have a pressure defining BC. Incompressible solvers only (default =nothing
)ncorrectors
number of non-orthogonality correction loops (default =0
)inner_loops
number to inner loops used in transient solver based on PISO algorithm (default =0
)
Output
This function returns a NamedTuple
for accessing the residuals (e.g. residuals.Ux
) with the following entries:
Ux
Vector of x-velocity residuals for each iteration.Uy
Vector of y-velocity residuals for each iteration.Uz
Vector of y-velocity residuals for each iteration.p
Vector of pressure residuals for each iteration.
XCALibre.Postprocess.boundary_average
— Methodfunction boundary_average(patch::Symbol, field, config; time=0)
# Extract mesh object
mesh = field.mesh
# Determine ID (index) of the boundary patch
ID = boundary_index(mesh.boundaries, patch)
@info "calculating average on patch: $patch at index $ID"
boundary = mesh.boundaries[ID]
(; IDs_range) = boundary
# Create face field of same type provided by user (scalar or vector)
sum = nothing
if typeof(field) <: VectorField
faceField = FaceVectorField(mesh)
sum = zeros(_get_float(mesh), 3) # create zero vector
else
faceField = FaceScalarField(mesh)
sum = zero(_get_float(mesh)) # create zero
end
# Interpolate CFD results to boundary
interpolate!(faceField, field, config)
correct_boundaries!(faceField, field, field.BCs, time, config)
# Calculate the average
for fID ∈ IDs_range
sum += faceField[fID]
end
ave = sum/length(IDs_range)
# return average
return ave
end
XCALibre.Postprocess.pressure_force
— Methodpressure_force(patch::Symbol, p::ScalarField, rho)
Function to calculate the pressure force acting on a given patch/boundary.
Input arguments
patch::Symbol
name of the boundary of interest (as aSymbol
)p::ScalarField
pressure fieldrho
density. Set to 1 for incompressible solvers
XCALibre.Postprocess.viscous_force
— Methodviscous_force(patch::Symbol, U::VectorField, rho, ν, νt)
Function to calculate the pressure force acting on a given patch/boundary.
Input arguments
patch::Symbol
name of the boundary of interest (as aSymbol
)U::VectorField
pressure fieldrho
density. Set to 1 for incompressible solversν
laminar viscosity of the fluidνt
eddy viscosity from turbulence models. Pass ConstantScalar(0) for laminar flows