MAST
frequency_domain_pressure_function.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 "numerics/utility.h"
27 #include "fluid/flight_condition.h"
28 #include "base/nonlinear_system.h"
29 
30 
31 // libMesh includes
32 #include "libmesh/numeric_vector.h"
33 #include "libmesh/system.h"
34 #include "libmesh/dof_map.h"
35 
36 
40 MAST::FieldFunction<Complex>("frequency_domain_pressure"),
41 _if_cp(false),
42 _system(sys),
43 _flt_cond(flt) {
44 
45 }
46 
47 
48 
49 
51 
52 }
53 
54 
55 
56 
57 void
59 init(const libMesh::NumericVector<Real>& steady_sol,
60  const libMesh::NumericVector<Real>& small_dist_sol_real,
61  const libMesh::NumericVector<Real>& small_dist_sol_imag) {
62 
64 
65  // first initialize the solution to the given vector
66  // steady state solution
67  _sol.reset(libMesh::NumericVector<Real>::build(sys.comm()).release());
68  _sol->init(steady_sol.size(), true, libMesh::SERIAL);
69 
70  // solution real part
71  _dsol_real.reset(libMesh::NumericVector<Real>::build(sys.comm()).release());
72  _dsol_real->init(steady_sol.size(), true, libMesh::SERIAL);
73 
74  // solution complex part
75  _dsol_imag.reset(libMesh::NumericVector<Real>::build(sys.comm()).release());
76  _dsol_imag->init(steady_sol.size(), true, libMesh::SERIAL);
77 
78 
79  // now localize the give solution to this objects's vector
80  steady_sol.localize(*_sol);
81  small_dist_sol_real.localize(*_dsol_real);
82  small_dist_sol_imag.localize(*_dsol_imag);
83 
84 
85  // if the mesh function has not been created so far, initialize it
86  _sol_function.reset(new libMesh::MeshFunction(sys.get_equation_systems(),
87  *_sol,
88  sys.get_dof_map(),
89  _system.vars()));
90  _sol_function->init();
91 
92 
93  _dsol_re_function.reset(new libMesh::MeshFunction(sys.get_equation_systems(),
94  *_dsol_real,
95  sys.get_dof_map(),
96  _system.vars()));
97  _dsol_re_function->init();
98 
99 
100  _dsol_im_function.reset(new libMesh::MeshFunction(sys.get_equation_systems(),
101  *_dsol_imag,
102  sys.get_dof_map(),
103  _system.vars()));
104  _dsol_im_function->init();
105 
106 }
107 
108 
109 
110 
111 
112 void
114 operator () (const libMesh::Point& p,
115  const Real t,
116  Complex& dpress) const {
117 
118 
119  libmesh_assert(_sol_function.get()); // should be initialized before this call
120 
121  dpress = 0.;
122 
123 
124  // get the nonlinear and linearized solution
126  v;
127 
129  sol = RealVectorX::Zero(_system.system().n_vars());
131  dsol = ComplexVectorX::Zero(_system.system().n_vars());
132 
133 
134  // first copy the real and imaginary solutions
135  (*_dsol_re_function)(p, 0., v);
136  MAST::copy(sol, v);
137  dsol.real() = sol;
138 
139 
140  // now the imaginary part
141  (*_dsol_im_function)(p, 0., v);
142  MAST::copy(sol, v);
143  dsol.imag() = sol;
144 
145 
146  // now the steady state function itself
147  (*_sol_function)(p, 0., v);
148  MAST::copy(sol, v);
149 
150 
153 
154  // now initialize the primitive variable contexts
155  p_sol.init(dynamic_cast<MAST::ConservativeFluidSystemInitialization&>(_system).dim(),
156  sol,
160  delta_p_sol.init(p_sol,
161  dsol,
163 
164  if (_if_cp)
165  dpress = delta_p_sol.c_pressure(_flt_cond.q0());
166  else
167  dpress = delta_p_sol.dp;
168 }
169 
virtual void operator()(const libMesh::Point &p, const Real t, Complex &dp) const
provides the complex pressure perturbation
MAST::SystemInitialization & _system
system associated with the mesh and solution vector
const std::vector< unsigned int > vars() const
MAST::NonlinearSystem & system()
Class defines basic operations and calculation of the small disturbance primitive variables...
This class implements a system for solution of nonlinear systems.
std::unique_ptr< libMesh::MeshFunction > _dsol_im_function
Class defines the conversion and some basic operations on primitive fluid variables used in calculati...
void init(const MAST::PrimitiveSolution &sol, const typename VectorType< ValType >::return_type &delta_sol, bool if_viscous)
Matrix< Complex, Dynamic, 1 > ComplexVectorX
libMesh::Real Real
libMesh::Complex Complex
std::unique_ptr< libMesh::NumericVector< Real > > _sol
steady part of solution
std::unique_ptr< libMesh::MeshFunction > _sol_function
mesh function that interpolates the solution
libMesh::DenseVector< Real > DenseRealVector
FrequencyDomainPressureFunction(MAST::SystemInitialization &sys, MAST::FlightCondition &flt)
Real q0() const
returns the flight dynamic pressure
This creates the base class for functions that have a saptial and temporal dependence, and provide sensitivity operations with respect to the functions and parameters.
std::unique_ptr< libMesh::NumericVector< Real > > _dsol_imag
imag part of small-disturbance solution
void copy(DenseRealMatrix &m1, const RealMatrixX &m2)
Definition: utility.h:167
Matrix< Real, Dynamic, 1 > RealVectorX
void init(const libMesh::NumericVector< Real > &steady_sol, const libMesh::NumericVector< Real > &small_dist_sol_real, const libMesh::NumericVector< Real > &small_dist_sol_imag)
initiate the mesh function for this solution
void init(const unsigned int dim, const RealVectorX &conservative_sol, const Real cp_val, const Real cv_val, bool if_viscous)
std::unique_ptr< libMesh::MeshFunction > _dsol_re_function
GasProperty gas_property
Ambient air properties.
bool _if_cp
the function will return cp instead of pressure if this option is true.
std::unique_ptr< libMesh::NumericVector< Real > > _dsol_real
real part of small-disturbance solution
MAST::FlightCondition & _flt_cond
flight condition