VM2D  1.12
Vortex methods for 2D flows simulation
VM2D::BoundaryConstLayerAver Class Reference

Класс, определяющий способ удовлетворения граничного условия на обтекаемом профиле More...

#include <Boundary2DConstLayerAver.h>

Inheritance diagram for VM2D::BoundaryConstLayerAver:
Collaboration diagram for VM2D::BoundaryConstLayerAver:

Public Member Functions

 BoundaryConstLayerAver (const World2D &W_, size_t numberInPassport_)
 Конструктор More...
 
virtual ~BoundaryConstLayerAver ()
 Деструктор More...
 
virtual void FillMatrixSelf (Eigen::MatrixXd &matr, Eigen::VectorXd &lastLine, Eigen::VectorXd &lactCol) override
 Генерация блока матрицы More...
 
virtual void FillIQSelf (std::pair< Eigen::MatrixXd, Eigen::MatrixXd > &IQ) override
 Генерация блока матрицы, состоящей из интегралов от (r-xi)/|r-xi|^2, влияние профиля самого на себя More...
 
virtual void FillMatrixFromOther (const Boundary &otherBoundary, Eigen::MatrixXd &matr) override
 Генерация блока матрицы влияния от другого профиля того же типа More...
 
virtual void FillIQFromOther (const Boundary &otherBoundary, std::pair< Eigen::MatrixXd, Eigen::MatrixXd > &IQ) override
 Генерация блока матрицы, состоящей из интегралов от (r-xi)/|r-xi|^2, влияние одного профиля на другой More...
 
virtual void SolutionToFreeVortexSheetAndVirtualVortex (const Eigen::VectorXd &sol) override
 Пересчет решения на интенсивность вихревого слоя и на рождаемые вихри на конкретном профиле More...
 
virtual void CalcConvVelocityToSetOfPointsFromSheets (const WakeDataBase &pointsDb, std::vector< Point2D > &velo) const override
 Вычисление конвективных скоростей в наборе точек, вызываемых наличием слоев вихрей и источников на профиле More...
 
virtual void CalcConvVelocityAtVirtualVortexes (std::vector< Point2D > &velo) const override
 Вычисление конвективной скорости в точках расположения виртуальных вихрей More...
 
virtual void ComputeAttachedSheetsIntensity () override
 Вычисление интенсивностей присоединенного вихревого слоя и присоединенного слоя источников More...
 
virtual void GetInfluenceFromVorticesToRectPanel (size_t panel, const Vortex2D *ptr, ptrdiff_t count, std::vector< double > &wakeRhs) const override
 Вычисление влияния части подряд идущих вихрей из вихревого следа на прямолинейную панель для правой части More...
 
virtual void GetInfluenceFromSourcesToRectPanel (size_t panel, const Vortex2D *ptr, ptrdiff_t count, std::vector< double > &wakeRhs) const override
 Вычисление влияния части подряд источников из области течения на прямолинейную панель для правой части More...
 
virtual void GetInfluenceFromVortexSheetAtRectPanelToVortex (size_t panel, const Vortex2D &vtx, Point2D &vel) const override
 Вычисление влияния вихревых слоев (свободный + присоединенный) конкретной прямолинейной панели на вихрь в области течения ///. More...
 
virtual void GetInfluenceFromSourceSheetAtRectPanelToVortex (size_t panel, const Vortex2D &ptr, Point2D &vel) const override
 Вычисление влияния слоя источников конкретной прямолинейной панели на вихрь в области течения More...
 
virtual void GetInfluenceFromVInfToRectPanel (std::vector< double > &vInfRhs) const override
 Вычисление влияния набегающего потока на прямолинейную панель для правой части More...
 
size_t GetUnknownsSize () const
 Возврат размерности вектора решения More...
 

Public Attributes

const Airfoilafl
 
int minVortexPerPanel
 Минимальное число вихрей, рождаемых на панели профиля и формирующих виртуальный вихревой след More...
 
std::vector< std::pair< int, int > > vortexBeginEnd
 Номера первого и последнего вихрей, рождаемых на каждой панели профиля (формируется после решения СЛАУ) More...
 
VirtualWake virtualWake
 Виртуальный вихревой след конкретного профиля More...
 
size_t sheetDim
 Размерность параметров каждого из слоев на каждой из панелей More...
 
Sheet sheets
 Слои на профиле More...
 
Sheet oldSheets
 Слои на профиле с предыдущего шага More...
 

Protected Attributes

const World2DW
 Константная ссылка на решаемую задачу More...
 
const size_t numberInPassport
 Номер профиля в паспорте More...
 

Detailed Description

Класс, определяющий способ удовлетворения граничного условия на обтекаемом профиле

Способ удовлетворения граничного условия:

  • генерация вихревого слоя постоянной интенсивности на панелях;
  • условие ортогональности невязки граничного условия константе (выполнение граничного условия в среднем по панели).
Author
Марчевский Илья Константинович
Сокол Ксения Сергеевна
Рятина Евгения Павловна
Колганова Александра Олеговна
\Version 1.12
\date 14 января 2024 г.

Definition at line 66 of file Boundary2DConstLayerAver.h.

