20 #ifndef LIBMESH_IMPLICIT_SYSTEM_H 21 #define LIBMESH_IMPLICIT_SYSTEM_H 24 #include "libmesh/explicit_system.h" 34 template <
typename T>
class LinearSolver;
35 class StaticCondensation;
65 const std::string &
name,
66 const unsigned int number);
98 virtual void clear ()
override;
117 virtual std::string
system_type ()
const override {
return "Implicit"; }
141 virtual std::pair<unsigned int, Real>
144 #ifdef LIBMESH_ENABLE_DEPRECATED 153 #endif // LIBMESH_ENABLE_DEPRECATED 166 { libmesh_not_implemented(); }
187 { libmesh_not_implemented(); }
196 virtual std::pair<unsigned int, Real>
207 virtual std::pair<unsigned int, Real>
220 virtual std::pair<unsigned int, Real>
235 virtual std::pair<unsigned int, Real>
357 template <
typename T>
381 #endif // LIBMESH_IMPLICIT_SYSTEM_H
This is the EquationSystems class.
virtual void qoi_parameter_hessian_vector_product(const QoISet &qoi_indices, const ParameterVector ¶meters, const ParameterVector &vector, SensitivityData &product) override
For each of the system's quantities of interest q in qoi[qoi_indices], and for a vector of parameters...
void create_static_condensation_system_matrix()
Create the static condensation system matrix.
Data structure for specifying which Parameters should be independent variables in a parameter sensiti...
virtual void create_static_condensation() override
Request that static condensation be performed for this system.
virtual std::pair< unsigned int, Real > get_linear_solve_parameters() const
Data structure for specifying which Quantities of Interest should be calculated in an adjoint or a pa...
virtual std::pair< unsigned int, Real > weighted_sensitivity_solve(const ParameterVector ¶meters, const ParameterVector &weights) override
Assembles & solves the linear system(s) (dR/du)*u_w = sum(w_p*-dR/dp), for those parameters p contain...
virtual std::string system_type() const override
virtual void clear() override
Clear all the data structures associated with the system.
virtual ~ImplicitSystem()
virtual LinearSolver< Number > * get_linear_solver() const
virtual std::pair< unsigned int, Real > sensitivity_solve(const ParameterVector ¶meters) override
Assembles & solves the linear system(s) (dR/du)*u_p = -dR/dp, for those parameters contained within p...
The libMesh namespace provides an interface to certain functionality in the library.
bool zero_out_matrix_and_rhs
By default, the system will zero out the matrix and the right hand side.
virtual std::pair< unsigned int, Real > adjoint_solve(const QoISet &qoi_indices=QoISet()) override
Assembles & solves the linear system (dR/du)^T*z = dq/du, for those quantities of interest q specifie...
StaticCondensation * _sc_system_matrix
The system matrix for static condensation problems.
Parameters parameters
Parameters for the system. If a parameter is not provided, it should be retrieved from the EquationSy...
ImplicitSystem & operator=(const ImplicitSystem &)=delete
void setup_static_condensation_preconditioner(T &solver)
Sets up the static condensation preconditioner for the supplied solver.
virtual void release_linear_solver(LinearSolver< Number > *) const
Currently a no-op.
virtual void forward_qoi_parameter_sensitivity(const QoISet &qoi_indices, const ParameterVector ¶meters, SensitivityData &sensitivities) override
Solves for the derivative of each of the system's quantities of interest q in qoi[qoi_indices] with r...
ImplicitSystem sys_type
The type of system.
unsigned int number() const
virtual void assembly(bool, bool, bool=false, bool=false)
Assembles a residual in rhs and/or a jacobian in matrix, as requested.
Data structure for holding completed parameter sensitivity calculations.
virtual void adjoint_qoi_parameter_sensitivity(const QoISet &qoi_indices, const ParameterVector ¶meters, SensitivityData &sensitivities) override
Solves for the derivative of each of the system's quantities of interest q in qoi[qoi_indices] with r...
StaticCondensation & get_static_condensation()
virtual void disable_cache() override
Avoids use of any cached data that might affect any solve result.
virtual void solve() override
For explicit systems, just assemble and solve the system A*x=b.
ExplicitSystem Parent
The type of the parent.
ImplicitSystem(EquationSystems &es, const std::string &name, const unsigned int number)
Constructor.
virtual void add_matrices() override
Adds the system matrix.
SparseMatrix< Number > * matrix
The system matrix.
virtual void assemble_residual_derivatives(const ParameterVector ¶meters) override
Residual parameter derivative function.
virtual void assemble() override
Prepares matrix and rhs for system assembly, then calls user assembly function.
const std::string & name() const
std::unique_ptr< LinearSolver< Number > > linear_solver
This class handles all the details of interfacing with various linear algebra packages like PETSc or ...
virtual void qoi_parameter_hessian(const QoISet &qoi_indices, const ParameterVector ¶meters, SensitivityData &hessian) override
For each of the system's quantities of interest q in qoi[qoi_indices], and for a vector of parameters...
virtual std::pair< unsigned int, Real > weighted_sensitivity_adjoint_solve(const ParameterVector ¶meters, const ParameterVector &weights, const QoISet &qoi_indices=QoISet()) override
Assembles & solves the linear system(s) (dR/du)^T*z_w = sum(w_p*(d^2q/dudp - d^2R/dudp*z)), for those parameters p contained within parameters, weighted by the values w_p found within weights.
const SparseMatrix< Number > & get_system_matrix() const
Manages consistently variables, degrees of freedom, and coefficient vectors for explicit systems...
Manages consistently variables, degrees of freedom, coefficient vectors, and matrices for implicit sy...