Распределение скоростей и давлений в движущейся жидкости

[tehnolog]

Здравствуйте, форумчане. Столкнулся с задачей по гидродинамике. Необходимо определить распределение скоростей и давлений в несжимаемой, идеальной жидкости, которая подается в сужающийся усеченный конус, при этом жидкость предварительно закручена, и до входа в конус имеет постоянную угловую скорость. В процессе прохождения по конусу жидкость не только увеличивает скорость поступательного движения, но и вращательного (по закону сохранения момента импульса)
Для решения поставленной задачи решил исходить из уравнений Эйлера для жидкости:

$$\frac{dV_x}{d\tau}=\frac{\partial{V_x}}{\partial{x}}\cdot{V_x}+\frac{\partial{V_x}}{\partial{y}}\cdot{V_y}+\frac{\partial{V_x}}{\partial{z}}\cdot{V_z}=-\frac{1}{\rho}\frac{\partial{p}}{\partial{x}}$$

$$\frac{dV_y}{d\tau}=\frac{\partial{V_y}}{\partial{x}}\cdot{V_x}+\frac{\partial{V_y}}{\partial{y}}\cdot{V_y}+\frac{\partial{V_y}}{\partial{z}}\cdot{V_z}=-\frac{1}{\rho}\frac{\partial{p}}{\partial{y}}\qquad\qquad (1)$$

$$\frac{dV_z}{d\tau}=\frac{\partial{V_z}}{\partial{x}}\cdot{V_x}+\frac{\partial{V_z}}{\partial{y}}\cdot{V_y}+\frac{\partial{V_z}}{\partial{z}}\cdot{V_z}=-\frac{1}{\rho}\frac{\partial{p}}{\partial{z}}$$