Constructor & Destructor Documentation

VM2D::BoundaryConstLayerAver::BoundaryConstLayerAver ( const World2D W_,
size_t  numberInPassport_ 
)
inline

Конструктор

Parameters
[in]W_константная ссылка на решаемую задачу
[in]numberInPassport_номер профиля в паспорте задачи

Definition at line 75 of file Boundary2DConstLayerAver.h.

75  :
76  Boundary(W_, numberInPassport_, 1) {};
Boundary(const World2D &W_, size_t numberInPassport_, int sheetDim_)
Конструктор
Definition: Boundary2D.cpp:53
virtual VM2D::BoundaryConstLayerAver::~BoundaryConstLayerAver ( )
inlinevirtual

Деструктор

Definition at line 79 of file Boundary2DConstLayerAver.h.

79 {};

Here is the call graph for this function:

Member Function Documentation

void BoundaryConstLayerAver::CalcConvVelocityAtVirtualVortexes ( std::vector< Point2D > &  velo) const
overridevirtual

Вычисление конвективной скорости в точках расположения виртуальных вихрей

Вычисление производится в точках расположения виртуальных вихрей. Скорости находятся как скорости соответствующих точек профиля плюс половина интенсивности присоединенного вихревого слоя, умноженная на касательную к профилю на данной панели.

Parameters
[out]veloссылка на заполняемый список скоростей
Warning
Массив velo заполняется путем присвоения, а не прибавления значений

Implements VM2D::Boundary.

Definition at line 257 of file Boundary2DConstLayerAver.cpp.

258 {
259  std::vector<Point2D>& Vel = velo;
260 
261  //Скорости виртуальных вихрей
262 
263  /*
264  std::stringstream ss;
265  ss << afl.numberInPassport << "-" << W.currentStep;
266  std::ofstream of(W.getPassport().dir + "dbg/tele_" + ss.str());
267  for (int i = 0; i < (int)Vel.size(); ++i)
268  of << afl.getV(virtualWake.aflPan[i].second) << " " << virtualWake.vecHalfGamma[i] << " " << W.getPassport().physicalProperties.V0() << std::endl;
269  of.close();
270  */
271 
272  //std::cout << virtualWake.vecHalfGamma.size() << std::endl;
273 #pragma omp parallel for default(none) shared(Vel)
274  for (int i = 0; i < (int)Vel.size(); ++i)
275  Vel[i] = afl.getV(virtualWake.aflPan[i].second) \
277  - W.getPassport().physicalProperties.V0(); // V0 потом прибавляется ко всем скоростям в функции MoveVortexes
278 
279 
280 }//CalcConvVelocityAtVirtualVortexes(...)
std::vector< Point2D > vecHalfGamma
Скорость вихрей виртуального следа конкретного профиля (равна Gamma/2) используется для расчета давле...
Definition: VirtualWake2D.h:70
PhysicalProperties physicalProperties
Структура с физическими свойствами задачи
Definition: Passport2D.h:296
const Point2D & getV(size_t q) const
Возврат константной ссылки на скорость вершины профиля
Definition: Airfoil2D.h:113
std::vector< std::pair< size_t, size_t > > aflPan
Пара чисел: номер профиля и номер панели, на которой рожден виртуальный вихрь
Definition: VirtualWake2D.h:73
const Passport & getPassport() const
Возврат константной ссылки на паспорт
Definition: World2D.h:222
VirtualWake virtualWake
Виртуальный вихревой след конкретного профиля
Definition: Boundary2D.h:85
Point2D V0() const
Функция скорости набегающего потока с учетом разгона
Definition: Passport2D.h:93
const World2D & W
Константная ссылка на решаемую задачу
Definition: Boundary2D.h:67
const Airfoil & afl
Definition: Boundary2D.h:76

Here is the call graph for this function:

Here is the caller graph for this function:

void BoundaryConstLayerAver::CalcConvVelocityToSetOfPointsFromSheets ( const WakeDataBase pointsDb,
std::vector< Point2D > &  velo 
) const
overridevirtual

Вычисление конвективных скоростей в наборе точек, вызываемых наличием слоев вихрей и источников на профиле

Вычисляет конвективные скорости в наборе точек, которые вызваны влиянием слоев вихрей и источников на профиле

Parameters
[in]pointsDbконстантная ссылка на базу данных вихрей, в точках которых вычисляются скорости
[out]veloссылка на вектор скоростей, которые приобретают точки из-за влияния слоев вихрей и источников на профиле
Warning
velo — накапливается!
Todo:
сделать вызов функции GetInfluenceFromVortexSheetAtRectPanelToVortex
Todo:
Тут надо разобраться, как должно быть...
Todo:
сделать if(move || deform)

Implements VM2D::Boundary.

Definition at line 158 of file Boundary2DConstLayerAver.cpp.

