Быстрый расчет ударных аэродинамических сил

[Petr76]

Здравствуйте. Возникла проблема в данном вопросе.

1. Общий случай.
Суть задачи - быстро рассчитать аэродинамические силы и моменты сил воздействия на абсолютно твердое герметичное тело произвольной формы. Следует игнорировать основные уравнения газодинамики и воспользоваться лишь кинетическим давлением. Иначе говоря, существует только простая зависимость давления, примененного к элементу поверхности, от скорости воздушного потока и его направления.
Для конкретной задачи эта зависимость описана так:

$$\vec p=\frac{-\rho u^2}2\left(a\vec n\cos ^2(\hat{\vec n,\vec u})+b\vec n\cos(\hat{\vec n,\vec u})+c\frac{\vec u}u\right)$$

где $\vec p$ - вектор давления, приложенного к элементу поверхности;
$\vec n$ - нормаль элемента поверхности;
$\rho$ - плотность окружающего газа (она постоянна для всей поверхности);
$\vec u$ - скорость элемента поверхности в пространстве (считается что окружающий газ всегда остается неподвижным);
$u$ - модуль этой скорости;
$a$, $b$ и $c$ - простые коэффициенты, характеризующие величину прямого, обратного и скользящего давления, об этом ниже.
Основное условие данного уравнения: $\vec p\cdot\vec u\leqslant 0$. Другими словами, сила, воздействующая на элемент поверхности, не должна увеличивать модуль скорости этого элемента (условие диссипативных сил).

Безразмерные коэффициенты уравнения задаются следующим образом:
$a,b,c\in[0,1]$ - ограничение, накладываемое их физическим смыслом.
$a+b+c=1$ - максимальное прямое давление должно быть правильным (когда $\vec u\uparrow\uparrow\vec n$).
Величина $c$ задает давление при скользящем воздушном потоке (когда $\vec u\perp\vec n$).
Выражение $b+c-a \in[0,1]$ характеризует величину обратного давления (когда $\vec u\uparrow\downarrow\vec n$).
К примеру, можно задать следующее: $a = 0.34;b = 0.64;c = 0.02$.

2. Случай поступательного движения.
При такой поставки задачи все довольно просто раскладывается. Зададим величины $\vec v$ - линейную скорость тела в пространстве и $\vec r$ - вектор от центра масс объекта до элемента его поверхности.
Очевидно, что поскольку тело не вращается, все элементы его поверхности движутся с одинаковой скоростью: $\vec u(\vec r)=\vec v$.
В этом случае, уравнение давления после замены косинусов будет выглядеть следующим образом:

$$\vec p=\frac{-\rho v^2}2\left(a\vec n\frac{(\vec n\cdot\vec v)^2}{v^2}+b\vec n\frac{\vec n\cdot\vec v}v+c\frac{\vec v}v\right)$$

Далее, избавляемся от отношений:

$$\vec p=\frac{-\rho}2\left(a\vec n(\vec n\cdot\vec v)^2+b v\vec n(\vec n\cdot\vec v)+c v\vec v\right)$$

Такое выражение легко разложить на 2 части - то, что постоянно и зависит только от формы объекта (величины $\vec n,a,b,c$) и то, что может меняться (величины $\vec v,\rho$). Таким образом, к примеру, для оси $x$ данное выражение можно представить в виде:

$$\vec{p_x}=\rho V\cdot P_x$$

где $V$ - 9 мерный вектор комбинаций компонент скорости тела;
$P_x$ - 9 мерный вектор, характеризующий воздействие скорости тела по оси $x$.
Тоже самое верно для размерностей $y$ и $z$. Вектор $V$ выглядит следующим образом:

$$V=\left\{\vec{v_x}^2; \vec{v_y}^2; \vec{v_z}^2; \vec{v_x}\vec{v_y}; \vec{v_x}\vec{v_z}; \vec{v_y}\vec{v_z}; \vec{v_x}v; \vec{v_y}v; \vec{v_z}v\right\}$$

Вектора $P_x,P_y,P_z$ я тут расписывать не буду, т.к. это все займет довольно много места. Обозначу лишь, что $\vec P=\left\{P_x;P_y;P_z\right\}$.
Далее, рассчитаем значение для момента давления, ведь потом нужно будет получить момент аэродинамических сил:

