Устойчивость первых двух собственных форм

Рассмотрим следующее уравнение: $$EI \frac{\partial^4 y(x, t)}{\partial x^4} + p \frac{\partial^2 y(x, t)}{\partial x^2} + m \frac{\partial^2 y(x, t)}{\partial t^2} = 0$$ вместе с граничными условиями $$\left . \frac{\partial y(x, t)}{\partial x} \right \vert_{x = 0} = \left . \frac{\partial^3 y(x, t)}{\partial x^3} \right \vert_{x = 0} = \left . \frac{\partial^2 y(x, t)}{\partial x^2} \right \vert_{x = l} = \left . \frac{\partial^3 y(x, t)}{\partial x^3} \right \vert_{x = l} =0$$

Заменой $x = l \phi$, $t = \sqrt{\frac{m}{EI}} l^2 \tau$ мы переходим к безразмерной форме уравнения: $$\frac{\partial^4 y(\phi, \tau)}{\partial \phi^4} + b \frac{\partial^2 y(\phi, \tau)}{\partial \phi^2} + \frac{\partial^2 y(\phi, \tau)}{\partial \tau^2} = 0,$$ где $b = \frac{pl^2}{EI}$ и граничные условия принимают вид: $$\left . \frac{\partial y(\phi, \tau)}{\partial \phi} \right \vert_{\phi = 0} = \left . \frac{\partial^3 y(\phi, \tau)}{\partial \phi^3} \right \vert_{\phi = 0} = \left . \frac{\partial^2 y(\phi, \tau)}{\partial \phi^2} \right \vert_{\phi = 1} = \left . \frac{\partial^3 y(\phi, \tau)}{\partial \phi^3} \right \vert_{\phi = 1} =0.$$

Теперь в качестве решения подставим $y(\phi, \tau) = X_1(\phi) T_1(\phi, \tau) + X_2(\phi) T_2 (\phi, \tau)$, где функции $X_1$, $X_2$, $T_1$ и $T_2$ связаны с первыми двумя собственными формами данного уравнения: $$\sum_{k=1}^{2} X_k^{IV}T_k + b \sum_{k=1}^{2}X_k^{II} T_k + \sum_{k=1}^{2}X_k \ddot{T}_k = 0$$

Умножим скалярно это соотношение на $X_1$ и $X_2$: $$\begin{equation} \left \lbrace \begin{array}{ccc} \ddot{T}_1 + T_1 \cdot \langle X_1^{IV} + b X_1^{II} ,\; X_1 \rangle + T_2 \cdot \langle X_2^{IV} + b X_2^{II} ,\; X_1 \rangle & = &0\\ \ddot{T}_2 + T_1 \cdot \langle X_1^{IV} + b X_1^{II} ,\; X_2 \rangle + T_2 \cdot \langle X_2^{IV} + b X_2^{II} ,\; X_2 \rangle & = &0 \end{array} \right . \end{equation}$$

Примем $a_{ij} = \langle X_j^{IV}+b X_j^{II} , X_i \; \rangle$. Тогда систему можно записать в виде: $$\ddot{T} = A T,$$ где $T = (T_1, T_2)$, а $A = \left ( \begin{array}{cc} -a_{11} & -a_{12}\\ -a_{21} & -a_{22} \end{array} \right )$. Вопросы асимптотического поведения (а, следовательно, и устойчивости) определяются решениями следующего характеристического уравнения: $$\det{(A - \lambda^2 I)} = 0$$

Теперь уточним, чему же именно равны коэффициенты матрицы $A$. Напомним, какому соотношению удовлетворяет собственная форма $X(\phi)$. В методе разделения переменных мы подставляем функцию вида $y(\phi, \tau) = X(\phi) T(\tau)$ в исходное уравнение и получаем $$X^{IV} T + b X^{II} T + X \ddot{T} = 0$$ $$\frac{X^{IV}(\phi) + b X^{II}(\phi)}{X(\phi)} = -\frac{\ddot{T}(\tau)}{T(\tau)} = \lambda,$$