159 {
160  std::vector<Point2D> selfVelo(pointsDb.vtx.size());
161 
162  double cft = IDPI;
163 
164 #pragma warning (push)
165 #pragma warning (disable: 4101)
166  //Локальные переменные для цикла
167  Point2D velI;
168  Point2D tempVel;
169  //double dst2eps, dst2;
170 #pragma warning (pop)
171 
172 #pragma omp parallel for default(none) shared(selfVelo, cft, pointsDb, std::cout) private(velI, tempVel)
173  for (int i = 0; i < pointsDb.vtx.size(); ++i)
174  {
176 
177  velI.toZero();
178 
179  const Point2D& posI = pointsDb.vtx[i].r();
180 
183  for (size_t j = 0; j < sheets.getSheetSize(); ++j)
184  {
185  Point2D dj = afl.getR(j + 1) - afl.getR(j);
186  Point2D tauj = dj.unit();
187 
188  Point2D s = posI - afl.getR(j);
189  Point2D p = posI - afl.getR(j + 1);
190 
191  double a = VMlib::Alpha(p, s);
192 
193  double lambda;
194  if ( (s.length2() > 1e-16) && (p.length2() > 1e-16) )
195  lambda = VMlib::Lambda(p, s);
196  else
197  lambda = 0.0;
198 
199  Point2D skos = -a * tauj.kcross() + lambda * tauj;
200 
201  velI += sheets.freeVortexSheet(j, 0) * skos.kcross();
202  velI += sheets.attachedVortexSheet(j, 0) * skos.kcross();
203  velI += sheets.attachedSourceSheet(j, 0) * skos;
204  }//for j
205 
206  velI *= cft;
207  selfVelo[i] = velI;
208  }//for i
209 
210  for (size_t i = 0; i < velo.size(); ++i)
211  velo[i] += selfVelo[i];
212 }//CalcConvVelocityToSetOfPointsFromSheets(...)
const double & freeVortexSheet(size_t n, size_t moment) const
Definition: Sheet2D.h:99
const double & attachedVortexSheet(size_t n, size_t moment) const
Definition: Sheet2D.h:104
const double IDPI
Число .
Definition: defs.h:76
const Point2D & getR(size_t q) const
Возврат константной ссылки на вершину профиля
Definition: Airfoil2D.h:101
double Alpha(const Point2D &p, const Point2D &s)
Вспомогательная функция вычисления угла между векторами
Definition: defs.cpp:259
auto length2() const -> typename std::remove_const< typename std::remove_reference< decltype(this->data[0])>::type >::type
Вычисление квадрата нормы (длины) вектора
Definition: numvector.h:383
double Lambda(const Point2D &p, const Point2D &s)
Вспомогательная функция вычисления логарифма отношения норм векторов
Definition: defs.cpp:265
std::vector< Vortex2D > vtx
Список вихревых элементов
size_t getSheetSize() const
Definition: Sheet2D.h:94
Класс, опеделяющий двумерный вектор
Definition: Point2D.h:75
auto unit(P newlen=1) const -> numvector< typename std::remove_const< decltype(this->data[0]*newlen)>::type, n >
Вычисление орта вектора или вектора заданной длины, коллинеарного данному
Definition: numvector.h:399
Sheet sheets
Слои на профиле
Definition: Boundary2D.h:95
const double & attachedSourceSheet(size_t n, size_t moment) const
Definition: Sheet2D.h:109
numvector< T, n > & toZero(P val=0)
Установка всех компонент вектора в константу (по умолчанию — нуль)
Definition: numvector.h:524
const Airfoil & afl
Definition: Boundary2D.h:76
numvector< T, 2 > kcross() const
Геометрический поворот двумерного вектора на 90 градусов
Definition: numvector.h:507

Here is the call graph for this function:

Here is the caller graph for this function:

void BoundaryConstLayerAver::ComputeAttachedSheetsIntensity ( )
overridevirtual

Вычисление интенсивностей присоединенного вихревого слоя и присоединенного слоя источников

Implements VM2D::Boundary.

Definition at line 284 of file Boundary2DConstLayerAver.cpp.

285 {
286 
287  for (size_t i = 0; i < sheets.getSheetSize(); ++i)
288  {
291  }
292 
293  for (size_t i = 0; i < sheets.getSheetSize(); ++i)
294  {
295  sheets.attachedVortexSheet(i, 0) = 0.5 * (afl.getV(i) + afl.getV(i + 1)) & afl.tau[i];
296  sheets.attachedSourceSheet(i, 0) = 0.5 * (afl.getV(i) + afl.getV(i + 1)) & afl.nrm[i];
297  }
298 }//ComputeAttachedSheetsIntensity()
const double & attachedVortexSheet(size_t n, size_t moment) const
Definition: Sheet2D.h:104
const Point2D & getV(size_t q) const
Возврат константной ссылки на скорость вершины профиля
Definition: Airfoil2D.h:113
size_t getSheetSize() const
Definition: Sheet2D.h:94
std::vector< Point2D > tau
Касательные к панелям профиля
Definition: Airfoil2D.h:182
Sheet sheets
Слои на профиле
Definition: Boundary2D.h:95
const double & attachedSourceSheet(size_t n, size_t moment) const
Definition: Sheet2D.h:109
Sheet oldSheets
Слои на профиле с предыдущего шага
Definition: Boundary2D.h:98
std::vector< Point2D > nrm
Нормали к панелям профиля
Definition: Airfoil2D.h:177
const Airfoil & afl
Definition: Boundary2D.h:76

