implicit_system.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #ifndef LIBMESH_IMPLICIT_SYSTEM_H
21 #define LIBMESH_IMPLICIT_SYSTEM_H
22 
23 // Local Includes
25 
26 // C++ includes
27 #include <cstddef>
28 
29 namespace libMesh
30 {
31 
32 // Forward declarations
33 template <typename T> class LinearSolver;
34 template <typename T> class SparseMatrix;
35 
49 {
50 public:
51 
57  const std::string & name,
58  const unsigned int number);
59 
63  virtual ~ImplicitSystem ();
64 
69 
73  sys_type & system () { return *this; }
74 
79 
84  virtual void clear () libmesh_override;
85 
90  virtual void reinit () libmesh_override;
91 
97  virtual void assemble () libmesh_override;
98 
103  virtual void disable_cache () libmesh_override;
104 
109  virtual std::string system_type () const libmesh_override { return "Implicit"; }
110 
125  virtual LinearSolver<Number> * get_linear_solver() const;
126 
132  virtual std::pair<unsigned int, Real>
134 
143  virtual void release_linear_solver(LinearSolver<Number> *) const;
144 
152  virtual void assembly (bool /* get_residual */,
153  bool /* get_jacobian */,
154  bool /* apply_heterogeneous_constraints */ = false,
155  bool /* apply_no_constraints */ = false)
156  { libmesh_not_implemented(); }
157 
169  virtual void assemble_residual_derivatives (const ParameterVector & parameters) libmesh_override;
170 
178  virtual std::pair<unsigned int, Real>
179  sensitivity_solve (const ParameterVector & parameters) libmesh_override;
180 
189  virtual std::pair<unsigned int, Real>
190  weighted_sensitivity_solve (const ParameterVector & parameters,
191  const ParameterVector & weights) libmesh_override;
192 
202  virtual std::pair<unsigned int, Real>
203  adjoint_solve (const QoISet & qoi_indices = QoISet()) libmesh_override;
204 
217  virtual std::pair<unsigned int, Real>
219  const ParameterVector & weights,
220  const QoISet & qoi_indices = QoISet()) libmesh_override;
221 
233  virtual void adjoint_qoi_parameter_sensitivity (const QoISet & qoi_indices,
234  const ParameterVector & parameters,
235  SensitivityData & sensitivities) libmesh_override;
236 
248  virtual void forward_qoi_parameter_sensitivity (const QoISet & qoi_indices,
249  const ParameterVector & parameters,
250  SensitivityData & sensitivities) libmesh_override;
251 
260  virtual void qoi_parameter_hessian(const QoISet & qoi_indices,
261  const ParameterVector & parameters,
262  SensitivityData & hessian) libmesh_override;
263 
274  virtual void qoi_parameter_hessian_vector_product(const QoISet & qoi_indices,
275  const ParameterVector & parameters,
276  const ParameterVector & vector,
277  SensitivityData & product) libmesh_override;
278 
282  typedef std::map<std::string, SparseMatrix<Number> *>::iterator matrices_iterator;
283  typedef std::map<std::string, SparseMatrix<Number> *>::const_iterator const_matrices_iterator;
284 
293  SparseMatrix<Number> & add_matrix (const std::string & mat_name);
294 
298  void remove_matrix(const std::string & mat_name);
299 
304  bool have_matrix (const std::string & mat_name) const;
305 
311  const SparseMatrix<Number> * request_matrix (const std::string & mat_name) const;
312 
318  SparseMatrix<Number> * request_matrix (const std::string & mat_name);
319 
328  const SparseMatrix<Number> & get_matrix (const std::string & mat_name) const;
329 
338  SparseMatrix<Number> & get_matrix (const std::string & mat_name);
339 
343  virtual unsigned int n_matrices () const libmesh_override;
344 
351 
358 
359 protected:
360 
365  virtual void init_data () libmesh_override;
366 
370  virtual void init_matrices ();
371 
372 private:
373 
378  void add_system_matrix ();
379 
383  std::map<std::string, SparseMatrix<Number> *> _matrices;
384 
389 };
390 
391 
392 
393 // ------------------------------------------------------------
394 // ImplicitSystem inline methods
395 inline
396 bool ImplicitSystem::have_matrix (const std::string & mat_name) const
397 {
398  return (_matrices.count(mat_name));
399 }
400 
401 
402 inline
403 unsigned int ImplicitSystem::n_matrices () const
404 {
405  return cast_int<unsigned int>(_matrices.size());
406 }
407 
408 } // namespace libMesh
409 
410 #endif // LIBMESH_IMPLICIT_SYSTEM_H
virtual LinearSolver< Number > * get_linear_solver() const
Manages multiples systems of equations.
virtual std::pair< unsigned int, Real > adjoint_solve(const QoISet &qoi_indices=QoISet()) libmesh_override
Specifies parameters for parameter sensitivity calculations.
virtual void disable_cache() libmesh_override
Used to specify quantities of interest in a simulation.
Definition: qoi_set.h:45
virtual void adjoint_qoi_parameter_sensitivity(const QoISet &qoi_indices, const ParameterVector &parameters, SensitivityData &sensitivities) libmesh_override
virtual void qoi_parameter_hessian_vector_product(const QoISet &qoi_indices, const ParameterVector &parameters, const ParameterVector &vector, SensitivityData &product) libmesh_override
SparseMatrix< Number > & add_matrix(const std::string &mat_name)
virtual std::pair< unsigned int, Real > sensitivity_solve(const ParameterVector &parameters) libmesh_override
std::map< std::string, SparseMatrix< Number > * > _matrices
virtual std::string system_type() const libmesh_override
const std::string & name() const
Definition: system.h:1996
virtual std::pair< unsigned int, Real > get_linear_solve_parameters() const
virtual void assembly(bool, bool, bool=false, bool=false)
Holds completed parameter sensitivity calculations.
virtual void qoi_parameter_hessian(const QoISet &qoi_indices, const ParameterVector &parameters, SensitivityData &hessian) libmesh_override
bool have_matrix(const std::string &mat_name) const
virtual void reinit() libmesh_override
const SparseMatrix< Number > & get_matrix(const std::string &mat_name) const
const SparseMatrix< Number > * request_matrix(const std::string &mat_name) const
ImplicitSystem(EquationSystems &es, const std::string &name, const unsigned int number)
virtual void init_data() libmesh_override
std::map< std::string, SparseMatrix< Number > * >::const_iterator const_matrices_iterator
unsigned int number() const
Definition: system.h:2004
void remove_matrix(const std::string &mat_name)
virtual std::pair< unsigned int, Real > weighted_sensitivity_adjoint_solve(const ParameterVector &parameters, const ParameterVector &weights, const QoISet &qoi_indices=QoISet()) libmesh_override
SparseMatrix< Number > * matrix
std::map< std::string, SparseMatrix< Number > * >::iterator matrices_iterator
virtual void release_linear_solver(LinearSolver< Number > *) const
virtual std::pair< unsigned int, Real > weighted_sensitivity_solve(const ParameterVector &parameters, const ParameterVector &weights) libmesh_override
virtual void forward_qoi_parameter_sensitivity(const QoISet &qoi_indices, const ParameterVector &parameters, SensitivityData &sensitivities) libmesh_override
virtual void assemble() libmesh_override
virtual void assemble_residual_derivatives(const ParameterVector &parameters) libmesh_override
virtual unsigned int n_matrices() const libmesh_override
virtual void init_matrices()
Used for solving explicit systems of equations.
virtual void clear() libmesh_override
Used for solving implicit systems of equations.