Каждому элементарному объему жидкости соответствуют определенные значения $x,y, z$ в декартовой системе координат, и $r, \beta, y$ в новой (не знаю, как правильно ее назвать (конической(?)) системе координат. Из рисунка можно видеть: $x=r\cdot\cos{\beta}; z=r\cdot\sin{\beta}$. Тогда:

$$V_x=\frac{dx}{d\tau}=\frac{dr}{d\tau}\cdot\cos{\beta}-r\cdot\sin{\beta}\cdot\frac{d\beta}{d\tau}$$

$$V_y=\frac{dy}{d\tau}\qquad\qquad (2)$$

$$V_z=\frac{dz}{d\tau}=\frac{dr}{d\tau}\cdot\sin{\beta}+r\cdot\cos{\beta}\cdot\frac{d\beta}{d\tau}$$

Где $\frac{dr}{d\tau}=V_r; \frac{d\beta}{d\tau}=\omega$. Находим частные производные компонент скоростей по координатам, как функции $r, \beta, y$:

$$\frac{\partial{V_x}}{\partial{x}}=\frac{\partial{V_x}}{\partial{r}}\cdot\frac{\partial{r}}{\partial{x}}+\frac{\partial{V_x}}{\partial{\beta}}\cdot\frac{\partial{\beta}}{\partial{x}}$$

$$\frac{\partial{r}}{\partial{x}}=\frac{\partial{(\sqrt{x^2+z^2})}}{\partial{x}}=\frac{x}{\sqrt{x^2+z^2}}=\cos{\beta}$$

$$\beta=\arccos\frac{x}{r}=\arccos\frac{x}{\sqrt{x^2+z^2}}$$

$$\frac{\partial\beta}{\partial x}=-\frac{1}{\sqrt{1-\frac{x^2}{x^2+z^2}}}\cdot\frac{\sqrt{x^2+z^2}- \frac{x^2}{\sqrt{x^2+z^2}} }{x^2+z^2}=-\frac{z}{r^2}=-\frac{\sin\beta}{r}$$

Тогда

$$\frac{\partial{V_x}}{\partial{x}}=\frac{\partial{V_r}}{\partial{r}}\cdot\cos^2\beta-(\omega+r\cdot\frac{\partial{\omega}}{\partial{r}})\cdot\sin\beta\cdot\cos\beta+\frac{V_r}{r}\cdot\sin^2\beta+\omega\cdot\cos\beta\cdot\sin\beta\qquad (3)$$

Аналогично:

$$\frac{\partial{V_x}}{\partial{z}}=\frac{\partial{V_x}}{\partial{r}}\cdot\frac{\partial{r}}{\partial{z}}+\frac{\partial{V_x}}{\partial{\beta}}\cdot\frac{\partial{\beta}}{\partial{z}}$$

$$\frac{\partial{V_x}}{\partial{z}}=\frac{\partial{V_r}}{\partial{r}}\cdot\cos\beta\cdot\sin\beta-(\omega+r\cdot\frac{\partial{\omega}}{\partial{r}})\cdot\sin^2\beta-\frac{V_r}{r}\cdot\sin\beta\cdot\cos\beta-\omega\cdot\cos^2\beta \qquad (4)$$

$$\frac{\partial V_x}{\partial y}=\frac{\partial V_r}{\partial y}\cdot\cos\beta-r\cdot\frac{\partial \omega}{\partial y}\cdot \sin\beta \qquad (5) $$

Используя первое из уравнений Эйлера (1), с учетом (2-5) найдем проекции ускорений на ось $x$. Для угла $\beta = 0$:
$\frac{1}{\rho}F_x=\frac{\partial V_r}{\partial r}\cdot V_r+\frac{\partial V_r}{\partial y}\cdot V_y-\omega^2r$ Это ускорение обусловлено центростремительной силой, которая в силу симметрии не зависит от угла $\beta$. Поэтому можно записать:

$$\frac{\partial V_r}{\partial r}\cdot V_r+\frac{\partial V_r}{\partial y}\cdot V_y-\omega^2r=-\frac{1}{\rho}\frac{\partial p}{\partial r}\qquad (6)$$

Для угла $\beta=\frac{\pi}{2}$:
$\frac{1}{\rho}F_x=2\omega \cdot V_r+r\cdot V_r\frac{\partial \omega}{\partial r}+r\cdot V_y\frac{\partial \omega}{\partial y}$ Очевидно, что это ускорение равно нулю (т.к. давление для постоянного $r$ и $y$ с изменением угла $\beta$ не меняется). Поэтому

$$2\omega V_r+r\cdot V_r\frac{\partial \omega}{\partial r}+r\cdot V_y\frac{\partial \omega}{\partial y}=0\qquad (7)$$

Далее, $\frac{\partial V_y}{\partial x}=\frac{\partial V_y}{\partial r}\cdot\frac{\partial r}{\partial x};\quad \frac{\partial V_y}{\partial z}=\frac{\partial V_y}{\partial r}\cdot\frac{\partial r}{\partial z}$ После подстановки соответствующих производных во второе уравнение Эйлера (1) упрощений и сокращений, получаем:

$$\frac{\partial V_y}{\partial y}\cdot V_y+\frac{\partial V_y}{\partial r}\cdot V_r=-\frac{1}{\rho}\frac{\partial p}{\partial y}\qquad (8)$$

Принимая, что давление является непрерывной функцией от $r$ и $y$, и существуют его частные производные до второго порядка справедливо соотношение: $\frac{\partial^2 p}{\partial r \partial y}=\frac{\partial^2 p}{\partial y \partial r}$ Тогда дифференцируя левую часть выражения (6) по $y$, а левую часть (8) по $r$, и приравнивая получаем:

$$\frac{\partial^2 V_r}{\partial r \partial y}\cdot V_r+\frac{\partial V_r}{\partial r}\cdot\frac{\partial V_r}{\partial y}+\frac{\partial^2V_r}{\partial y^2}\cdot V_y+\frac{\partial V_r}{\partial y}\cdot\frac{\partial V_y}{\partial y}-2\omega r\cdot\frac{\partial \omega}{\partial y}=$$

$$=\frac{\partial^2 V_y}{\partial y \partial r}\cdot V_y+\frac{\partial V_y}{\partial y}\cdot\frac{\partial V_y}{\partial r}+\frac{\partial^2V_y}{\partial r^2}\cdot V_r+\frac{\partial V_y}{\partial r}\cdot\frac{\partial V_r}{\partial r}\quad (9)$$

Уже есть два дифференциальных уравнения и три неизвестных функции, подлежащие определению: $V_y=f_1(r, y) ; V_r=f_2(r, y); \omega=f_3(r, y)$ Есть еще одно условие описываемое уравнением: $\operatorname{div}V = \frac{\partial V_x}{\partial x}+\frac{\partial V_y}{\partial y}+\frac{\partial V_z}{\partial z}=0$, которое после соответствующих подстановок для нашего случая приобретает вид:

$$\frac{\partial V_r}{\partial r}+\frac{\partial V_y}{\partial y}+\frac{V_r}{r}=0\qquad (10)$$

От смешанных производных в уравнении (9) можно избавиться, используя уравнение (10). Возникают следующие вопросы: 1) Хватает ли уравнений 7), 9) и 10) для определения искомых функций, при надлежащем выборе граничных условий; 2) Как выбрать дополнительные граничные условия, кроме этих: $V_r(r, y_0) = 0;\quad \omega(r, y_0)=\omega_0$; 3) Можно ли решить эту систему уравнений средствами маткад 15, если нет, то что можно применить?