Here is the call graph for this function:

Here is the caller graph for this function:

void BoundaryConstLayerAver::FillIQFromOther ( const Boundary otherBoundary,
std::pair< Eigen::MatrixXd, Eigen::MatrixXd > &  IQ 
)
overridevirtual

Генерация блока матрицы, состоящей из интегралов от (r-xi)/|r-xi|^2, влияние одного профиля на другой

Генерирует блок матрицы влияния от другого профиля того же типа

Parameters
[in]otherBoundaryконстантная ссылка на граничное условие на втором профиле
[out]IQссылка на пару матриц, выражающих взаимные влияния (касательные и нормальные) панелей профиля
Todo:
Пока считается, что граничные условия одинаковые

Implements VM2D::Boundary.

Definition at line 151 of file Boundary2DConstLayerAver.cpp.

152 {
153  afl.calcIQ(1, otherBoundary.afl, IQ);
154 }//FillIQFromOther(...)
virtual void calcIQ(size_t p, const Airfoil &otherAirfoil, std::pair< Eigen::MatrixXd, Eigen::MatrixXd > &matrPair) const =0
Вычисление коэффициентов матрицы, состоящей из интегралов от (r-xi)/|r-xi|^2.
const Airfoil & afl
Definition: Boundary2D.h:76

Here is the call graph for this function:

Here is the caller graph for this function:

void BoundaryConstLayerAver::FillIQSelf ( std::pair< Eigen::MatrixXd, Eigen::MatrixXd > &  IQ)
overridevirtual

Генерация блока матрицы, состоящей из интегралов от (r-xi)/|r-xi|^2, влияние профиля самого на себя

Parameters
[out]IQссылка на генерируемую матрицу

Implements VM2D::Boundary.

Definition at line 137 of file Boundary2DConstLayerAver.cpp.

138 {
139  afl.calcIQ(1, afl, IQ);
140 }//FillIQSelf(...)
virtual void calcIQ(size_t p, const Airfoil &otherAirfoil, std::pair< Eigen::MatrixXd, Eigen::MatrixXd > &matrPair) const =0
Вычисление коэффициентов матрицы, состоящей из интегралов от (r-xi)/|r-xi|^2.
const Airfoil & afl
Definition: Boundary2D.h:76

Here is the call graph for this function:

Here is the caller graph for this function:

void BoundaryConstLayerAver::FillMatrixFromOther ( const Boundary otherBoundary,
Eigen::MatrixXd &  matr 
)
overridevirtual

Генерация блока матрицы влияния от другого профиля того же типа

Генерирует блок матрицы влияния от другого профиля того же типа

Parameters
[in]otherBoundaryконстантная ссылка на граничное условие на втором профиле
[out]matrссылка на генерируемый блок матрицы
Todo:
Пока считается, что граничные условия одинаковые

Implements VM2D::Boundary.

Definition at line 143 of file Boundary2DConstLayerAver.cpp.

144 {
145  for (size_t i = 0; i < afl.getNumberOfPanels(); ++i)
146  for (size_t j = 0; j < otherBoundary.afl.getNumberOfPanels(); ++j)
147  matr(i, j) = afl.getA(1, i, otherBoundary.afl, j)[0];
148 }//FillMatrixFromOther(...)
virtual std::vector< double > getA(size_t p, size_t i, const Airfoil &airfoil, size_t j) const =0
Вычисление коэффициентов матрицы A для расчета влияния панели на панель
size_t getNumberOfPanels() const
Возврат количества панелей на профиле
Definition: Airfoil2D.h:139
const Airfoil & afl
Definition: Boundary2D.h:76

Here is the call graph for this function:

Here is the caller graph for this function:

void BoundaryConstLayerAver::FillMatrixSelf ( Eigen::MatrixXd &  matr,
Eigen::VectorXd &  lastLine,
Eigen::VectorXd &  lactCol 
)
overridevirtual

Генерация блока матрицы

Генерирует следующие компоненты матрицы:

  • диагональный блок матрицы — влияние данного профиля на самого себя;
  • нижнюю строку для матрицы для данного профиля;
  • правый столбец матрицы для данного профиля.
Parameters
[out]matrссылка на генерируемую матрицу
[out]lastLineссылка на нижнюю строку
[out]lactColссылка на правый столбец

Implements VM2D::Boundary.

Definition at line 121 of file Boundary2DConstLayerAver.cpp.

122 {
123  size_t np = afl.getNumberOfPanels();
124 
125  for (size_t i = 0; i < np; ++i)
126  {
127  lactCol(i) = 1.0;
128  lastLine(i) = afl.len[i];
129  }
130 
131  for (size_t i = 0; i < np; ++i)
132  for (size_t j = 0; j < np; ++j)
133  matr(i, j) = afl.getA(1, i, afl, j)[0];
134 
135 }//FillMatrixSelf(...)
virtual std::vector< double > getA(size_t p, size_t i, const Airfoil &airfoil, size_t j) const =0
Вычисление коэффициентов матрицы A для расчета влияния панели на панель
size_t getNumberOfPanels() const
Возврат количества панелей на профиле
Definition: Airfoil2D.h:139
std::vector< double > len
Длины панелей профиля
Definition: Airfoil2D.h:185
const Airfoil & afl
Definition: Boundary2D.h:76

