MAST
flutter_solver_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 #ifndef __mast__flutter_solver_base_h__
21 #define __mast__flutter_solver_base_h__
22 
23 // C++ includes
24 #include <string>
25 #include <fstream>
26 #include <iomanip>
27 
28 
29 // MAST includes
30 #include "base/mast_data_types.h"
31 
32 
33 // libMesh includes
34 #include "libmesh/numeric_vector.h"
35 
36 
37 
38 namespace MAST {
39 
40  // Forward declerations
41  class FunctionBase;
42  class FlutterModel;
43  class FlutterRootBase;
44  class FlutterSolutionBase;
45  class FlutterRootCrossoverBase;
46  class StructuralFluidInteractionAssembly;
47  template <typename ValType> class BasisMatrix;
48 
49 
51 
52  public:
53 
58 
59 
60  virtual ~FlutterSolverBase();
61 
62 
67  class SteadySolver {
68  public:
69 
70 
72 
73  virtual ~SteadySolver() {}
74 
79  virtual const libMesh::NumericVector<Real>&
80  solve() = 0;
81 
82 
86  virtual const libMesh::NumericVector<Real>&
87  solution() const = 0;
88 
89  };
90 
91 
96 
97 
102 
103 
107  virtual void clear();
108 
109 
113  virtual void clear_assembly_object();
114 
115 
119  void initialize(std::vector<libMesh::NumericVector<Real>*>& basis);
120 
121 
122 
123  void set_output_file(const std::string& nm) {
124 
125  if (!_output)
126  _output = new std::ofstream;
127 
128  _output->close();
129  _output->open(nm.c_str(), std::ofstream::out);
130  }
131 
132 
136  virtual void print_sorted_roots() = 0;
137 
138 
143  virtual void print_crossover_points() = 0;
144 
145 
146  protected:
147 
148 
154 
155 
159  std::vector<libMesh::NumericVector<Real>*>* _basis_vectors;
160 
161 
165  std::ofstream* _output;
166 
167 
172 
173  };
174 }
175 
176 
177 #endif // __mast__flutter_solver_base_h__
void set_output_file(const std::string &nm)
virtual const libMesh::NumericVector< Real > & solution() const =0
virtual void print_crossover_points()=0
Prints the crossover points output.
std::vector< libMesh::NumericVector< Real > * > * _basis_vectors
basis vector used to define the reduced order model
void attach_steady_solver(MAST::FlutterSolverBase::SteadySolver &solver)
attaches the steady solution object
virtual void print_sorted_roots()=0
Prints the sorted roots to the output.
MAST::FlutterSolverBase::SteadySolver * _steady_solver
object provides the steady state solution.
void attach_assembly(MAST::StructuralFluidInteractionAssembly &assembly)
attaches the assembly object to this solver.
abstract class defines the interface to provide the steady-state solution
std::ofstream * _output
file to which the result will be written
FlutterSolverBase()
defalut constructor
void initialize(std::vector< libMesh::NumericVector< Real > * > &basis)
initializes the data structres for a flutter solution.
virtual void clear_assembly_object()
clears the assembly object
virtual const libMesh::NumericVector< Real > & solve()=0
solves for the steady state solution, and
MAST::StructuralFluidInteractionAssembly * _assembly
structural assembly that provides the assembly of the system matrices.
virtual void clear()
clears the solution and other data from this solver