Time-dependent Solvers
This section describes how to solve time-dependent (non-stationary) PDEs using the high-level API.
Approaches to Time-dependent Problems
- Manual approach:
- Add custom time-derivative terms to the problem (e.g., a mass matrix as a
BilinearOperator
andLinearOperator
s for previous time steps). - If more than one previous time step is needed (e.g., for BDF2 or multi-step methods), you must manage the storage and update of previous solutions manually.
- Add custom time-derivative terms to the problem (e.g., a mass matrix as a
- Automatic approach:
- Reframe the
ProblemDescription
as an ODE problem and solve it using DifferentialEquations.jl via theExtendableFEMDiffEQExt.jl
extension.
- Reframe the
Several time-dependent examples are available, including both approaches. See, for example, Example103 (Burgers' equation) and Example205 (Heat equation).
Using SciMLBase.ODEProblem and DifferentialEquations.jl
You can reframe the spatial part of your PDE as the right-hand side of an ODEProblem
. The ProblemDescription
then describes the spatial operators and right-hand side:
\[\begin{aligned} M u_t(t) & = b(u(t)) - A(u(t)) u(t) \end{aligned}\]
where:
A
andb
are the assembled (linearized) spatial operator and right-hand side operators in theProblemDescription
(note:A
comes with a minus sign).M
is the mass matrix, which can be customized (as long as it remains constant).- Operators in the
ProblemDescription
may depend on time (e.g., if their kernels useqpinfo.time
) and will be reassembled at each time step. - To avoid unnecessary reassembly, use
store = true
for operators that do not change in time, orconstant_matrix = true
in theSolverConfiguration
to skip full matrix reassembly.
ExtendableFEM.generate_ODEProblem
— Functionfunction generate_ODEProblem(
PD::ProblemDescription,
FES,
tspan;
mass_matrix = nothing)
kwargs...)
Reframes the ProblemDescription inside the SolverConfiguration into an ODEProblem, for DifferentialEquations.jl where tspan is the desired time interval.
If no mass matrix is provided the standard mass matrix for the respective finite element space(s) for all unknowns is assembled.
Additional keyword arguments:
constant_matrix
: matrix is constant (skips reassembly and refactorization in solver). Default: falseconstant_rhs
: right-hand side is constant (skips reassembly). Default: falseinit
: initial solution (otherwise starts with a zero vector). Default: nothinginitialized
: linear system in solver configuration is already assembled (turns true after first assembly). Default: falsesametol
: tolerance to identify two solution vectors to be identical (and to skip reassemblies called by DifferentialEquations.jl). Default: 1.0e-15verbosity
: verbosity level. Default: 0
ExtendableFEM.eval_jacobian!
— MethodProvides the jacobi matrix calculation function for DifferentialEquations.jl/ODEProblem.
ExtendableFEM.eval_rhs!
— MethodProvides the rhs function for DifferentialEquations.jl/ODEProblem.
ExtendableFEM.generate_ODEProblem
— Methodfunction generate_ODEProblem(
SC::SolverConfiguration,
tspan;
mass_matrix = nothing)
kwargs...)
Reframes the ProblemDescription inside the SolverConfiguration into a SciMLBase.ODEProblem, for DifferentialEquations.jl where tspan is the desired time interval.
If no mass matrix is provided the standard mass matrix for the respective finite element space(s) for all unknowns is assembled.
Keyword arguments:
constant_matrix
: matrix is constant (skips reassembly and refactorization in solver). Default: falseconstant_rhs
: right-hand side is constant (skips reassembly). Default: falseinit
: initial solution (otherwise starts with a zero vector). Default: nothinginitialized
: linear system in solver configuration is already assembled (turns true after first assembly). Default: falsesametol
: tolerance to identify two solution vectors to be identical (and to skip reassemblies called by DifferentialEquations.jl). Default: 1.0e-15verbosity
: verbosity level. Default: 0
ExtendableFEM.jac_prototype
— MethodProvides the system matrix as prototype for the jacobian.
When using DifferentialEquations.jl, set autodiff=false
in the solver options, as automatic differentiation of the generated ODEProblem with respect to time is not currently supported.
Example: 2D Heat Equation (extracted from Example205)
The following ProblemDescription
yields the space discretization of the heat equation (including homogeneous boundary conditions; equivalent to the Poisson equation):
PD = ProblemDescription("Heat Equation")
u = Unknown("u"; name = "temperature")
assign_unknown!(PD, u)
assign_operator!(PD, BilinearOperator([grad(u)]; store = true, kwargs...))
assign_operator!(PD, HomogeneousBoundaryData(u))
Given a finite element space FES
and an initial FEVector
sol
for the unknown, the ODEProblem
for a time interval (0, T)
can be generated and solved as follows:
prob = generate_ODEProblem(PD, FES, (0, T); init = sol)
DifferentialEquations.solve(prob, ImplicitEuler(autodiff = false), dt = 1e-3, dtmin = 1e-6, adaptive = true)