Here is the call graph for this function:

Here is the caller graph for this function:

void BoundaryConstLayerAver::GetInfluenceFromSourceSheetAtRectPanelToVortex ( size_t  panel,
const Vortex2D vtx,
Point2D vel 
) const
overridevirtual

Вычисление влияния слоя источников конкретной прямолинейной панели на вихрь в области течения

Parameters
[in]panelномер панели профиля, от которой считается влияние
[in]vtxссылка на вихрь
[out]velссылка на вектор полученной скорости

Implements VM2D::Boundary.

Definition at line 353 of file Boundary2DConstLayerAver.cpp.

354 {
355  vel.toZero();
356 
357  const Point2D& posI = ptr.r();
358 
359  Point2D dj = afl.getR(panel + 1) - afl.getR(panel);
360  Point2D tauj = dj.unit();
361 
362  Point2D s = posI - afl.getR(panel);
363  Point2D p = posI - afl.getR(panel + 1);
364 
365  double a = VMlib::Alpha(p, s);
366 
367  double lambda;
368  if ((s.length2() > 1e-16) && (p.length2() > 1e-16))
369  lambda = VMlib::Lambda(p, s);
370  else
371  lambda = 0.0;
372 
373  vel += sheets.attachedSourceSheet(panel, 0) * (-a * tauj.kcross() + lambda * tauj);
374 }// GetInfluenceFromSourceSheetAtRectPanelToVortex(...)
const Point2D & getR(size_t q) const
Возврат константной ссылки на вершину профиля
Definition: Airfoil2D.h:101
double Alpha(const Point2D &p, const Point2D &s)
Вспомогательная функция вычисления угла между векторами
Definition: defs.cpp:259
auto length2() const -> typename std::remove_const< typename std::remove_reference< decltype(this->data[0])>::type >::type
Вычисление квадрата нормы (длины) вектора
Definition: numvector.h:383
double Lambda(const Point2D &p, const Point2D &s)
Вспомогательная функция вычисления логарифма отношения норм векторов
Definition: defs.cpp:265
Класс, опеделяющий двумерный вектор
Definition: Point2D.h:75
auto unit(P newlen=1) const -> numvector< typename std::remove_const< decltype(this->data[0]*newlen)>::type, n >
Вычисление орта вектора или вектора заданной длины, коллинеарного данному
Definition: numvector.h:399
Sheet sheets
Слои на профиле
Definition: Boundary2D.h:95
const double & attachedSourceSheet(size_t n, size_t moment) const
Definition: Sheet2D.h:109
numvector< T, n > & toZero(P val=0)
Установка всех компонент вектора в константу (по умолчанию — нуль)
Definition: numvector.h:524
const Airfoil & afl
Definition: Boundary2D.h:76
numvector< T, 2 > kcross() const
Геометрический поворот двумерного вектора на 90 градусов
Definition: numvector.h:507

Here is the call graph for this function:

Here is the caller graph for this function:

void BoundaryConstLayerAver::GetInfluenceFromSourcesToRectPanel ( size_t  panel,
const Vortex2D ptr,
ptrdiff_t  count,
std::vector< double > &  wakeRhs 
) const
overridevirtual

Вычисление влияния части подряд источников из области течения на прямолинейную панель для правой части

Вычисляет влияния части подряд идущих источников из области течения на прямолинейную панель для правой части

Parameters
[in]panelномер панели профиля, на которую считается влияние
[in]ptrуказатель на начало диапазона источников
[in]countдлина диапазона источников
[out]wakeRhsссылка на вектор полученных влияние для правой части СЛАУ

Implements VM2D::Boundary.

Definition at line 329 of file Boundary2DConstLayerAver.cpp.

330 {
331  double& velI = wakeRhs[0];
332 
333  const Point2D& posI0 = afl.getR(panel);
334  const Point2D& posI1 = afl.getR(panel + 1);
335 
336  for (size_t it = 0; it != count; ++it)
337  {
338  const Vortex2D& vt = ptr[it];
339  const Point2D& posJ = vt.r();
340  const double& gamJ = vt.g();
341 
342  Point2D s = posJ - posI0;
343  Point2D p = posJ - posI1;
344 
345  double lambda = VMlib::Lambda(p, s);
346 
347  velI -= gamJ * lambda;
348  }
349 }// GetInfluenceFromSourcesToRectPanel(...)
const Point2D & getR(size_t q) const
Возврат константной ссылки на вершину профиля
Definition: Airfoil2D.h:101
HD Point2D & r()
Функция для доступа к радиус-вектору вихря
Definition: Vortex2D.h:85
HD double & g()
Функция для доступа к циркуляции вихря
Definition: Vortex2D.h:93
double Lambda(const Point2D &p, const Point2D &s)
Вспомогательная функция вычисления логарифма отношения норм векторов
Definition: defs.cpp:265
Класс, опеделяющий двумерный вектор
Definition: Point2D.h:75
Класс, опеделяющий двумерный вихревой элемент
Definition: Vortex2D.h:51
const Airfoil & afl
Definition: Boundary2D.h:76

