MAST
time_domain_flutter_solution.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 // C++ includes
21 #include <iomanip>
22 
23 
24 // MAST includes
28 
29 
32 { }
33 
34 
35 
37 { }
38 
39 
40 
41 
42 void
44  const Real v_ref,
45  const MAST::LAPACK_DGGEV& eig_sol) {
46 
47  // make sure that it hasn't already been initialized
48  libmesh_assert(!_roots.size());
49 
50  _ref_val = v_ref;
51 
52  // iterate over the roots and initialize the vector
53  _Amat = eig_sol.A();
54  _Bmat = eig_sol.B();
55  const ComplexMatrixX
56  &VR = eig_sol.right_eigenvectors(),
57  &VL = eig_sol.left_eigenvectors();
58  const ComplexVectorX
59  &num = eig_sol.alphas();
60 
61  const RealVectorX
62  &den = eig_sol.betas();
63 
64  // use only half the modes that have a positive frequency.
65  unsigned int nvals = (int)_Bmat.rows();
66 
67  _roots.resize(nvals);
68  for (unsigned int i=0; i<nvals; i++) {
69 
71  root->init(v_ref,
72  num(i),
73  den(i),
74  _Bmat,
75  VR.col(i),
76  VL.col(i));
77 
78  _roots[i] = root;
79  }
80 }
81 
82 
83 
84 unsigned int
87 
88  unsigned int n=0;
89 
90  // iterate over all the roots and find the ones with g >= 0
91  std::vector<MAST::FlutterRootBase*>::const_iterator
92  it = _roots.begin(),
93  end = _roots.end();
94 
95  for ( ; it != end; it++) {
96 
97  // only look at upper half of the complex domain
98  if ((**it).root.imag() >= 0. &&
99  (**it).root.real() >= tol)
100  n++;
101  }
102 
103  return n;
104 }
105 
106 
107 
108 
109 
112 
113  // If there is an unstable root, then find one with the lowest velocity
114  // and send it back.
115  //
116  // Otherwise, find the root that is closest to the g=0 condition
117 
119  *unstable = nullptr,
120  *max_g = nullptr;
121 
122  // iterate over all the roots and find the ones with g >= 0
123  std::vector<MAST::FlutterRootBase*>::const_iterator
124  it = _roots.begin(),
125  end = _roots.end();
126 
127  for ( ; it != end; it++) {
128  // only look at the upper half of the complex domain.
129  if ((**it).root.imag() >= 0.) {
130 
131  // look for the lowest g root
132  if (!max_g)
133  max_g = *it;
134  else if ((**it).root.real() >= max_g->root.real())
135  max_g = *it;
136 
137  // now look for the unstable root with lowest velocity
138  // only roots evaluated at V > 0 are considered, since in the absensce
139  // of aerodynamics and structural damping, the roots would all
140  // have zero damping.
141 
142  if ((**it).root.real() >= tol &&
143  (**it).V > 0.) {
144 
145  if (!unstable)
146  unstable = *it;
147  else if ( (**it).V <= unstable->V )
148  unstable = *it;
149  }
150  }
151 
152  }
153 
154  if (unstable)
155  return unstable;
156  else
157  return max_g;
158 }
159 
160 
161 
162 
163 
164 
165 void
167 
168  const unsigned int nvals = this->n_roots();
169  libmesh_assert_equal_to(nvals, sol.n_roots());
170 
171  // two roots with highest modal_participation dot product are sorted
172  // in the same serial order
173  for (unsigned int i=0; i<nvals-1; i++)
174  {
175  const MAST::FlutterRootBase& r = sol.get_root(i);
176  Real max_val = 0.;
177  Complex val = 0.;
178  unsigned int max_val_root = nvals-1;
179  for (unsigned int j=i; j<nvals; j++) {
180 
181  // use a combination of both the eigenvectors from both
182  // roots
183  val = .5*(r.eig_vec_left.dot(_Bmat*_roots[j]->eig_vec_right) +
184  _roots[j]->eig_vec_left.dot(_Bmat*r.eig_vec_right));
185  //_roots[j]->modal_participation.dot(r.modal_participation);
186  // scale by the eigenvalue separation with the assumption that
187  // the roots will be closer to each other than any other
188  // root at two consecutive eigenvalues. In other words,
189  // we are penalizing the dot product with the eigenvalue
190  // distance
191  val /= abs(r.root-_roots[j]->root);
192  if (abs(val) > max_val) {
193  max_val = abs(val);
194  max_val_root = j;
195  }
196  }
197 
198  // now we should have the one with highest dot product
199  if (i != max_val_root)
200  std::swap(_roots[i], _roots[max_val_root]);
201  }
202 }
203 
204 
205 
206 
207 
208 
209 void
211 
212  // make sure that some roots exist for this solution
213  libmesh_assert(this->n_roots() > 0);
214 
215  const unsigned int nvals = this->n_roots();
216  unsigned int n_participation_vals =
217  (unsigned int)this->get_root(0).modal_participation.size();
218  libmesh_assert(nvals);
219 
220  // first write the reference values of the root
221  output << " Flutter Root " << std::endl;
222 
223  // now write the root
224  output
225  << std::setw(5) << "#"
226  << std::setw(15) << "V_ref"
227  << std::setw(15) << "Re"
228  << std::setw(15) << "Im"
229  << std::setw(3) << " | ";
230 
231  // output the headers for flutter mode participation
232  for (unsigned int i=0; i<n_participation_vals; i++)
233  output
234  << std::setw(2) << " "
235  << std::setw(5) << "Mode "
236  << std::setw(5) << i
237  << std::setw(2) << " ";
238 
239  // output the headers for flutter mode
240  for (unsigned int i=0; i<nvals; i++)
241  output
242  << std::setw(10) << "| "
243  << std::setw(5) << "Mode "
244  << std::setw(5) << i
245  << std::setw(10) << " |";
246  output << std::endl;
247 
248  for (unsigned int i=0; i<nvals; i++)
249  {
250  const MAST::FlutterRootBase& root = this->get_root(i);
251 
252  // flutter root details
253  output
254  << std::setw(5) << i
255  << std::setw(15) << root.V
256  << std::setw(15) << std::real(root.root)
257  << std::setw(15) << std::imag(root.root)
258  << std::setw(3) << " | ";
259 
260  // now write the modal participation
261  for (unsigned int j=0; j<n_participation_vals; j++)
262  output
263  << std::setw(12) << root.modal_participation(j)
264  << std::setw(2) << " ";
265 
266  // now write the flutter mode
267  for (unsigned int j=0; j<nvals; j++)
268  {
269  output
270  << std::setw(2) << "| "
271  << std::setw(12) << std::real(root.eig_vec_right(j))
272  << std::setw(2) << " "
273  << std::setw(12) << std::imag(root.eig_vec_right(j))
274  << std::setw(2) << " |";
275  }
276  output << std::endl;
277  }
278 }
279 
280 
const MAST::FlutterRootBase & get_root(const unsigned int i) const
ComplexVectorX eig_vec_left
std::vector< MAST::FlutterRootBase * > _roots
virtual void print(std::ostream &output)
prints the data and modes from this solution
unsigned int n_unstable_roots_in_upper_complex_half(Real tol) const
number of unstable roots in this solution.
void init(const Real v_ref_val, const Complex num, const Complex den, const RealMatrixX &Bmat, const ComplexVectorX &evec_right, const ComplexVectorX &evec_left)
initializes the data
Matrix< Complex, Dynamic, 1 > ComplexVectorX
libMesh::Real Real
libMesh::Complex Complex
const RealVectorX & betas() const
virtual ~TimeDomainFlutterSolution()
delete the flutter root objects
const ComplexMatrixX & right_eigenvectors() const
void init(const MAST::TimeDomainFlutterSolver &solver, const Real v_ref, const MAST::LAPACK_DGGEV &eig_sol)
initializes the root
unsigned int n_roots() const
number of roots in this solution
const RealMatrixX & B() const
const ComplexVectorX & alphas() const
const ComplexMatrixX & left_eigenvectors() const
Matrix< Complex, Dynamic, Dynamic > ComplexMatrixX
MAST::FlutterRootBase * get_critical_root(Real tol)
Matrix< Real, Dynamic, 1 > RealVectorX
RealMatrixX _Amat
Matrix used for scaling of eigenvectors, and sorting of roots.
const RealMatrixX & A() const
Real _ref_val
Reference value of the sweeping parameter for which this solution was obtained.
virtual void sort(const MAST::FlutterSolutionBase &sol)
sort this root with respect to the given solution from a previous eigen solution. ...
This implements a solver for a single parameter instability problem, for example a flutter solver whe...
RealVectorX modal_participation
ComplexVectorX eig_vec_right
right and left eigenvevtors