20 #ifndef LIBMESH_DIFF_SYSTEM_H 21 #define LIBMESH_DIFF_SYSTEM_H 24 #include "libmesh/diff_context.h" 25 #include "libmesh/diff_physics.h" 26 #include "libmesh/diff_qoi.h" 27 #include "libmesh/implicit_system.h" 28 #include "libmesh/time_solver.h" 40 template <
typename T>
class NumericVector;
64 const std::string &
name,
65 const unsigned int number);
92 virtual void clear ()
override;
98 virtual void reinit ()
override;
123 virtual std::pair<unsigned int, Real>
136 virtual void assembly (
bool get_residual,
138 bool apply_heterogeneous_constraints =
false,
139 bool apply_no_constraints =
false)
override = 0;
146 virtual void solve ()
override;
152 virtual std::pair<unsigned int, Real>
160 libmesh_not_implemented();
162 return std::unique_ptr<DifferentiablePhysics>(
nullptr);
168 virtual std::unique_ptr<DifferentiableQoI>
clone()
override 170 libmesh_not_implemented();
172 return std::unique_ptr<DifferentiableQoI>(
nullptr);
414 #ifdef LIBMESH_ENABLE_DIRICHLET 439 std::stack<std::unique_ptr<DifferentiablePhysics>,
447 std::stack<std::unique_ptr<DifferentiableQoI>,
448 std::vector<std::unique_ptr<DifferentiableQoI>>>
_diff_qoi;
457 libmesh_assert_equal_to (&(
time_solver->system()),
this);
465 libmesh_assert_equal_to (&(
time_solver->system()),
this);
472 #endif // LIBMESH_DIFF_SYSTEM_H void attach_physics(DifferentiablePhysics *physics_in)
Attach external Physics object.
const DifferentiableQoI * get_qoi() const
This is the EquationSystems class.
virtual void assemble() override
Prepares matrix and rhs for matrix assembly.
bool _constrain_in_solver
_constrain_in_solver defaults to true; if false then we apply constraints only via residual terms in ...
This class provides all data required for a physics package (e.g.
virtual void assembly(bool get_residual, bool get_jacobian, bool apply_heterogeneous_constraints=false, bool apply_no_constraints=false) override=0
Assembles a residual in rhs and/or a jacobian in matrix, as requested.
This is a generic class that defines a solver to handle time integration of DifferentiableSystems.
virtual void set_constrain_in_solver(bool enable)
set_constrain_in_solver to false to apply constraints only via residual terms in the systems to be so...
Data structure for specifying which Quantities of Interest should be calculated in an adjoint or a pa...
bool print_jacobian_norms
Set print_jacobian_norms to true to print |J| whenever it is assembled.
std::unique_ptr< TimeSolver > time_solver
A pointer to the solver object we're going to use.
DifferentiableQoI * get_qoi()
void add_second_order_dot_vars()
Helper function to add "velocity" variables that are cousins to second order-in-time variables in the...
bool postprocess_sides
If postprocess_sides is true (it is false by default), the postprocessing loop will loop over all sid...
DifferentiablePhysics * get_physics()
void add_dot_var_dirichlet_bcs(unsigned int var_idx, unsigned int dot_var_idx)
Helper function to and Dirichlet boundary conditions to "dot" variable cousins of second order variab...
The libMesh namespace provides an interface to certain functionality in the library.
virtual std::unique_ptr< DiffContext > build_context()
Builds a DiffContext object with enough information to do evaluations on each element.
bool print_jacobians
Set print_jacobians to true to print J whenever it is assembled.
void pop_physics()
Pop a physics object off of our stack.
bool print_element_residuals
Set print_element_residuals to true to print each R_elem contribution.
virtual void reinit() override
Reinitializes the member data fields associated with the system, so that, e.g., assemble() may be use...
bool print_element_solutions
Set print_element_solutions to true to print each U_elem input.
This class provides a specific system class.
virtual bool get_constrain_in_solver()
bool print_solution_norms
Set print_residual_norms to true to print |U| whenever it is used in an assembly() call...
void swap_physics(DifferentiablePhysics *&swap_physics)
Swap current physics object with external object.
unsigned int number() const
unsigned int get_second_order_dot_var(unsigned int var) const
For a given second order (in time) variable var, this method will return the index to the correspondi...
bool print_residual_norms
Set print_residual_norms to true to print |F| whenever it is assembled.
virtual void side_postprocess(DiffContext &)
Does any work that needs to be done on side of elem in a postprocessing loop.
Real deltat
For time-dependent problems, this is the amount delta t to advance the solution in time...
std::stack< std::unique_ptr< DifferentiablePhysics >, std::vector< std::unique_ptr< DifferentiablePhysics > > > _diff_physics
Stack of pointers to objects to use for physics assembly evaluations.
std::stack< std::unique_ptr< DifferentiableQoI >, std::vector< std::unique_ptr< DifferentiableQoI > > > _diff_qoi
Pointer to object to use for quantity of interest assembly evaluations.
bool print_residuals
Set print_residuals to true to print F whenever it is assembled.
virtual std::unique_ptr< DifferentiablePhysics > clone_physics() override
We don't allow systems to be attached to each other.
DifferentiableSystem sys_type
The type of system.
void attach_qoi(DifferentiableQoI *qoi_in)
Attach external QoI object.
virtual void element_postprocess(DiffContext &)
Does any work that needs to be done on elem in a postprocessing loop.
bool have_second_order_scalar_vars() const
Check for any second order vars that are also belong to FEFamily::SCALAR.
bool print_solutions
Set print_solutions to true to print U whenever it is used in an assembly() call. ...
const DifferentiablePhysics * get_physics() const
This class provides a specific system class.
DifferentiableSystem(EquationSystems &es, const std::string &name, const unsigned int number)
Constructor.
This class provides a specific system class.
virtual void init_data() override
Initializes the member data fields associated with the system, so that, e.g., assemble() may be used...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual LinearSolver< Number > * get_linear_solver() const override
virtual std::pair< unsigned int, Real > get_linear_solve_parameters() const override
ImplicitSystem Parent
The type of the parent.
virtual std::unique_ptr< DifferentiableQoI > clone() override
We don't allow systems to be attached to each other.
bool print_element_jacobians
Set print_element_jacobians to true to print each J_elem contribution.
void set_time_solver(std::unique_ptr< TimeSolver > _time_solver)
Sets the time_solver FIXME: This code is a little dangerous as it transfers ownership from the TimeSo...
virtual std::unique_ptr< DifferentiablePhysics > clone_physics()=0
Copy of this object.
void push_physics(DifferentiablePhysics &new_physics)
Push a clone of a new physics object onto our stack, overriding the current physics until the new phy...
virtual void solve() override
Invokes the solver associated with the system.
virtual void clear() override
Clear all the data structures associated with the system.
const std::string & name() const
bool have_first_order_scalar_vars() const
Check for any first order vars that are also belong to FEFamily::SCALAR.
virtual void postprocess()
Executes a postprocessing loop over all elements, and if postprocess_sides is true over all sides...
virtual std::pair< unsigned int, Real > adjoint_solve(const QoISet &qoi_indices=QoISet()) override
This function sets the _is_adjoint boolean member of TimeSolver to true and then calls the adjoint_so...
virtual ~DifferentiableSystem()
DifferentiableSystem & operator=(const DifferentiableSystem &)=delete
TimeSolver & get_time_solver()
Manages consistently variables, degrees of freedom, coefficient vectors, and matrices for implicit sy...