Here is the call graph for this function:

Here is the caller graph for this function:

void BoundaryConstLayerAver::GetInfluenceFromVInfToRectPanel ( std::vector< double > &  vInfRhs) const
overridevirtual

Вычисление влияния набегающего потока на прямолинейную панель для правой части

Вычисляет влияния набегающего потока на прямолинейную панель для правой части

Parameters
[out]vInfRhsссылка на вектор полученных влияние для правой части СЛАУ

Implements VM2D::Boundary.

Definition at line 405 of file Boundary2DConstLayerAver.cpp.

406 {
407  size_t np = afl.getNumberOfPanels();
408  vInfRhs.resize(np);
409 
410 #pragma omp parallel for default(none) shared(vInfRhs, np)
411  for (int i = 0; i < np; ++i)
412  vInfRhs[i] = afl.tau[i] & W.getPassport().physicalProperties.V0();
413 
414 }
PhysicalProperties physicalProperties
Структура с физическими свойствами задачи
Definition: Passport2D.h:296
size_t getNumberOfPanels() const
Возврат количества панелей на профиле
Definition: Airfoil2D.h:139
std::vector< Point2D > tau
Касательные к панелям профиля
Definition: Airfoil2D.h:182
const Passport & getPassport() const
Возврат константной ссылки на паспорт
Definition: World2D.h:222
Point2D V0() const
Функция скорости набегающего потока с учетом разгона
Definition: Passport2D.h:93
const World2D & W
Константная ссылка на решаемую задачу
Definition: Boundary2D.h:67
const Airfoil & afl
Definition: Boundary2D.h:76

Here is the call graph for this function:

Here is the caller graph for this function:

void BoundaryConstLayerAver::GetInfluenceFromVortexSheetAtRectPanelToVortex ( size_t  panel,
const Vortex2D vtx,
Point2D vel 
) const
overridevirtual

Вычисление влияния вихревых слоев (свободный + присоединенный) конкретной прямолинейной панели на вихрь в области течения ///.

Parameters
[in]panelномер панели профиля, от которой считается влияние
[in]vtxссылка на вихрь
[out]velссылка на вектор полученной скорости

Implements VM2D::Boundary.

Definition at line 377 of file Boundary2DConstLayerAver.cpp.

378 {
379  vel.toZero();
380 
381  const Point2D& posI = ptr.r();
382 
383  Point2D dj = afl.getR(panel + 1) - afl.getR(panel);
384  Point2D tauj = dj.unit();
385 
386  Point2D s = posI - afl.getR(panel);
387  Point2D p = posI - afl.getR(panel + 1);
388  double a = VMlib::Alpha(p, s);
389 
390  double lambda;
391  if ((s.length2() > 1e-16) && (p.length2() > 1e-16))
392  lambda = VMlib::Lambda(p, s);
393  else
394  lambda = 0.0;
395 
396  Point2D skos = -a * tauj.kcross() + lambda * tauj;
397 
398  vel += sheets.freeVortexSheet(panel, 0) * skos.kcross();
399  vel += sheets.attachedVortexSheet(panel, 0) * skos.kcross();
400 
401 }// GetInfluenceFromVortexSheetAtRectPanelToVortex(...)
const double & freeVortexSheet(size_t n, size_t moment) const
Definition: Sheet2D.h:99
const double & attachedVortexSheet(size_t n, size_t moment) const
Definition: Sheet2D.h:104
const Point2D & getR(size_t q) const
Возврат константной ссылки на вершину профиля
Definition: Airfoil2D.h:101
double Alpha(const Point2D &p, const Point2D &s)
Вспомогательная функция вычисления угла между векторами
Definition: defs.cpp:259
auto length2() const -> typename std::remove_const< typename std::remove_reference< decltype(this->data[0])>::type >::type
Вычисление квадрата нормы (длины) вектора
Definition: numvector.h:383
double Lambda(const Point2D &p, const Point2D &s)
Вспомогательная функция вычисления логарифма отношения норм векторов
Definition: defs.cpp:265
Класс, опеделяющий двумерный вектор
Definition: Point2D.h:75
auto unit(P newlen=1) const -> numvector< typename std::remove_const< decltype(this->data[0]*newlen)>::type, n >
Вычисление орта вектора или вектора заданной длины, коллинеарного данному
Definition: numvector.h:399
Sheet sheets
Слои на профиле
Definition: Boundary2D.h:95
numvector< T, n > & toZero(P val=0)
Установка всех компонент вектора в константу (по умолчанию — нуль)
Definition: numvector.h:524
const Airfoil & afl
Definition: Boundary2D.h:76
numvector< T, 2 > kcross() const
Геометрический поворот двумерного вектора на 90 градусов
Definition: numvector.h:507

Here is the call graph for this function:

Here is the caller graph for this function:

void BoundaryConstLayerAver::GetInfluenceFromVorticesToRectPanel ( size_t  panel,
const Vortex2D ptr,
ptrdiff_t  count,
std::vector< double > &  wakeRhs 
) const
overridevirtual

