implicit_system.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2014 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 
36 
37 
46 // ------------------------------------------------------------
47 // ImplicitSystem class definition
48 
50 {
51 public:
52 
58  const std::string& name,
59  const unsigned int number);
60 
64  virtual ~ImplicitSystem ();
65 
70 
74  sys_type & system () { return *this; }
75 
80 
85  virtual void clear ();
86 
91  virtual void reinit ();
92 
98  virtual void assemble ();
99 
100 // /**
101 // * Assembles & solves the linear system Ax=b.
102 // */
103 // virtual void solve ();
104 
105 
110  virtual void disable_cache ();
111 
116  virtual std::string system_type () const { return "Implicit"; }
117 
126  virtual LinearSolver<Number> *get_linear_solver() const;
127 
133  virtual std::pair<unsigned int, Real>
135 
140  virtual void release_linear_solver(LinearSolver<Number> *) const;
141 
149  virtual void assembly(bool /* get_residual */ , bool /* get_jacobian */)
150  { libmesh_error(); }
151 
163  virtual void assemble_residual_derivatives (const ParameterVector& parameters);
164 
172  virtual std::pair<unsigned int, Real>
173  sensitivity_solve (const ParameterVector& parameters);
174 
183  virtual std::pair<unsigned int, Real>
184  weighted_sensitivity_solve (const ParameterVector& parameters,
185  const ParameterVector& weights);
186 
196  virtual std::pair<unsigned int, Real>
197  adjoint_solve (const QoISet& qoi_indices = QoISet());
198 
211  virtual std::pair<unsigned int, Real>
213  const ParameterVector& weights,
214  const QoISet& qoi_indices = QoISet());
215 
227  virtual void adjoint_qoi_parameter_sensitivity (const QoISet& qoi_indices,
228  const ParameterVector& parameters,
229  SensitivityData& sensitivities);
230 
242  virtual void forward_qoi_parameter_sensitivity (const QoISet& qoi_indices,
243  const ParameterVector& parameters,
244  SensitivityData& sensitivities);
245 
254  virtual void qoi_parameter_hessian(const QoISet& qoi_indices,
255  const ParameterVector& parameters,
256  SensitivityData& hessian);
257 
268  virtual void qoi_parameter_hessian_vector_product(const QoISet& qoi_indices,
269  const ParameterVector& parameters,
270  const ParameterVector& vector,
271  SensitivityData& product);
272 
276  typedef std::map<std::string, SparseMatrix<Number>* >::iterator matrices_iterator;
277  typedef std::map<std::string, SparseMatrix<Number>* >::const_iterator const_matrices_iterator;
278 
287  SparseMatrix<Number> & add_matrix (const std::string& mat_name);
288 
293  bool have_matrix (const std::string& mat_name) const;
294 
300  const SparseMatrix<Number> * request_matrix (const std::string& mat_name) const;
301 
307  SparseMatrix<Number> * request_matrix (const std::string& mat_name);
308 
315  const SparseMatrix<Number> & get_matrix (const std::string& mat_name) const;
316 
323  SparseMatrix<Number> & get_matrix (const std::string& mat_name);
324 
328  virtual unsigned int n_matrices () const;
329 
336 
337 
338 
339 protected:
340 
345  virtual void init_data ();
346 
350  virtual void init_matrices ();
351 
352 
353 
354 private:
355 
360  void add_system_matrix ();
361 
365  std::map<std::string, SparseMatrix<Number>* > _matrices;
366 
371 };
372 
373 
374 
375 // ------------------------------------------------------------
376 // ImplicitSystem inline methods
377 inline
378 bool ImplicitSystem::have_matrix (const std::string& mat_name) const
379 {
380  return (_matrices.count(mat_name));
381 }
382 
383 
384 inline
385 unsigned int ImplicitSystem::n_matrices () const
386 {
387  return libmesh_cast_int<unsigned int>(_matrices.size());
388 }
389 
390 } // namespace libMesh
391 
392 #endif // LIBMESH_IMPLICIT_SYSTEM_H

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:05 UTC

Hosted By:
SourceForge.net Logo