MAST
Multidisciplinary-design Adaptation and Sensitivity Toolkit (MAST)
structural_nonlinear_assembly.cpp
Go to the documentation of this file.
1 /*
2  * MAST: Multidisciplinary-design Adaptation and Sensitivity Toolkit
3  * Copyright (C) 2013-2020 Manav Bhatia and MAST authors
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 // MAST includes
28 #include "mesh/geom_elem.h"
29 
30 
33 _incompatible_sol_assembly(nullptr) {
34 
35 }
36 
37 
38 
41 
42 }
43 
44 
45 
46 void
48 
49  unsigned int
50  n = (unsigned int)sol.size();
51 
53  zero = RealVectorX::Zero(n);
54 
56  _physics_elem->set_velocity (zero); // set to zero vector for a quasi-steady analysis
57  _physics_elem->set_acceleration(zero); // set to zero vector for a quasi-steady analysis
58 
59 
61 
62  // set the incompatible mode solution if required by the
63  // element
64 
66  dynamic_cast<MAST::StructuralElementBase&>(*_physics_elem);
67 
68  if (s_elem.if_incompatible_modes())
70  }
71 }
72 
73 
74 void
76 set_elem_data(unsigned int dim,
77  const libMesh::Elem& ref_elem,
78  MAST::GeomElem& elem) const {
79 
80  libmesh_assert(!_physics_elem);
81 
82  if (dim == 1) {
83 
85  dynamic_cast<const MAST::ElementPropertyCard1D&>(_discipline->get_property_card(ref_elem));
86 
87  elem.set_local_y_vector(p.y_vector());
88  }
89 }
90 
91 
92 
93 void
95 init(const MAST::GeomElem& elem) {
96 
97  libmesh_assert(!_physics_elem);
98  libmesh_assert(_system);
99  libmesh_assert(_assembly);
100 
102  dynamic_cast<const MAST::ElementPropertyCardBase&>
104 
105  _physics_elem =
106  MAST::build_structural_element(*_system, elem, p).release();
107 }
108 
109 
110 
111 
112 void
114 elem_calculations(bool if_jac,
115  RealVectorX& vec,
116  RealMatrixX& mat) {
117 
118  libmesh_assert(_physics_elem);
119 
121  dynamic_cast<MAST::StructuralElementBase&>(*_physics_elem);
122 
123  vec.setZero();
124  mat.setZero();
126  dummy = RealMatrixX::Zero(mat.rows(), mat.cols());
127 
128  e.internal_residual(if_jac, vec, mat);
129  e.side_external_residual(if_jac,
130  vec,
131  dummy,
132  mat,
134  e.volume_external_residual(if_jac,
135  vec,
136  dummy,
137  mat,
139 }
140 
141 
142 
143 void
146 
147  libmesh_assert(_physics_elem);
148 
150  dynamic_cast<MAST::StructuralElementBase&>(*_physics_elem);
151 
152  vec.setZero();
154  mat = RealMatrixX::Zero(vec.size(), vec.size());
155 
156  e.linearized_internal_residual(false, vec, mat);
158  vec,
159  mat,
160  mat,
163  vec,
164  mat,
165  mat,
167 }
168 
169 
170 
171 
172 void
175  RealVectorX& vec) {
176 
177  libmesh_assert(_physics_elem);
178 
180  dynamic_cast<MAST::StructuralElementBase&>(*_physics_elem);
181 
182  vec.setZero();
184  dummy = RealMatrixX::Zero(vec.size(), vec.size());
185 
186  e.internal_residual_sensitivity(f, false, vec, dummy);
188  vec,
189  dummy,
190  dummy,
193  vec,
194  dummy,
195  dummy,
197 }
198 
199 
200 
201 void
204  RealVectorX& vec) {
205 
206  libmesh_assert(_physics_elem);
207  libmesh_assert(f.is_topology_parameter());
208 
209  std::pair<const MAST::FieldFunction<RealVectorX>*, unsigned int>
210  val = this->get_elem_boundary_velocity_data();
211 
212  if (val.first) {
213 
215  dynamic_cast<MAST::StructuralElementBase&>(*_physics_elem);
216 
217  vec.setZero();
219  dummy = RealMatrixX::Zero(vec.size(), vec.size());
220 
222  val.second,
223  *val.first,
224  false,
225  vec,
226  dummy);
228  val.second,
229  *val.first,
231  false,
232  vec,
233  dummy);
234  /*e.side_external_residual_sensitivity(f, false,
235  vec,
236  dummy,
237  dummy,
238  _discipline->side_loads());*/
239  }
240 }
241 
242 
243 
244 void
248  RealVectorX& vec) {
249 
250  libmesh_assert(_physics_elem);
251  libmesh_assert(f.is_topology_parameter());
252 
254  &elem = dynamic_cast<const MAST::LevelSetIntersectedElem&>(_physics_elem->elem());
255 
256  // sensitivity only exists at the boundary. So, we proceed with calculation
257  // only if this element has an intersection in the interior, or with a side.
258  if (elem.if_elem_has_level_set_boundary() &&
259  elem.if_subelem_has_side_on_level_set_boundary()) {
260 
262  dynamic_cast<MAST::StructuralElementBase&>(*_physics_elem);
263 
264  vec.setZero();
266  dummy = RealMatrixX::Zero(vec.size(), vec.size());
267 
269  elem.get_subelem_side_on_level_set_boundary(),
270  vel,
271  false,
272  vec,
273  dummy);
275  elem.get_subelem_side_on_level_set_boundary(),
276  vel,
278  false,
279  vec,
280  dummy);
281  /*e.side_external_residual_sensitivity(f, false,
282  vec,
283  dummy,
284  dummy,
285  _discipline->side_loads());*/
286  }
287 }
288 
289 
290 
291 void
294 
295  libmesh_assert(_physics_elem);
296 
298  dynamic_cast<MAST::StructuralElementBase&>(*_physics_elem);
299  m.setZero();
300 
302 }
303 
304 
virtual void internal_residual_boundary_velocity(const MAST::FunctionBase &p, const unsigned int s, const MAST::FieldFunction< RealVectorX > &vel_f, bool request_jacobian, RealVectorX &f, RealMatrixX &jac)=0
calculates the term on side s: .
const MAST::ElementPropertyCardBase & get_property_card(const libMesh::Elem &elem) const
get property card for the specified element
virtual void set_velocity(const RealVectorX &vec, bool if_sens=false)
stores vec as velocity for element level calculations, or its sensitivity if if_sens is true...
Definition: elem_base.cpp:98
virtual void elem_sensitivity_calculations(const MAST::FunctionBase &f, RealVectorX &vec)
performs the element sensitivity calculations over elem, and returns the element residual sensitivity...
RealVectorX & y_vector()
returns value of the property val.
void set_local_y_vector(const RealVectorX &y_vec)
for 1D elements the transformed coordinate system attached to the element defines the local x-axis al...
Definition: geom_elem.cpp:119
StructuralNonlinearAssemblyElemOperations()
constructor associates this assembly object with the system
bool linearized_volume_external_residual(bool request_jacobian, RealVectorX &f, RealMatrixX &jac_xdot, RealMatrixX &jac, std::multimap< libMesh::subdomain_id_type, MAST::BoundaryConditionBase *> &bc)
volume external force contribution to system residual.
virtual bool internal_residual_sensitivity(const MAST::FunctionBase &p, bool request_jacobian, RealVectorX &f, RealMatrixX &jac)=0
sensitivity of the internal force contribution to system residual
MAST::PhysicsDisciplineBase * _discipline
virtual void set_elem_data(unsigned int dim, const libMesh::Elem &ref_elem, MAST::GeomElem &elem) const
sets the structural element y-vector if 1D element is used.
void set_elem_incompatible_sol(MAST::StructuralElementBase &elem)
virtual std::pair< const MAST::FieldFunction< RealVectorX > *, unsigned int > get_elem_boundary_velocity_data()
searches through the side load data and populates the data with the boundary id and velocity function...
virtual bool if_incompatible_modes() const =0
virtual bool is_topology_parameter() const
Definition: function_base.h:97
const MAST::SideBCMapType & side_loads() const
virtual void elem_linearized_jacobian_solution_product(RealVectorX &vec)
performs the element calculations over elem, and returns the element vector quantity in vec...
Matrix< Real, Dynamic, Dynamic > RealMatrixX
bool side_external_residual(bool request_jacobian, RealVectorX &f, RealMatrixX &jac_xdot, RealMatrixX &jac, std::multimap< libMesh::boundary_id_type, MAST::BoundaryConditionBase *> &bc)
side external force contribution to system residual.
This class inherits from MAST::GeomElem and provides an interface to initialize FE objects on sub-ele...
Matrix< Real, Dynamic, 1 > RealVectorX
bool volume_external_residual(bool request_jacobian, RealVectorX &f, RealMatrixX &jac_xdot, RealMatrixX &jac, std::multimap< libMesh::subdomain_id_type, MAST::BoundaryConditionBase *> &bc)
volume external force contribution to system residual.
std::unique_ptr< MAST::StructuralElementBase > build_structural_element(MAST::SystemInitialization &sys, const MAST::GeomElem &elem, const MAST::ElementPropertyCardBase &p)
builds the structural element for the specified element type
virtual ~StructuralNonlinearAssemblyElemOperations()
destructor resets the association of this assembly object with the system
This class acts as a wrapper around libMesh::Elem for the purpose of providing a uniform interface fo...
Definition: geom_elem.h:59
virtual void elem_second_derivative_dot_solution_assembly(RealMatrixX &mat)
calculates over elem, and returns the matrix in vec .
const MAST::VolumeBCMapType & volume_loads() const
virtual void set_elem_solution(const RealVectorX &sol)
sets the element solution(s) before calculations
const MAST::GeomElem & elem() const
Definition: elem_base.h:109
virtual void init(const MAST::GeomElem &elem)
initializes the object for the geometric element elem.
virtual void set_solution(const RealVectorX &vec, bool if_sens=false)
stores vec as solution for element level calculations, or its sensitivity if if_sens is true...
Definition: elem_base.cpp:60
virtual void elem_calculations(bool if_jac, RealVectorX &vec, RealMatrixX &mat)
performs the element calculations over elem, and returns the element vector and matrix quantities in ...
bool side_external_residual_sensitivity(const MAST::FunctionBase &p, bool request_jacobian, RealVectorX &f, RealMatrixX &jac_xdot, RealMatrixX &jac, std::multimap< libMesh::boundary_id_type, MAST::BoundaryConditionBase *> &bc)
sensitivity of the side external force contribution to system residual
bool linearized_side_external_residual(bool request_jacobian, RealVectorX &f, RealMatrixX &jac_xdot, RealMatrixX &jac, std::multimap< libMesh::boundary_id_type, MAST::BoundaryConditionBase *> &bc)
side external force contribution to system residual.
virtual bool internal_residual_jac_dot_state_sensitivity(RealMatrixX &jac)=0
calculates d[J]/d{x} .
virtual bool linearized_internal_residual(bool request_jacobian, RealVectorX &f, RealMatrixX &jac)
internal force contribution to system residual of the linearized problem
bool volume_external_residual_sensitivity(const MAST::FunctionBase &p, bool request_jacobian, RealVectorX &f, RealMatrixX &jac_xdot, RealMatrixX &jac, std::multimap< libMesh::subdomain_id_type, MAST::BoundaryConditionBase *> &bc)
sensitivity of the volume external force contribution to system residual
virtual bool internal_residual(bool request_jacobian, RealVectorX &f, RealMatrixX &jac)=0
internal force contribution to system residual
virtual void elem_topology_sensitivity_calculations(const MAST::FunctionBase &f, RealVectorX &vec)
performs the element topology sensitivity calculations over elem, and returns the element residual se...
void volume_external_residual_boundary_velocity(const MAST::FunctionBase &p, const unsigned int s, const MAST::FieldFunction< RealVectorX > &vel_f, std::multimap< libMesh::subdomain_id_type, MAST::BoundaryConditionBase *> &bc, bool request_jacobian, RealVectorX &f, RealMatrixX &jac)
boundary velocity contribution of volume external force.
virtual void set_acceleration(const RealVectorX &vec, bool if_sens=false)
stores vec as acceleration for element level calculations, or its sensitivity if if_sens is true...
Definition: elem_base.cpp:122
MAST::SystemInitialization * _system