Вычисление влияния части подряд идущих вихрей из вихревого следа на прямолинейную панель для правой части

Вычисляет влияния части подряд идущих вихрей из вихревого следа на прямолинейную панель для правой части

Parameters
[in]panelномер панели профиля, на которую считается влияние
[in]ptrуказатель на начало диапазона вихрей
[in]countдлина диапазона вихрей
[out]wakeRhsссылка на вектор полученных влияние для правой части СЛАУ

Implements VM2D::Boundary.

Definition at line 302 of file Boundary2DConstLayerAver.cpp.

303 {
304  double& velI = wakeRhs[0];
305 
306  const Point2D& posI0 = afl.getR(panel);
307  const Point2D& posI1 = afl.getR(panel + 1);
308 
309 
310  for (size_t it = 0; it != count; ++it)
311  {
312  const Vortex2D& vt = ptr[it];
313 
314  const Point2D& posJ = vt.r();
315  const double& gamJ = vt.g();
316 
317  Point2D s = posJ - posI0;
318  Point2D p = posJ - posI1;
319 
320  double alpha = VMlib::Alpha(p, s);
321 
322  velI -= gamJ * alpha;
323  }
324 }// GetInfluenceFromVorticesToRectPanel(...)
const Point2D & getR(size_t q) const
Возврат константной ссылки на вершину профиля
Definition: Airfoil2D.h:101
HD Point2D & r()
Функция для доступа к радиус-вектору вихря
Definition: Vortex2D.h:85
HD double & g()
Функция для доступа к циркуляции вихря
Definition: Vortex2D.h:93
double Alpha(const Point2D &p, const Point2D &s)
Вспомогательная функция вычисления угла между векторами
Definition: defs.cpp:259
Класс, опеделяющий двумерный вектор
Definition: Point2D.h:75
Класс, опеделяющий двумерный вихревой элемент
Definition: Vortex2D.h:51
const Airfoil & afl
Definition: Boundary2D.h:76

Here is the call graph for this function:

Here is the caller graph for this function:

size_t Boundary::GetUnknownsSize ( ) const
inherited

Возврат размерности вектора решения

(без учета регуляризирующей переменной)

Returns
размерность вектора решения

Definition at line 71 of file Boundary2D.cpp.

72 {
73  return sheetDim * afl.getNumberOfPanels();
74 }
size_t getNumberOfPanels() const
Возврат количества панелей на профиле
Definition: Airfoil2D.h:139
size_t sheetDim
Размерность параметров каждого из слоев на каждой из панелей
Definition: Boundary2D.h:92
const Airfoil & afl
Definition: Boundary2D.h:76

Here is the call graph for this function:

Here is the caller graph for this function:

void BoundaryConstLayerAver::SolutionToFreeVortexSheetAndVirtualVortex ( const Eigen::VectorXd &  sol)
overridevirtual

Пересчет решения на интенсивность вихревого слоя и на рождаемые вихри на конкретном профиле

1) Приводит решение к интенсивности вихревого слоя и записывает его в sheets.freeVortexSheet:

  • если неизвестное — интенсивность вихря, то он "размазывается" по панели;
  • если неизвестное — интенсивность слоя, то она передается непосредственно.

2) Приводит интенсивность вихревого слоя к рождаемым вихрям, а также вычисляет их положения

Parameters
[in]solвектор решения СЛАУ

Implements VM2D::Boundary.

Definition at line 57 of file Boundary2DConstLayerAver.cpp.

