MAST
Multidisciplinary-design Adaptation and Sensitivity Toolkit (MAST)
integrated_force_output.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 // MAST includes
24 #include "base/assembly_base.h"
25 #include "base/elem_base.h"
27 #include "base/nonlinear_system.h"
28 #include "mesh/geom_elem.h"
29 
30 
31 // libMesh includes
32 #include "libmesh/boundary_info.h"
33 #include "libmesh/parallel.h"
34 
35 
38 _n_vec (nvec),
39 _force (0.),
40 _force_sens (0.) {
41 
42  _n_vec /= _n_vec.norm();
43 }
44 
45 
46 
48 
49 }
50 
51 
52 
53 void
55 
56  libmesh_assert(!_physics_elem);
57  libmesh_assert(_system);
58  libmesh_assert(_assembly);
59 
60  const MAST::FlightCondition& p =
62  (_assembly->discipline()).flight_condition();
63 
66 }
67 
68 
69 
70 void
72 
73  _force = 0.;
74  _force_sens = 0.;
75 }
76 
77 
78 void
80 
81  // nothing to be done here
82  _force_sens = 0.;
83 }
84 
85 
86 Real
88 
89  Real val = _force;
90 
91  if (!_skip_comm_sum)
92  _system->system().comm().sum(val);
93 
94  return val;
95 }
96 
97 
98 
99 Real
101 
102  Real val = _force_sens;
103 
104  if (!_skip_comm_sum)
105  _system->system().comm().sum(val);
106 
107  return val;
108 }
109 
110 
111 
112 void
114 
115  libmesh_assert(_physics_elem);
116 
118  dynamic_cast<MAST::ConservativeFluidElementBase&>(*_physics_elem);
119 
120  const MAST::GeomElem&
121  elem = _physics_elem->elem();
122 
124  f = RealVectorX::Zero(3);
125 
126  for (unsigned short int n=0; n<elem.n_sides_quadrature_elem(); n++)
127  if (this->if_evaluate_for_boundary(elem, n)) {
128 
129  e.side_integrated_force(n, f);
130  _force += f.dot(_n_vec);
131  }
132 }
133 
134 
135 
136 void
138 
139  libmesh_assert(_physics_elem);
140 
142  dynamic_cast<MAST::ConservativeFluidElementBase&>(*_physics_elem);
143 
144  const MAST::GeomElem&
145  elem = _physics_elem->elem();
146 
148  df = RealVectorX::Zero(3);
149 
150  for (unsigned short int n=0; n<elem.n_sides_quadrature_elem(); n++)
151  if (this->if_evaluate_for_boundary(elem, n)) {
152 
154  _force_sens += df.dot(_n_vec);
155  }
156 }
157 
158 
159 
160 void
162 
163  libmesh_assert(_physics_elem);
164 
166  dynamic_cast<MAST::ConservativeFluidElementBase&>(*_physics_elem);
167 
168  const MAST::GeomElem&
169  elem = _physics_elem->elem();
170 
172  f = RealVectorX::Zero(3);
173 
175  dfdX = RealMatrixX::Zero(3, dq_dX.size());
176 
177  for (unsigned short int n=0; n<elem.n_sides_quadrature_elem(); n++)
178  if (this->if_evaluate_for_boundary(elem, n)) {
179 
180  e.side_integrated_force(n, f, &dfdX);
181  dq_dX += _n_vec.transpose() * dfdX;
182  }
183 }
184 
185 
186 
MAST::NonlinearSystem & system()
This class provides the necessary functionality for spatial discretization of the conservative fluid ...
virtual Real output_sensitivity_total(const MAST::FunctionBase &p)
Real _force
integrated value of the force
virtual void side_integrated_force_sensitivity(const MAST::FunctionBase &p, const unsigned int s, RealVectorX &f)
This provides the base class for definitin of element level contribution of output quantity in an ana...
libMesh::Real Real
virtual void side_integrated_force(const unsigned int s, RealVectorX &f, RealMatrixX *dfdX=nullptr)
surface integrated force
virtual void zero_for_sensitivity()
zeroes the output quantity values stored inside this object so that assembly process can begin...
unsigned int n_sides_quadrature_elem() const
number of sides on quadrature element.
Definition: geom_elem.cpp:101
virtual void output_derivative_for_elem(RealVectorX &dq_dX)
returns the output quantity derivative with respect to state vector in dq_dX.
bool _skip_comm_sum
If an output has contrinutions only from local processor then the user can request that the global co...
Matrix< Real, Dynamic, Dynamic > RealMatrixX
Matrix< Real, Dynamic, 1 > RealVectorX
virtual void evaluate_sensitivity(const MAST::FunctionBase &p)
this evaluates all relevant sensitivity components on the 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 init(const MAST::GeomElem &elem)
initialize for the element.
const MAST::GeomElem & elem() const
Definition: elem_base.h:109
Real _force_sens
integrated value of the sensitivity of force
IntegratedForceOutput(const RealVectorX &nvec)
RealVectorX _n_vec
vector along which the force is measured
const MAST::PhysicsDisciplineBase & discipline() const
virtual void zero_for_analysis()
zeroes the output quantity values stored inside this object so that assembly process can begin...
virtual void evaluate()
this is the abstract interface to be implemented by derived classes.
virtual bool if_evaluate_for_boundary(const MAST::GeomElem &elem, const unsigned int s) const
checks to see if the specified side of the element needs evaluation of the output contribution...
MAST::SystemInitialization * _system