System
Types and Constants
VoronoiFVM.AbstractSystem
— Typeabstract type AbstractSystem
Abstract type for finite volume system structure
VoronoiFVM.DenseSystem
— Typemutable struct DenseSystem{Tv, Ti, Tm} <: VoronoiFVM.AbstractSystem{Tv,Ti,Tm}
Structure holding data for finite volume system solution. Information on species distribution is kept in dense matrices, and the solution array is of type Array{2}.
Unlike in the SparseSystem, the system matrix handles exactly those degrees of freedom which correspond to unknowns, and dummy degrees of freedom where unknowns are not defined. Handling of the sparse matrix structures for the bookeeping of the unknowns has less overhead, but additional dummy equations are added to the system matrix.
grid::Any
Grid
physics::VoronoiFVM.Physics
Physics data
boundary_values::Array{Tv,2} where Tv
Array of boundary values
boundary_factors::Array{Tv,2} where Tv
Array of boundary factors
region_species::Array{Int8,2}
Full matrix containing species numbers for inner regions
bregion_species::Array{Int8,2}
Full matrix containing species numbers for boundary regions
node_dof::Array{Int8,2}
Full matrix containing degree of freedom numbers for each node
matrix::ExtendableSparse.ExtendableSparseMatrix{Tv,Tm} where Tm where Tv
Jacobi matrix for nonlinear problem
species_homogeneous::Bool
Flag which says if the number of unknowns per node is constant
update::Array{Tv,2} where Tv
Solution vector holding Newton update
residual::Array{Tv,2} where Tv
Solution vector holding Newton residual
cellnodefactors::Array{Tv,2} where Tv
celledgefactors::Array{Tv,2} where Tv
bfacenodefactors::Array{Tv,2} where Tv
generic_matrix::SparseArrays.SparseMatrixCSC
generic_matrix_colors::Array{T,1} where T
VoronoiFVM.DenseSystem
— MethodDenseSystem(grid::Any, physics::VoronoiFVM.Physics; matrixindextype) -> VoronoiFVM.DenseSystem{_A,_B,Int32} where _B where _A
Constructor for DenseSystem.
VoronoiFVM.SparseSolutionArray
— Typestruct SparseSolutionArray{Tv, Ti} <: AbstractArray{Tv,2}
Struct holding solution information for SparseSystem. Solution is stored in a sparse matrix structure.
This class plays well with the abstract array interface.
node_dof::SparseArrays.SparseMatrixCSC{Tv,Ti} where Ti where Tv
Sparse matrix holding actual data.
VoronoiFVM.SparseSystem
— Typemutable struct SparseSystem{Tv, Ti, Tm} <: VoronoiFVM.AbstractSystem{Tv,Ti,Tm}
Structure holding data for finite volume system solution. Information on species distribution is kept in sparse matrices, and the solution array is of type SparseSolutionArray, i.e. effectively it is a sparse matrix.
Unlike in the DenseSystem, the system matrix handles exactly those degrees of freedom which correspond to unknowns. However, handling of the sparse matrix structures for the bookkeeping of the unknowns creates overhead.
grid::Any
Grid
physics::VoronoiFVM.Physics
Physics data
boundary_values::Array{Tv,2} where Tv
Array of boundary values
boundary_factors::Array{Tv,2} where Tv
Array of boundary factors
region_species::SparseArrays.SparseMatrixCSC{Int8,Ti} where Ti
Sparse matrix containing species numbers for inner regions
bregion_species::SparseArrays.SparseMatrixCSC{Int8,Ti} where Ti
Sparse matrix containing species numbers for boundary regions
node_dof::SparseArrays.SparseMatrixCSC{Int8,Tm} where Tm
Sparse matrix containing degree of freedom numbers for each node
matrix::ExtendableSparse.ExtendableSparseMatrix{Tv,Tm} where Tm where Tv
Jacobi matrix for nonlinear problem
species_homogeneous::Bool
Flag which says if the number of unknowns per node is constant
update::VoronoiFVM.SparseSolutionArray{Tv,Tm} where Tm where Tv
Solution vector holding Newton update
residual::VoronoiFVM.SparseSolutionArray{Tv,Tm} where Tm where Tv
Solution vector holding Newton residual
cellnodefactors::Array{Tv,2} where Tv
celledgefactors::Array{Tv,2} where Tv
bfacenodefactors::Array{Tv,2} where Tv
generic_matrix::SparseArrays.SparseMatrixCSC
generic_matrix_colors::Array{T,1} where T
VoronoiFVM.SparseSystem
— MethodSparseSystem(grid::Any, physics::VoronoiFVM.Physics; matrixindextype) -> VoronoiFVM.SparseSystem{_A,_B,Int32} where _B where _A
Constructor for SparseSystem.
VoronoiFVM.BNode
— Typemutable struct BNode{Tv, Ti} <: VoronoiFVM.AbstractGeometryItem{Tv,Ti}
Structure holding local boundary node information.
index::Any
Index in grid
region::Any
Boundary region number
nspec::Any
Number of species defined in node
coord::Array{Tv,2} where Tv
Grid coordinates
VoronoiFVM.Edge
— Typemutable struct Edge{Tv, Ti} <: VoronoiFVM.AbstractGeometryItem{Tv,Ti}
Structure holding local edge information.
index::Any
Index in grid
node::Array{Ti,1} where Ti
Index
region::Any
Inner region number corresponding to edge
nspec::Any
Number of species defined in edge
icell::Any
Number of discretization cell the edge is invoked from
coord::Array{Tv,2} where Tv
Grid coordinates
VoronoiFVM.Node
— Typemutable struct Node{Tv, Ti} <: VoronoiFVM.AbstractGeometryItem{Tv,Ti}
Structure holding local node information.
index::Any
Index in grid
region::Any
Inner region number
nspec::Any
Number of species defined in node
icell::Any
Number of discretization cell the node is invoked from
coord::Array{Tv,2} where Tv
Grid coordinates
VoronoiFVM.AssemblyError
— Typestruct AssemblyError <: Exception
Exception thrown if error occured during assembly (e.g. domain error)
VoronoiFVM.ConvergenceError
— Typestruct ConvergenceError <: Exception
Exception thrown if Newton's method convergence fails.
VoronoiFVM.EmbeddingError
— Typestruct EmbeddingError <: Exception
Exception thrown if embedding fails
VoronoiFVM.FactorizationError
— Typestruct FactorizationError <: Exception
Exception thrown if error occured during factorization.
VoronoiFVM.NewtonControl
— Typemutable struct NewtonControl
Control parameters for Newton method.
Newton's method solves $F(u)=0$ by the iterative procedure $u_{i+1}=u_{i} - d_i F'(u_i)^{-1}F(u_i)$ starting with some inital value $u_0$, where $d_i$ is the damping.
tol_absolute::Float64
Tolerance (in terms of norm of Newton update): terminate if $\Delta_i=||u_{i+1}-u_i||_\infty <$
tol_absolute
.Default value:
1.0e-10
.
tol_relative::Float64
Tolerance (relative to the first residual): terminate if $\Delta_i/\Delta_0<$
tol_relative
.Default value:
1.0e-10
.
tol_round::Float64
Tolerance for roundoff error detection: terminate if $|\;||u_{i+1}||_1 - ||u_{i}||_1\;|/ ||u_{i}||_1<$
tol_round
occuredmax_round
times in a row.Default value:
1.0e-10
.
tol_mono::Float64
Tolerance for monotonicity test: terminat with error if $\Delta_i/\Delta_{i-1}>$
1/tol_mono
.Default value: 1.0e-3
damp_initial::Float64
Initial damping parameter $d_0$.
Default value:
1.0
.To handle convergence problems, set this to a value less than 1.
damp_growth::Float64
Damping parameter growth factor: $d_{i+1}=\min(d_i\cdot$
max_growth
$,1)$Default value:
1.2
Generally it should be set to a value between 1 and 2.
max_iterations::Int32
Maximum number of iterations.
Default value: 100
max_lureuse::Int32
Maximum number of reuses of lu factorization. It this value is 0, linear systems are solved by a sparse direct solver, and it's LU factorization is called in every Newton step.
Otherwise, a BICGstab iterative method is used for linear system solution with an LU factorization as preconditioner which is updated only every
max_lureuse
Newton step.Default value: 0.
max_round::Int32
Maximum number of consecutive iterations within roundoff error tolerance
Default value: 1000 (effectively switching of this criterion).
tol_linear::Float64
Tolerance of iterative linear solver.
Default value:
1.0e-4
.
verbose::Bool
Verbosity flag
Default value:
false
.
handle_exceptions::Bool
Handle exceptions in
embed!
andevolve!
methods. Iftrue
, exceptions in Newton solves are catched, the embedding resp. time step is lowered, and solution is retried.Default value:
false
Δp::Float64
Initial parameter step for
embed!
method.Default value:
1.0
Δp_max::Float64
Maximal parameter step for
embed!
methodDefault value:
1.0
Δp_min::Float64
Minimal parameter step for
embed!
method.Default value:
1.0e-3
Δt::Float64
Time step for
evolve!
method.Default value: 0.1
Δt_max::Float64
Maximal time step for
evolve!
method.Default value: 1
Δt_min::Float64
Minimal time step for
evolve!
method.Default value:
1.0e-3
Δt_grow::Float64
Maximal step size growth for
evolve!
method.Default: 1.2
Δu_opt::Float64
Optimal size of update for
evolve!
method. The algorithm keeps this value approximately constant.Default: 0.1
edge_cutoff::Float64
Edge cutoff for rectangular triangles.
Default value:
0.0
.
umfpack_pivot_tolerance::Float64
Pivot tolerance for umfpack
Default value: 0.1.
VoronoiFVM.NewtonControl
— MethodNewtonControl(this::Any) -> Any
Default constructor
VoronoiFVM.TestFunctionFactory
— Typemutable struct TestFunctionFactory{Tv}
Data structure containing DenseSystem used to calculate test functions for boundary flux calculations.
system::VoronoiFVM.AbstractSystem{Tv,Ti,Tm} where Tm<:Integer where Ti<:Integer where Tv
Original system
tfsystem::VoronoiFVM.DenseSystem{Tv,Ti,Tm} where Tm where Ti where Tv
Test function system
VoronoiFVM.TestFunctionFactory
— MethodTestFunctionFactory(system::VoronoiFVM.AbstractSystem{Tv,Ti,Tm} where Tm<:Integer where Ti<:Integer) -> VoronoiFVM.TestFunctionFactory{_A} where _A
Constructor for TestFunctionFactory from System
VoronoiFVM.Dirichlet
— ConstantConstant to be used as boundary condition factor to mark Dirichlet boundary conditons.
Methods
VoronoiFVM.System
— MethodSystem(grid::Any, physics::VoronoiFVM.Physics; unknown_storage, matrixindextype) -> Union{VoronoiFVM.DenseSystem{_A,_B,_C} where _C where _B where _A, VoronoiFVM.SparseSystem{_A,_B,_C} where _C where _B where _A}
Create Finite Volume System.
grid
: 1D/2D/3D discretization gridphysics
: Physics struct containing node and edge callbacksunknown_storage
: string or symbol::dense
: solution vector is annspecies
xnnodes
dense matrix:sparse
: solution vector is annspecies
xnnodes
sparse matrix
VoronoiFVM.boundary_dirichlet!
— Methodboundary_dirichlet!(this::VoronoiFVM.AbstractSystem, ispec::Integer, ibc::Integer, val::Any) -> Any
Set Dirichlet boundary conditon for species ispec at boundary ibc.
VoronoiFVM.boundary_neumann!
— Methodboundary_neumann!(this::VoronoiFVM.AbstractSystem, ispec::Integer, ibc::Integer, val::Any) -> Any
Set Neumann boundary conditon for species ispec at boundary ibc.
VoronoiFVM.boundary_robin!
— Methodboundary_robin!(this::VoronoiFVM.AbstractSystem, ispec::Integer, ibc::Integer, fac::Any, val::Any) -> Any
Set Robin boundary conditon for species ispec at boundary ibc.
VoronoiFVM.data
— Methoddata(this::VoronoiFVM.AbstractSystem) -> Any
Retrieve user data record.
VoronoiFVM.enable_boundary_species!
— Methodenable_boundary_species!(this::VoronoiFVM.AbstractSystem, ispec::Integer, bregions::AbstractArray{T,1} where T)
Add species to a list of boundary regions. Species numbers for bulk and boundary species have to be distinct.
VoronoiFVM.enable_species!
— Methodenable_species!(this::VoronoiFVM.AbstractSystem, ispec::Integer, regions::AbstractArray{T,1} where T)
Add species to a list of bulk regions. Species numbers for bulk and boundary species have to be distinct.
VoronoiFVM.is_boundary_species
— Methodis_boundary_species(this::VoronoiFVM.AbstractSystem, ispec::Integer) -> Bool
Check if species number corresponds to boundary species.
VoronoiFVM.is_bulk_species
— Methodis_bulk_species(this::VoronoiFVM.AbstractSystem, ispec::Integer) -> Bool
Check if species number corresponds bulk species.
VoronoiFVM.isdof
— Methodisdof(this::VoronoiFVM.AbstractSystem, ispec::Any, inode::Any) -> Bool
Check if degree of freedom is defined.
VoronoiFVM.num_species
— Methodnum_species(a::AbstractArray) -> Any
Number of species (size of first dimension) of solution array.
VoronoiFVM.num_species
— Methodnum_species(this::VoronoiFVM.AbstractSystem) -> Any
Number of species in system
VoronoiFVM.update_grid!
— MethodUpdate grid (e.g. after rescaling of coordinates).
VoronoiFVM.Grid
— FunctionGrid=ExtendableGrids.simplexgrid
Re-Export of ExtendableGrids.simplexgrid
Base.reshape
— Methodreshape(v, sys)
Reshape vector to fit as solution to system.
VoronoiFVM.dof
— Methoddof(a::Array{Tv,2}, ispec::Integer, K::Integer) -> Any
Get number of degree of freedom.
VoronoiFVM.num_dof
— Methodnum_dof(this::VoronoiFVM.DenseSystem) -> Int64
Number of degrees of freedom for system.
VoronoiFVM.unknowns
— Methodunknowns(Tu::Type, sys::VoronoiFVM.DenseSystem{Tv,Ti,Tm} where Tm where Ti; inival) -> Any
Create a solution vector for dense system. If inival is not specified, the entries of the returned vector are undefined.
VoronoiFVM.unknowns
— MethodCreate a solution vector for dense system. If inival is not specified, the entries of the returned vector are undefined.
VoronoiFVM.values
— Methodvalues(a::Array) -> Union{Base.ReshapedArray{_A,1,_B,Tuple{}} where _B<:AbstractArray where _A, Array{T,M} where M where T}
Array of values in solution array.
Base.copy
— MethodCreate a copy of solution array
Base.getindex
— Methodgetindex(a::VoronoiFVM.SparseSolutionArray, ispec::Integer, inode::Integer) -> Any
Accessor for solution array.
Base.reshape
— Methodreshape(v, sys)
Reshape vector to fit as solution to system.
Base.setindex!
— Methodsetindex!(a::VoronoiFVM.SparseSolutionArray, v::Any, ispec::Integer, inode::Integer) -> Union{Nothing, VoronoiFVM.SparseSolutionArray}
Accessor for solution array.
Base.size
— Methodsize(a::VoronoiFVM.SparseSolutionArray) -> Tuple{Int64,Int64}
Return size of solution array.
VoronoiFVM.dof
— Methoddof(a, i, j)
Get number of degree of freedom. Return 0 if species is not defined in node.
VoronoiFVM.getdof
— Methodgetdof(a::VoronoiFVM.SparseSolutionArray, i::Integer) -> Any
Return value for degree of freedom.
VoronoiFVM.num_dof
— Methodnum_dof(this::VoronoiFVM.SparseSystem) -> Any
Number of degrees of freedom for system.
VoronoiFVM.setdof!
— Methodsetdof!(a::VoronoiFVM.SparseSolutionArray, v::Any, i::Integer) -> Any
Set value for degree of freedom.
VoronoiFVM.unknowns
— Methodunknowns(Tu, sys; inival)
Create a solution vector for sparse system with given type. If inival is not specified, the entries of the returned vector are undefined.
VoronoiFVM.unknowns
— MethodCreate a solution vector for sparse system. The entries of the returned vector are undefined.
VoronoiFVM.values
— Methodvalues(a::VoronoiFVM.SparseSolutionArray) -> Array{Tv,1} where Tv
Array of values in solution array.
VoronoiFVM.meas
— Methodmeas(edge::VoronoiFVM.Edge) -> Any
Calculate the length of an edge.
VoronoiFVM.num_species
— Methodnum_species(edge::VoronoiFVM.Edge) -> Any
Return number of species for edge
VoronoiFVM.viewK
— MethodviewK(nspec::Any, u::AbstractArray) -> Any
Solution view on first edge node
VoronoiFVM.viewK
— MethodviewK(edge::VoronoiFVM.Edge, u::AbstractArray) -> Any
Solution view on first edge node
VoronoiFVM.viewL
— MethodviewL(nspec::Any, u::AbstractArray) -> Any
Solution view on second edge node
VoronoiFVM.viewL
— MethodviewL(edge::VoronoiFVM.Edge, u::AbstractArray) -> Any
Solution view on second edge node
VoronoiFVM.embed!
— Methodembed!(solution, xinival, this; control, pre, post)
Solution method for instance of abstract system.
Perform solution via parameter embedding, calling solve! for each value of the parameter p from interval (0,1). The user is responsible for the interpretation of the parameter. The optional pre()
callback can be used to communicate its value to the system. The optional post()
callback method can be used to perform various postprocessing steps.
If control.handle_error
is true, solve!
throws an error, and stepsize control.Δp
is lowered, and solve!
is called again with a smaller parameter value. If control.Δp<control.Δp_min
, embed!
is aborted with error.
VoronoiFVM.eval_and_assemble
— Methodeval_and_assemble(system, U, UOld, F, tstep; edge_cutoff)
Main assembly method.
Evaluate solution with result in right hand side F and assemble matrix into system.matrix.
VoronoiFVM.evolve!
— Methodevolve!(solution::AbstractArray{Tv,2}, xinival::AbstractArray{Tv,2}, this::VoronoiFVM.AbstractSystem{Tv,Ti,Tm} where Tm<:Integer where Ti<:Integer, times::AbstractArray{T,1} where T; control, pre, post, delta) -> AbstractArray{Tv,2}
Time dependent solver for abstract system
Use implicit Euler method with + damped Newton's method to solve time dependent problem.
VoronoiFVM.integrate
— Methodintegrate(this, F, U; boundary)
Integrate function F
of solution vector over domain or boundary The result contains the integral for each species separately.
VoronoiFVM.solve!
— MethodSolution method for instance of abstract system.
Perform solution of stationary system (if tstep==Inf
) or one tine step of implicit Euler time step system using Newton's method with damping. Initial damping is chosen according control.damp_initial
.
VoronoiFVM.value
— FunctionExtract value from dual number. Use to debug physics callbacks. Re-exported from ForwardDiff.jl
VoronoiFVM.integrate
— Methodintegrate(this, tf, U, Uold, tstep)
Calculate test function integral for transient solution.
VoronoiFVM.integrate
— Methodintegrate(this, tf, U)
Calculate test function integral for steady state solution.
VoronoiFVM.integrate_stdy
— MethodSteady state part of test function integral.
VoronoiFVM.integrate_tran
— Methodintegrate_tran(this, tf, U)
Calculate transient part of test function integral.
VoronoiFVM.testfunction
— Methodtestfunction(factory::VoronoiFVM.TestFunctionFactory{Tv}, bc0::Any, bc1::Any) -> Any
Create testfunction which has Dirichlet zero boundary conditions for boundary regions in bc0 and Dirichlet one boundary conditions for boundary regions in bc1.