58 {
59  Vortex2D virtVort;
60  Point2D midNorm;
61 
62  size_t np = afl.getNumberOfPanels();
63 
65 
67 
68  //Очистка и резервирование памяти
69  virtualWake.vecHalfGamma.clear();
70  virtualWake.vecHalfGamma.reserve(np * nVortPerPan);
71 
72  //Очистка и резервирование памяти
73  virtualWake.aflPan.clear();
74  virtualWake.aflPan.reserve(np * nVortPerPan);
75 
76  //Резервирование памяти
77  virtualWake.vtx.clear();
78  virtualWake.vtx.reserve(np * nVortPerPan);
79 
80  //Очистка и резервирование памяти
81  vortexBeginEnd.clear();
82  vortexBeginEnd.reserve(np);
83 
85 
86  std::pair<int, int> pair;
87 
88  for (size_t i = 0; i < np; ++i)
89  {
90  midNorm = afl.nrm[i] * delta;
91 
92  size_t NEWnVortPerPan = (size_t)std::max((int)std::ceil(fabs(sol(i) * afl.len[i]) / maxG), nVortPerPan);
93 
94  pair.first = (int)virtualWake.vtx.size();
95 
96 
97  Point2D dr = 1.0 / NEWnVortPerPan * (afl.getR(i + 1) - afl.getR(i));
98 
99  for (size_t j = 0; j < NEWnVortPerPan; ++j)
100  {
101  virtVort.r() = afl.getR(i) + dr * (j * 1.0 + 0.5) + midNorm;
102  virtVort.g() = sol(i) * afl.len[i] / NEWnVortPerPan;
103  virtualWake.vtx.push_back(virtVort);
104 
105  virtualWake.vecHalfGamma.push_back(0.5 * sol(i) * afl.tau[i]);
106  virtualWake.aflPan.push_back({ numberInPassport, i });
107  }
108 
109  pair.second = (int)virtualWake.vtx.size();
110  vortexBeginEnd.push_back(pair);
111  }
112 
113 
114  for (size_t j = 0; j < np; ++j)
115  sheets.freeVortexSheet(j, 0) = sol(j);
116 
117 }//SolutionToFreeVortexSheetAndVirtualVortex(...)
const double & freeVortexSheet(size_t n, size_t moment) const
Definition: Sheet2D.h:99
double delta
Расстояние, на которое рождаемый вихрь отодвигается от профиля
Definition: Passport2D.h:151
double maxGamma
Максимально допустимая циркуляция вихря
Definition: Passport2D.h:157
std::vector< Point2D > vecHalfGamma
Скорость вихрей виртуального следа конкретного профиля (равна Gamma/2) используется для расчета давле...
Definition: VirtualWake2D.h:70
const Point2D & getR(size_t q) const
Возврат константной ссылки на вершину профиля
Definition: Airfoil2D.h:101
HD Point2D & r()
Функция для доступа к радиус-вектору вихря
Definition: Vortex2D.h:85
HD double & g()
Функция для доступа к циркуляции вихря
Definition: Vortex2D.h:93
size_t getNumberOfPanels() const
Возврат количества панелей на профиле
Definition: Airfoil2D.h:139
std::vector< std::pair< size_t, size_t > > aflPan
Пара чисел: номер профиля и номер панели, на которой рожден виртуальный вихрь
Definition: VirtualWake2D.h:73
const size_t numberInPassport
Номер профиля в паспорте
Definition: Boundary2D.h:70
std::vector< Vortex2D > vtx
Список вихревых элементов
Класс, опеделяющий двумерный вектор
Definition: Point2D.h:75
std::vector< Point2D > tau
Касательные к панелям профиля
Definition: Airfoil2D.h:182
Sheet sheets
Слои на профиле
Definition: Boundary2D.h:95
const Passport & getPassport() const
Возврат константной ссылки на паспорт
Definition: World2D.h:222
Класс, опеделяющий двумерный вихревой элемент
Definition: Vortex2D.h:51
int minVortexPerPanel
Минимальное число вихрей, рождаемых на каждой панели профииля
Definition: Passport2D.h:154
VirtualWake virtualWake
Виртуальный вихревой след конкретного профиля
Definition: Boundary2D.h:85
std::vector< Point2D > nrm
Нормали к панелям профиля
Definition: Airfoil2D.h:177
WakeDiscretizationProperties wakeDiscretizationProperties
Структура с параметрами дискретизации вихревого следа
Definition: Passport2D.h:299
std::vector< double > len
Длины панелей профиля
Definition: Airfoil2D.h:185
std::vector< std::pair< int, int > > vortexBeginEnd
Номера первого и последнего вихрей, рождаемых на каждой панели профиля (формируется после решения СЛА...
Definition: Boundary2D.h:82
const World2D & W
Константная ссылка на решаемую задачу
Definition: Boundary2D.h:67
const Airfoil & afl
Definition: Boundary2D.h:76

Here is the call graph for this function:

Here is the caller graph for this function:

Member Data Documentation

const Airfoil& VM2D::Boundary::afl
inherited

Константная ссылка на профиль
инициализируется автоматом в конструкторе

Definition at line 76 of file Boundary2D.h.

int VM2D::Boundary::minVortexPerPanel
inherited

Минимальное число вихрей, рождаемых на панели профиля и формирующих виртуальный вихревой след

Definition at line 79 of file Boundary2D.h.

const size_t VM2D::Boundary::numberInPassport
protectedinherited

Номер профиля в паспорте

Definition at line 70 of file Boundary2D.h.

Sheet VM2D::Boundary::oldSheets
inherited

Слои на профиле с предыдущего шага

Definition at line 98 of file Boundary2D.h.

size_t VM2D::Boundary::sheetDim
inherited

Размерность параметров каждого из слоев на каждой из панелей

Указывает, сколькими числами задается интенсивность каждого из слоев на каждой панели:

  • 1 — одно число — задается только среднее значение;
  • 2 — два числа — задается среднее значение и "наклон".

Definition at line 92 of file Boundary2D.h.

Sheet VM2D::Boundary::sheets
inherited

Слои на профиле

Definition at line 95 of file Boundary2D.h.

VirtualWake VM2D::Boundary::virtualWake
inherited

Виртуальный вихревой след конкретного профиля

Definition at line 85 of file Boundary2D.h.

std::vector<std::pair<int, int> > VM2D::Boundary::vortexBeginEnd
inherited

Номера первого и последнего вихрей, рождаемых на каждой панели профиля (формируется после решения СЛАУ)

Definition at line 82 of file Boundary2D.h.

const World2D& VM2D::Boundary::W
protectedinherited

Константная ссылка на решаемую задачу

Definition at line 67 of file Boundary2D.h.


The documentation for this class was generated from the following files: