MAST
assembly_base.h
Go to the documentation of this file.
1 /*
2  * MAST: Multidisciplinary-design Adaptation and Sensitivity Toolkit
3  * Copyright (C) 2013-2019 Manav Bhatia
4  *
5  * This library is free software; you can redistribute it and/or
6  * modify it under the terms of the GNU Lesser General Public
7  * License as published by the Free Software Foundation; either
8  * version 2.1 of the License, or (at your option) any later version.
9  *
10  * This library is distributed in the hope that it will be useful,
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13  * Lesser General Public License for more details.
14  *
15  * You should have received a copy of the GNU Lesser General Public
16  * License along with this library; if not, write to the Free Software
17  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
18  */
19 
20 #ifndef __mast__system_assembly_base_h__
21 #define __mast__system_assembly_base_h__
22 
23 // C++ includes
24 #include <map>
25 #include <memory>
26 
27 
28 // MAST includes
29 #include "base/mast_data_types.h"
30 
31 
32 // libMesh includes
33 #include "libmesh/system.h"
34 #include "libmesh/nonlinear_implicit_system.h"
35 #include "libmesh/sparse_matrix.h"
36 
37 
38 namespace MAST {
39 
40  // Forward declerations
41  class PhysicsDisciplineBase;
42  class SystemInitialization;
43  class ElementBase;
44  class MeshFieldFunction;
45  class NonlinearSystem;
46  class FEBase;
47  class AssemblyElemOperations;
48  class OutputAssemblyElemOperations;
49  class FunctionBase;
50 
51  class AssemblyBase:
52  public libMesh::NonlinearImplicitSystem::ComputeResidualandJacobian {
53  public:
54 
59  AssemblyBase();
60 
64  virtual ~AssemblyBase();
65 
66 
67  class SolverMonitor {
68  public:
70  virtual ~SolverMonitor(){}
71  virtual void init(MAST::AssemblyBase& assembly) = 0;
72  virtual void clear() = 0;
73  };
74 
87  public:
88  ElemParameterDependence(bool o_flag): override_flag(o_flag) {}
90  virtual bool if_elem_depends_on_parameter(const libMesh::Elem& e,
91  const MAST::FunctionBase& p) const = 0;
98  };
99 
104 
110 
116 
117 
122 
123 
127  virtual void
130 
137  void
139 
140 
141  void
143 
144 
148  virtual void
150 
151 
156  virtual void
158 
159 
164  virtual void
166 
167 
172  const MAST::NonlinearSystem& system() const;
173 
179 
184 
185 
191 
197 
201  void clear_solver_monitor();
202 
208 
209 
214 
219  virtual void
220  residual_and_jacobian (const libMesh::NumericVector<Real>& X,
221  libMesh::NumericVector<Real>* R,
222  libMesh::SparseMatrix<Real>* J,
223  libMesh::NonlinearImplicitSystem& S) {
224  libmesh_assert(false); // implement in the derived class
225  }
226 
238  virtual bool
240  libMesh::NumericVector<Real>& sensitivity_rhs) {
241  libmesh_assert(false); // implemented in the derived class
242  }
243 
247  virtual void
248  calculate_output(const libMesh::NumericVector<Real>& X,
250 
251 
255  virtual void
256  calculate_output_derivative(const libMesh::NumericVector<Real>& X,
258  libMesh::NumericVector<Real>& dq_dX);
259 
260 
270  virtual void
271  calculate_output_direct_sensitivity(const libMesh::NumericVector<Real>& X,
272  const libMesh::NumericVector<Real>* dXdp,
273  const MAST::FunctionBase& p,
275 
276 
282  virtual Real
283  calculate_output_adjoint_sensitivity(const libMesh::NumericVector<Real>& X,
284  const libMesh::NumericVector<Real>& dq_dX,
285  const MAST::FunctionBase& p,
288  const bool include_partial_sens = true);
289 
290 
296  std::unique_ptr<libMesh::NumericVector<Real> >
297  build_localized_vector(const libMesh::System& sys,
298  const libMesh::NumericVector<Real>& global) const;
299 
300 
301  protected:
302 
307 
308 
313 
314 
319 
324 
330 
337  };
338 
339 }
340 
341 #endif // __mast__system_assembly_base_h__
342 
MAST::AssemblyElemOperations * _elem_ops
provides assembly elem operations for use by this class
AssemblyBase()
constructor takes a reference to the discipline that provides the boundary conditions, volume loads, properties, etc.
void clear_solver_monitor()
clears the monitor object
This class implements a system for solution of nonlinear systems.
bool override_flag
if true, assume zero solution sensitivity when elem does not dependent on parameter.
Definition: assembly_base.h:97
virtual void clear_elem_operation_object()
clears the association of this object with the assembly element operation object. ...
virtual bool sensitivity_assemble(const MAST::FunctionBase &f, libMesh::NumericVector< Real > &sensitivity_rhs)
Assembly function.
This provides a wrapper FieldFunction compatible class that interpolates the solution using libMesh&#39;s...
void attach_solution_function(MAST::MeshFieldFunction &f)
tells the assembly object that this function is will need to be initialized before each residual eval...
virtual void calculate_output_derivative(const libMesh::NumericVector< Real > &X, MAST::OutputAssemblyElemOperations &output, libMesh::NumericVector< Real > &dq_dX)
calculates
void attach_elem_parameter_dependence_object(MAST::AssemblyBase::ElemParameterDependence &dep)
This object, if provided by user, will be used to reduce unnecessary computations in sensitivity anal...
virtual Real calculate_output_adjoint_sensitivity(const libMesh::NumericVector< Real > &X, const libMesh::NumericVector< Real > &dq_dX, const MAST::FunctionBase &p, MAST::AssemblyElemOperations &elem_ops, MAST::OutputAssemblyElemOperations &output, const bool include_partial_sens=true)
Evaluates the total sensitivity of output wrt p using the adjoint solution provided in dq_dX for a li...
void set_solver_monitor(MAST::AssemblyBase::SolverMonitor &monitor)
attaches the solver monitor, which is a user provided routine that is called each time ...
This provides the base class for definitin of element level contribution of output quantity in an ana...
bool close_matrix
flag to control the closing fo the Jacobian after assembly
virtual void set_elem_operation_object(MAST::AssemblyElemOperations &elem_ops)
attaches a element operation to this object, and associated this with the element operation object...
libMesh::Real Real
MAST::SystemInitialization * _system
System for which this assembly is performed.
virtual void init(MAST::AssemblyBase &assembly)=0
void detach_solution_function()
removes the attachment of the solution function
MAST::PhysicsDisciplineBase * _discipline
PhysicsDisciplineBase object for which this class is assembling.
virtual void clear_discipline_and_system()
clears association with a system to this discipline
virtual void calculate_output(const libMesh::NumericVector< Real > &X, MAST::OutputAssemblyElemOperations &output)
calculates the value of quantity .
virtual void calculate_output_direct_sensitivity(const libMesh::NumericVector< Real > &X, const libMesh::NumericVector< Real > *dXdp, const MAST::FunctionBase &p, MAST::OutputAssemblyElemOperations &output)
evaluates the sensitivity of the outputs in the attached discipline with respect to the parametrs in ...
void clear_elem_parameter_dependence_object()
std::unique_ptr< libMesh::NumericVector< Real > > build_localized_vector(const libMesh::System &sys, const libMesh::NumericVector< Real > &global) const
localizes the parallel vector so that the local copy stores all values necessary for calculation of t...
MAST::AssemblyBase::SolverMonitor * get_solver_monitor()
Inherited objects from this class can be provided by the user provide assessment of whether or not an...
Definition: assembly_base.h:86
MAST::SystemInitialization & system_init()
MAST::AssemblyBase::ElemParameterDependence * _param_dependence
If provided by user, this object is used by sensitiivty analysis to check for whether or the current ...
virtual void set_discipline_and_system(MAST::PhysicsDisciplineBase &discipline, MAST::SystemInitialization &system)
attaches a system to this discipline
virtual ~AssemblyBase()
virtual destructor
MAST::AssemblyElemOperations & get_elem_ops()
virtual void residual_and_jacobian(const libMesh::NumericVector< Real > &X, libMesh::NumericVector< Real > *R, libMesh::SparseMatrix< Real > *J, libMesh::NonlinearImplicitSystem &S)
function that assembles the matrices and vectors quantities for nonlinear solution ...
MAST::AssemblyBase::SolverMonitor * _solver_monitor
User provided solver monitor is attached to the linear nonlinear solvers, if provided.
const MAST::NonlinearSystem & system() const
MAST::MeshFieldFunction * _sol_function
system solution that will be initialized before each solution
const MAST::PhysicsDisciplineBase & discipline() const