$$\vec{lp}=(\rho V\cdot\vec P)\times\vec r=\rho V\cdot(\vec P\times\vec r)=\rho V\cdot\vec{LP}$$

где $\vec{LP}$ - очередной вектор 9-мерных векторов (Конечно, можно назвать это матрицей, но тогда теряется наглядность).
Затем, для давления и его момента выполняется интегрирование по всей площади тела:

$$\vec f=\iint\limits_{S}\rho V\cdot\vec P\,d\vec r = \rho V\cdot\vec F$$

$$\vec l=\iint\limits_{S}\rho V\cdot\vec{LP}\,d\vec r = \rho V\cdot\vec L$$

где $\vec F$ - проинтегрированный по поверхности вектор$\vec P$;
$\vec L$ - проинтегрированный по поверхности вектор$\vec{LP}$;
$\vec f$ - общая аэродинамическая сила, воздействующая на тело и $\vec l$ - момент этой силы.

Таким образом получается, что эту силу и момент силы для тела любой формы можно охарактеризовать всего 54-мя числами, что можно довольно быстро подсчитать. Этот алгоритм уже был осуществлен в виде конечной функции в программе. Его работа была полностью проверена и подтверждена. Возможно, тут есть ошибки в оформлении формул, но последовательность действий ясна.

3. Случай вращающегося тела, с неподвижным центром масс.
А вот тут у меня и возникла основная проблема, ведь теперь появилась зависимость скорости элемента поверхности от расстояния: $\vec u(\vec r)=\vec w\times\vec r$. Здесь $\vec w$ - вектор угловой скорости, измеряемый в радианах. Также, поскольку центр масс не движется, в формуле отсутствует линейная скорость тела $v$.

Я начал расписывать уравнение давления:

$$\vec p=\frac{-\rho(\vec w\times\vec r)^2}2\left(a\vec n\frac{(\vec n\cdot(\vec w\times\vec r))^2}{(\vec w\times\vec r)^2}+b\vec n\frac{\vec n\cdot(\vec w\times\vec r)}{|\vec w\times\vec r|}+c\frac{\vec w\times\vec r}{|\vec w\times\vec r|}\right)$$

Сократил его до следующего:

$$\vec p=\frac{-\rho}2\Bigl(a\vec n\bigl(\vec n\cdot(\vec w\times\vec r)\bigr)^2+|\vec w\times\vec r|\Bigl(b\vec n\bigl(\vec n\cdot(\vec w\times\vec r)\bigr)+c(\vec w\times\vec r)\Bigr)\Bigr)$$

И дальше просто разложить это выражение на компоненты не получается, т.к. в формуле присутствует модуль от угловой скорости $|\vec w\times\vec r|$, который нельзя просто так выразить как произведение постоянных и изменяющихся величин.
Расписать модуль в виде квадратного корня и вложить внутрь него все постоянные величины также невозможно, т.к. при этом теряется знак выражения, что мешает дальнейшему интегрированию.

Хотелось бы узнать в каком направлении стоит двигаться, ведь очевидно, что преобразование угловой скорости в силы существует, и его параметры постоянны. Я могу найти эти параметры косвенно, например путем аппроксимации уравнения силы по некоторым ключевым направлениям угловой скорости, а ее квадрат вынести за пределы искомых векторов, но это лишь приближенный метод расчета и он довольно медленный, ведь придется много раз интегрировать по всей площади. Может, стоить упростить исходную формулу давления, заменив ее такой, где не появиться модуль? Я пытался подобрать такое выражение но мне так и не удалось соблюсти условие диссипации сил.

[amon]

Petr76 в сообщении #1125301 писал(а):

где $\vec p$ - вектор давления

Я, видать, отстал от жизни и пропустил тот момент, когда давление стало вектором...

[manul91]

amon в сообщении #1125303 писал(а):

Я, видать, отстал от жизни и пропустил тот момент, когда давление стало вектором...

м.б. модуль вектора $\vec{p}$ это давление, а направление $\vec{p}$ ортогонально элементу поверхности dS? (в рассчеты не вникал..)