откуда $X^{IV}(\phi) + b X^{II}(\phi) = \lambda X(\phi)$. Но тогда это означает, что $a_{ij} = \langle X^{IV}_j + b X^{II}_j, X_i \rangle = \langle \lambda_j X_j, X_i \rangle = \lambda_j \delta_{ij}$, где $\delta_{ij}$ — дельта Кронекера, и на самом деле $$A = \left ( \begin{array}{cc} -\lambda_1 & 0\\ 0 & -\lambda_2 \end{array} \right ).$$

Система уравнений распадается на два независимых уравнения $\ddot{T}_1 = -\lambda_1 T_1$ и $\ddot{T}_2 = -\lambda_2 T_2$. Вопрос асимптотики решений для них очевиден: если $\lambda_i > 0$, то в решении получаются гармонические функции, которые ограничены при всех $\tau \in \mathbb{R}$; если же $\lambda_i < 0$, то решение есть суперпозиция $\cosh(\sqrt{-\lambda_i} \tau)$ и $\sinh(\sqrt{-\lambda_i} \tau)$, которая ограничена только в случае тривиального решения.

Всё, к чему теперь сводится задача — определить знак собственных чисел $\lambda_1$ и $\lambda_2$. Проанализируем какими могут быть корни характеристического уравнения в зависимости от $\lambda$: $$p^4 + b p^2 - \lambda = 0,$$ $$\left (p^2 + \frac{b}{2} \right)^2 = \lambda + \frac{b^2}{4}.$$ Из физических соображений $b > 0$, поэтому можно выделить следующие случаи:

• если $\lambda < -\frac{b^2}{4}$, то корни имеют вид $\pm r e^{\pm i \theta}$, $r>0$, $\theta \in \lbrack 0, \frac{\pi}{2} \rbrack$

• если $-\frac{b^2}{4} < \lambda < 0$, то корни имеют вид $\pm i \omega$ и $\pm i \theta$, $\omega > 0$, $\theta > 0$

• если $\lambda > 0$, то корни имеют вид $\pm i\omega$ и $\pm \alpha$, $\alpha > 0$, $\beta > 0$

Для устойчивости достаточно выяснить, возможно ли разрешить задачу на собственные числа при $\lambda > 0$. Если определитель системы для решения граничных условий всегда будет отличен от нуля, то это точно свидетельствует о том, что все собственные числа отрицательны, а значит система неустойчива.

Численное определение области устойчивости

В качестве аппроксимации области устойчивости воспользуемся следующим методом. Вместо настоящих первых двух мод будем брать их полиномиальные аппроксимации $\hat{X}_i$, а матрица $A$ будет состоять из элементов $a_{ij} = \langle \hat{X}^{IV}_j + b \hat{X}^{II}_j, \hat{X}_i \rangle$. Зная собственные числа матрицы $A$ для системы $\ddot{T}=AT$, мы можем ответить на вопрос об устойчивости. А именно, если матрица $A$ имеет положительное собственное число $p$, то есть и положительное собственное число $\lambda = \sqrt{p}$ в системе дифференциальных уравнений $\ddot{T} = AT$, что приведёт к неустойчивости решений.

Так как полиномиальная аппроксимация зависит только от граничных условий, то просто повторим поиск первых нескольких собственных форм:

Так как первое решение вообще почти никакого отношения к реальным модам не имеет, то мы его можем отбросить и вычислить производные только для $\hat{X}_1$ и $\hat{X}_2$:

В итоге это приведёт к следующей аппроксимации матрицы $A$ (которая сейчас зависит от параметра $b$):

Изобразим на графике след и определитель матрицы $A$ в зависимости от параметра $b$:

Мы наблюдаем, что до некоторого критического значения $b$ собственные числа матрицы $A(b)$ либо одного знака, либо комплексные, но их сумма имеет тот же знак, что и след матрицы $A(b)$, а он отрицателен. Это гарантирует отрицательность вещественных частей собственных чисел матрицы $A(b)$.

Смена устойчивости происходит при $b \approx 9.92$, после чего определитель меняет знак и собственные числа точно разных знаков, что приводит к неустойчивости.

Однако, есть ещё одна смена устойчивости при $b \approx 40.6$ — после этого значения двухмодовое решение снова становится устойчивым: