MAST
compliance_output.cpp
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 // MAST includes
24 #include "base/assembly_base.h"
26 #include "base/nonlinear_system.h"
32 #include "mesh/geom_elem.h"
33 
34 
37 _compliance (0.),
38 _dcompliance_dp (0.) {
39 
40 }
41 
42 
43 
44 
46 
47 }
48 
49 
50 
51 
52 void
54 
55  _compliance = 0.;
56  _dcompliance_dp = 0.;
57 }
58 
59 
60 
61 void
63 
64  _compliance = 0.;
65  _dcompliance_dp = 0.;
66 }
67 
68 
69 void
71 
72  // make sure that this has not been initialized ana calculated for all elems
73  libmesh_assert(_physics_elem);
74 
76 
78  dynamic_cast<MAST::StructuralElementBase&>(*_physics_elem);
79 
81  vec = RealVectorX::Zero(e.sol().size());
82 
84  dummy = RealMatrixX::Zero(vec.size(), vec.size());
85 
86  e.side_external_residual(false,
87  vec,
88  dummy,
89  dummy,
92  vec,
93  dummy,
94  dummy,
96 
97  // compute the contribution of this element to compliance
98  _compliance -= vec.dot(e.sol());
99  }
100 }
101 
102 
103 
104 void
106 
107  // make sure that this has not been initialized ana calculated for all elems
108  libmesh_assert(_physics_elem);
109 
111 
113  dynamic_cast<MAST::StructuralElementBase&>(*_physics_elem);
114 
116  vec = RealVectorX::Zero(e.sol().size());
117 
119  dummy = RealMatrixX::Zero(vec.size(), vec.size());
120 
122  vec,
123  dummy,
124  dummy,
127  vec,
128  dummy,
129  dummy,
131 
132  // ask for the values
133  _dcompliance_dp -= vec.dot(e.sol());
134 
135  vec.setZero();
136  e.side_external_residual(false,
137  vec,
138  dummy,
139  dummy,
141  e.volume_external_residual(false,
142  vec,
143  dummy,
144  dummy,
146 
147  // ask for the values
148  _dcompliance_dp -= vec.dot(e.sol(true));
149  }
150 }
151 
152 
153 
154 void
158 
159  // the primal data should have been calculated
160  libmesh_assert(_physics_elem);
161  libmesh_assert(f.is_topology_parameter());
162 
164  &elem = dynamic_cast<const MAST::LevelSetIntersectedElem&>(_physics_elem->elem());
165 
166  // sensitivity only exists at the boundary. So, we proceed with calculation
167  // only if this element has an intersection in the interior, or with a side.
168 
169  if (this->if_evaluate_for_element(elem) &&
170  elem.if_elem_has_level_set_boundary() &&
171  elem.if_subelem_has_side_on_level_set_boundary()) {
172 
174  dynamic_cast<MAST::StructuralElementBase&>(*_physics_elem);
175 
177  vec = RealVectorX::Zero(e.sol().size());
178 
180  dummy = RealMatrixX::Zero(vec.size(), vec.size());
181 
183  elem.get_subelem_side_on_level_set_boundary(),
184  vel,
186  false,
187  vec,
188  dummy);
189 
190  // compute the contribution of this element to compliance
191  _dcompliance_dp -= vec.dot(e.sol());
192  }
193 }
194 
195 
196 
197 Real
199 
200  Real val = _compliance;
201  _system->system().comm().sum(val);
202  return val;
203 }
204 
205 
206 
207 Real
209 
210  Real val = _dcompliance_dp;
211  _system->system().comm().sum(val);
212  return val;
213 }
214 
215 
216 
217 void
219 
220  // make sure that this has not been initialized ana calculated for all elems
221  libmesh_assert(_physics_elem);
222 
223  // since compliance = -X^T F, derivative wrt X = -F.
224 
226 
228  dynamic_cast<MAST::StructuralElementBase&>(*_physics_elem);
229 
230  dq_dX.setZero();
231 
233  dummy = RealMatrixX::Zero(dq_dX.size(), dq_dX.size());
234 
235  e.side_external_residual(false,
236  dq_dX,
237  dummy,
238  dummy,
240  e.volume_external_residual(false,
241  dq_dX,
242  dummy,
243  dummy,
245 
246  dq_dX *= -1.;
247  }
248 }
249 
250 
251 
252 void
254 set_elem_data(unsigned int dim,
255  const libMesh::Elem& ref_elem,
256  MAST::GeomElem& elem) const {
257 
258  libmesh_assert(!_physics_elem);
259 
260  if (dim == 1) {
261 
262  const MAST::ElementPropertyCard1D& p =
263  dynamic_cast<const MAST::ElementPropertyCard1D&>(_discipline->get_property_card(ref_elem));
264 
265  elem.set_local_y_vector(p.y_vector());
266  }
267 }
268 
269 
270 void
272 
273  libmesh_assert(!_physics_elem);
274  libmesh_assert(_assembly);
275  libmesh_assert(_system);
276 
278  dynamic_cast<const MAST::ElementPropertyCardBase&>(_discipline->get_property_card(elem));
279 
280  _physics_elem =
281  MAST::build_structural_element(*_system, *_assembly, elem, p).release();
282 }
283 
284 
virtual void zero_for_sensitivity()
zeroes the output quantity values stored inside this object so that assembly process can begin...
const MAST::VolumeBCMapType & volume_loads() const
MAST::NonlinearSystem & system()
const RealVectorX & sol(bool if_sens=false) const
Definition: elem_base.cpp:53
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.
virtual void init(const MAST::GeomElem &elem)
initialize for the element.
virtual bool if_evaluate_for_element(const MAST::GeomElem &elem) const
checks to see if the object has been told about the subset of elements and if the specified element i...
std::unique_ptr< MAST::StructuralElementBase > build_structural_element(MAST::SystemInitialization &sys, MAST::AssemblyBase &assembly, const MAST::GeomElem &elem, const MAST::ElementPropertyCardBase &p)
builds the structural element for the specified element type
virtual void zero_for_analysis()
zeroes the output quantity values stored inside this object so that assembly process can begin...
RealVectorX & y_vector()
returns value of the property val.
const MAST::GeomElem & elem() const
Definition: elem_base.h:117
This provides the base class for definitin of element level contribution of output quantity in an ana...
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
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.
libMesh::Real Real
MAST::PhysicsDisciplineBase * _discipline
virtual Real output_sensitivity_total(const MAST::FunctionBase &p)
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.
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
virtual void evaluate_topology_sensitivity(const MAST::FunctionBase &f, const MAST::FieldFunction< RealVectorX > &vel)
This evaluates the contribution to the topology sensitivity on the boundary.
virtual void evaluate_sensitivity(const MAST::FunctionBase &f)
this evaluates all relevant stress sensitivity components on the element to evaluate the p-averaged q...
Matrix< Real, Dynamic, Dynamic > RealMatrixX
ComplianceOutput()
default constructor
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.
This class inherits from MAST::GeomElem and provides an interface to initialize FE objects on sub-ele...
Matrix< Real, Dynamic, 1 > RealVectorX
const MAST::ElementPropertyCardBase & get_property_card(const libMesh::Elem &elem) const
get property card for the specified element
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 evaluate()
this evaluates all relevant stress components on the element to evaluate the p-averaged quantity...
virtual void output_derivative_for_elem(RealVectorX &dq_dX)
calculates the derivative of p-norm von Mises stress for the norm identified using set_p_val()...
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
const MAST::SideBCMapType & side_loads() const
virtual bool is_topology_parameter() const
Definition: function_base.h:97
MAST::SystemInitialization * _system