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 
48 {
49 public:
50 
56  const std::string & name,
57  const unsigned int number);
58 
62  virtual ~ImplicitSystem ();
63 
68 
72  sys_type & system () { return *this; }
73 
78 
83  virtual void clear () libmesh_override;
84 
89  virtual void reinit () libmesh_override;
90 
96  virtual void assemble () libmesh_override;
97 
102  virtual void disable_cache () libmesh_override;
103 
108  virtual std::string system_type () const libmesh_override { return "Implicit"; }
109 
124  virtual LinearSolver<Number> * get_linear_solver() const;
125 
131  virtual std::pair<unsigned int, Real>
133 
142  virtual void release_linear_solver(LinearSolver<Number> *) const;
143 
151  virtual void assembly (bool /* get_residual */,
152  bool /* get_jacobian */,
153  bool /* apply_heterogeneous_constraints */ = false)
154  { libmesh_not_implemented(); }
155 
167  virtual void assemble_residual_derivatives (const ParameterVector & parameters) libmesh_override;
168 
176  virtual std::pair<unsigned int, Real>
177  sensitivity_solve (const ParameterVector & parameters) libmesh_override;
178 
187  virtual std::pair<unsigned int, Real>
188  weighted_sensitivity_solve (const ParameterVector & parameters,
189  const ParameterVector & weights) libmesh_override;
190 
200  virtual std::pair<unsigned int, Real>
201  adjoint_solve (const QoISet & qoi_indices = QoISet()) libmesh_override;
202 
215  virtual std::pair<unsigned int, Real>
217  const ParameterVector & weights,
218  const QoISet & qoi_indices = QoISet()) libmesh_override;
219 
231  virtual void adjoint_qoi_parameter_sensitivity (const QoISet & qoi_indices,
232  const ParameterVector & parameters,
233  SensitivityData & sensitivities) libmesh_override;
234 
246  virtual void forward_qoi_parameter_sensitivity (const QoISet & qoi_indices,
247  const ParameterVector & parameters,
248  SensitivityData & sensitivities) libmesh_override;
249 
258  virtual void qoi_parameter_hessian(const QoISet & qoi_indices,
259  const ParameterVector & parameters,
260  SensitivityData & hessian) libmesh_override;
261 
272  virtual void qoi_parameter_hessian_vector_product(const QoISet & qoi_indices,
273  const ParameterVector & parameters,
274  const ParameterVector & vector,
275  SensitivityData & product) libmesh_override;
276 
280  typedef std::map<std::string, SparseMatrix<Number> *>::iterator matrices_iterator;
281  typedef std::map<std::string, SparseMatrix<Number> *>::const_iterator const_matrices_iterator;
282 
291  SparseMatrix<Number> & add_matrix (const std::string & mat_name);
292 
297  bool have_matrix (const std::string & mat_name) const;
298 
304  const SparseMatrix<Number> * request_matrix (const std::string & mat_name) const;
305 
311  SparseMatrix<Number> * request_matrix (const std::string & mat_name);
312 
319  const SparseMatrix<Number> & get_matrix (const std::string & mat_name) const;
320 
327  SparseMatrix<Number> & get_matrix (const std::string & mat_name);
328 
332  virtual unsigned int n_matrices () const libmesh_override;
333 
340 
347 
348 
349 protected:
350 
355  virtual void init_data () libmesh_override;
356 
360  virtual void init_matrices ();
361 
362 
363 
364 private:
365 
370  void add_system_matrix ();
371 
375  std::map<std::string, SparseMatrix<Number> *> _matrices;
376 
381 };
382 
383 
384 
385 // ------------------------------------------------------------
386 // ImplicitSystem inline methods
387 inline
388 bool ImplicitSystem::have_matrix (const std::string & mat_name) const
389 {
390  return (_matrices.count(mat_name));
391 }
392 
393 
394 inline
395 unsigned int ImplicitSystem::n_matrices () const
396 {
397  return cast_int<unsigned int>(_matrices.size());
398 }
399 
400 } // namespace libMesh
401 
402 #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:1987
virtual std::pair< unsigned int, Real > get_linear_solve_parameters() const
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
virtual void assembly(bool, bool, bool=false)
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:1995
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.