MAST
level_set_elem_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 
21 #ifndef __mast__level_set_elem_base__
22 #define __mast__level_set_elem_base__
23 
24 // MAST includes
25 #include "base/elem_base.h"
26 
27 namespace MAST {
28 
29  // Forward declerations
30  class ElementPropertyCardBase;
31  class BoundaryConditionBase;
32  class FEMOperatorMatrix;
33  template <typename ValType> class FieldFunction;
34 
35 
37  public MAST::ElementBase
38  {
39  public:
45  const MAST::GeomElem& elem);
46 
47 
48  virtual ~LevelSetElementBase();
49 
50 
51  void
53 
54  libmesh_assert(!_phi_vel);
55  _phi_vel = &vel;
56  }
57 
63  void set_propagation_mode(bool f) {
64  _if_propagation = f;
65  }
66 
67 
74 
75 
79  virtual bool
80  internal_residual (bool request_jacobian,
81  RealVectorX& f,
82  RealMatrixX& jac);
83 
84 
88  virtual bool
89  velocity_residual (bool request_jacobian,
90  RealVectorX& f,
91  RealMatrixX& jac_xdot,
92  RealMatrixX& jac);
93 
97  bool
98  side_external_residual (bool request_jacobian,
99  RealVectorX& f,
100  RealMatrixX& jac,
101  std::multimap<libMesh::boundary_id_type, MAST::BoundaryConditionBase*>& bc);
102 
106  bool
107  volume_external_residual (bool request_jacobian,
108  RealVectorX& f,
109  RealMatrixX& jac,
110  std::multimap<libMesh::subdomain_id_type, MAST::BoundaryConditionBase*>& bc);
111 
115  virtual bool
116  internal_residual_sensitivity (bool request_jacobian,
117  RealVectorX& f,
118  RealMatrixX& jac);
122  virtual bool
123  velocity_residual_sensitivity (bool request_jacobian,
124  RealVectorX& f,
125  RealMatrixX& jac);
126 
130  bool
131  side_external_residual_sensitivity (bool request_jacobian,
132  RealVectorX& f,
133  RealMatrixX& jac,
134  std::multimap<libMesh::boundary_id_type, MAST::BoundaryConditionBase*>& bc);
135 
139  bool
140  volume_external_residual_sensitivity (bool request_jacobian,
141  RealVectorX& f,
142  RealMatrixX& jac,
143  std::multimap<libMesh::subdomain_id_type, MAST::BoundaryConditionBase*>& bc);
144 
149  Real
150  volume();
151 
165  Real
166  perimeter();
167 
173  Real
175 
176 
182  Real
183  volume_boundary_velocity_on_side(unsigned int s);
184 
185  protected:
186 
190  void _velocity_and_source(const unsigned int qp,
191  const libMesh::Point& p,
192  const Real t,
193  const MAST::FEMOperatorMatrix& Bmat,
194  const std::vector<MAST::FEMOperatorMatrix>& dBmat,
195  RealVectorX& vel,
196  Real& source);
197 
201  void _tau(const MAST::FEBase& fe,
202  unsigned int qp,
203  const MAST::FEMOperatorMatrix& Bmat,
204  const std::vector<MAST::FEMOperatorMatrix>& dBmat,
205  const RealVectorX& vel,
206  RealMatrixX& tau);
207 
208 
209 
210 
211  void
212  _dc_operator(const MAST::FEBase& fe,
213  const unsigned int qp,
214  const std::vector<MAST::FEMOperatorMatrix>& dB_mat,
215  const RealVectorX& vel,
216  Real& dc);
217 
218  void
219  _calculate_dxidX (const MAST::FEBase& fe,
220  const unsigned int qp,
221  RealMatrixX& dxi_dX);
222 
235  void _initialize_fem_operators(const unsigned int qp,
236  const MAST::FEBase& fe,
238  std::vector<MAST::FEMOperatorMatrix>& dBmat);
239 
240 
245 
251 
256  };
257 
258 }
259 
260 
261 #endif // __mast__level_set_elem_base__
262 
virtual bool internal_residual(bool request_jacobian, RealVectorX &f, RealMatrixX &jac)
internal force contribution to system residual
void _calculate_dxidX(const MAST::FEBase &fe, const unsigned int qp, RealMatrixX &dxi_dX)
const RealVectorX & sol(bool if_sens=false) const
Definition: elem_base.cpp:53
const MAST::FieldFunction< Real > * _phi_vel
element property
virtual bool velocity_residual(bool request_jacobian, RealVectorX &f, RealMatrixX &jac_xdot, RealMatrixX &jac)
inertial force contribution to system residual
void _tau(const MAST::FEBase &fe, unsigned int qp, const MAST::FEMOperatorMatrix &Bmat, const std::vector< MAST::FEMOperatorMatrix > &dBmat, const RealVectorX &vel, RealMatrixX &tau)
initializes the tau operator
const MAST::GeomElem & elem() const
Definition: elem_base.h:117
Real perimeter()
Approximates the integral of the Dirac delta function to approximate the perimeter.
Real perimeter_sensitivity()
computes the partial derivative of the integral of the Dirac delta function using the solution and se...
libMesh::Real Real
RealVectorX _ref_sol
reference solution for reinitialization of the level set
Real volume_boundary_velocity_on_side(unsigned int s)
MAST::AssemblyBase & assembly()
Definition: elem_base.h:103
void set_reference_solution_for_initialization(const RealVectorX &sol)
For reinitialization to , the solution before initialization is used to calculate the source and velo...
bool volume_external_residual_sensitivity(bool request_jacobian, RealVectorX &f, RealMatrixX &jac, std::multimap< libMesh::subdomain_id_type, MAST::BoundaryConditionBase * > &bc)
sensitivity of the volume external force contribution to system residual
void set_propagation_mode(bool f)
This can operate in one of two modes: propagation of level set given Vn, or reinitialization of level...
Matrix< Real, Dynamic, Dynamic > RealMatrixX
bool side_external_residual_sensitivity(bool request_jacobian, RealVectorX &f, RealMatrixX &jac, std::multimap< libMesh::boundary_id_type, MAST::BoundaryConditionBase * > &bc)
sensitivity of the side external force contribution to system residual
virtual bool velocity_residual_sensitivity(bool request_jacobian, RealVectorX &f, RealMatrixX &jac)
sensitivity of the damping force contribution to system residual
Matrix< Real, Dynamic, 1 > RealVectorX
void set_velocity_function(const MAST::FieldFunction< Real > &vel)
This class acts as a wrapper around libMesh::Elem for the purpose of providing a uniform interface fo...
Definition: geom_elem.h:59
bool _if_propagation
this can operate in one of two modes: propagation of level set given Vn, or reinitialization of level...
void _initialize_fem_operators(const unsigned int qp, const MAST::FEBase &fe, MAST::FEMOperatorMatrix &Bmat, std::vector< MAST::FEMOperatorMatrix > &dBmat)
When mass = false, initializes the FEM operator matrix to the shape functions as ...
bool volume_external_residual(bool request_jacobian, RealVectorX &f, RealMatrixX &jac, std::multimap< libMesh::subdomain_id_type, MAST::BoundaryConditionBase * > &bc)
volume external force contribution to system residual
void _velocity_and_source(const unsigned int qp, const libMesh::Point &p, const Real t, const MAST::FEMOperatorMatrix &Bmat, const std::vector< MAST::FEMOperatorMatrix > &dBmat, RealVectorX &vel, Real &source)
calculates the velocity at the quadrature point
LevelSetElementBase(MAST::SystemInitialization &sys, MAST::AssemblyBase &assembly, const MAST::GeomElem &elem)
Constructor.
virtual bool internal_residual_sensitivity(bool request_jacobian, RealVectorX &f, RealMatrixX &jac)
sensitivity of the internal force contribution to system residual
bool side_external_residual(bool request_jacobian, RealVectorX &f, RealMatrixX &jac, std::multimap< libMesh::boundary_id_type, MAST::BoundaryConditionBase * > &bc)
side external force contribution to system residual
void _dc_operator(const MAST::FEBase &fe, const unsigned int qp, const std::vector< MAST::FEMOperatorMatrix > &dB_mat, const RealVectorX &vel, Real &dc)
This is the base class for elements that implement calculation of finite element quantities over the ...
Definition: elem_base.h:72