Обсужденный в п. 4.2 метод определения собственных частот колеблющихся систем обычно используется только в тех случаях, когда найти корни характеристического уравнения не представляет труда. Здесь также возможно применение различных численных методов *, но они обычно эффективнее в случае систем с большим числом \”степеней свободы. Обсуждаемый в данном параграфе подход иногда называют методом степенных рядов или методом Стодолы-Вианелло, но, как правило, его именуют просто итерационным методом. Этот подход удобно применять для работы с матрицами невысокого порядка, используя при расчетах логарифмическую линейку или настольный калькулятор, но решения больших задач следует программировать, чтобы проводить вычисления на цифровых ЭВМ.

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

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

Для того чтобы построить итерационный процесс, не связывая его с конкретным приложением к колебаниям, начнем с рассмотрения задачи на собственные значения, представленной в следующей стандартной форме:
\[
\mathbf{A X}_{M_{i}}=\lambda_{i} \mathbf{X}_{M_{i}},
\]

совпадающей с уравнением (э) в п. 4.2. Первое приближение для одного из собственных значений $\lambda_{i}$ можно найти, подставив пробный собственный вектор $(\mathbf{X})_{1}$ в левую и правую части уравнения (4.99) и решив результирующее соотношение по $\lambda_{i}$. Для этого обозначим через вектор (Y) $)_{1}$ стоящее в левой части соотношения произведение матрицы А и выбранного пробного значения собственного вектора-столбца ( $\mathbf{X})_{1}$ :
\[
(\mathbf{Y})_{\mathbf{1}}=\mathbf{A}(\mathbf{X})_{\mathbf{1}} .
\]

Если вектор $(\mathbf{X})_{1}$ не является точным значением собственного вектора, то после подстановки его в уравнение (4.99) последнее будет удовлетворяться только приближенно:
\[
(\mathbf{Y})_{1} \approx \lambda_{i}(\mathbf{X})_{1} .
\]

Первое приближение для собственного значения $\left(\lambda_{i}\right)_{1}$ можно получить, разделив любую из $n$ компонент вектора (Y) $)_{1}$ на соответствующую компоненту вектора $(\mathbf{X})_{1}$. Отметим, что если бы удалось установить точное значение собственного вектора, то все подобные отношения были бы равны между собой. Обозначим все возможные значения, получающиеся при таком делении, в виде
\[
\left(\lambda_{i}\right)_{1}=\left(y_{j}\right)_{1} /\left(x_{j}\right)_{1},
\]

где $1 \leqslant j \leqslant n$.
Прежде чем перейти ко второму шагу итерации, обычно проводят нормирование вектора (Y) $)_{1}$ тем или иным способом, например, разделив все его компоненты на первую или последнюю компоненту. В общем случае разделим вектор ( $\mathbf{Y})_{1}$ на произвольную постоянную и результат возьмем в качестве второго приближения для вектора
\[
(\mathbf{X})_{2}=(\mathbf{Y})_{1} / b_{1} .
\]

Этот вектор умножаем слева на матрицу $\mathbf{A}$, что дает вектор
\[
(\mathbf{Y})_{2}=\mathbf{A}(\mathbf{X})_{2} .
\]

Далее вычисляем второе приближение для собственного значения
\[
\left(\lambda_{i}\right)_{2}=\left(y_{j}\right)_{2} /\left(x_{j}\right)_{2} .
\]

Затем вектор (Y) $)_{2}$ масштабируем путем деления на произвольную постоянную $b_{2}$, что дает третье приближение для вектора
\[
(\mathbf{X})_{3}=(\mathbf{Y})_{2} / b_{2} .
\]

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

На $k$-м шаге итерации описанным выше процедурам соответствуют следующие рекуррентные уравнения:
\[
\begin{array}{c}
(\mathbf{Y})_{k}=\mathbf{A}(\mathbf{X})_{k} ; \\
\left(\lambda_{i}\right)_{k}=\left(y_{j}\right)_{k} /\left(x_{j}\right)_{k} ; \\
(\mathbf{X})_{k+1}=(\mathbf{Y})_{k} / b_{k},
\end{array}
\]

где $b_{k}$ – произвольно выбранное число. Ниже будет показано, что в результате повторяющейся процедуры вычислений по этим уравнениям процесс сходится к численно наибольшим собственному значению и собственному вектору. Поэтому для ускорения сходимости первое приближение для вектора ( $\mathbf{X})_{1}$ следует выбирать таким, чтобы оно численно соответствовало наибольшему собственному значению.

Для того чтобы доказать, что процесс итераций сходится к наибольшему собственному значению, разложим первое приближение для вектора ( $\mathbf{X})_{1}$ по действительным значениям собственных векторов системы
\[
(\mathbf{X})_{1}=\sum_{i=1}^{n} a_{i} \mathbf{X}_{M i}
\]

где $a_{1}, a_{2}, \ldots, a_{i}, \ldots, a_{n}$ – масштабирующие множители. Возможность подобного представления зависит от существования линейно независимых (хотя и неизвестных) собственных векторов, каждый из которых удовлетворяет рассматриваемым колеблющимся системам. Здесь также предполагается, что соответствующие собственные значения располагаются в порядке убывания, при этом имеется только одно наибольшее собственное значение
\[
\lambda_{1}>\lambda_{2} \geqslant \cdots \geqslant \lambda_{i} \geqslant \cdots \geqslant \lambda_{n} .
\]

Подставляя представление (з) в соотношение (а) и используя уравнение (4.99), можно получить следующее разложение для вектора (Y) $)_{1}$ по точным значениям собственных векторов:
\[
(\mathbf{Y})_{1}=\sum_{i=1}^{n} a_{i} \mathbf{A} \mathbf{X}_{M i}=\sum_{i=1}^{n} a_{i} \lambda_{i} \mathbf{X}_{M i} .
\]

На втором шаге итерации имеем
\[
(\mathbf{Y})_{2}=\mathbf{A}(\mathbf{X})_{2}=\mathbf{A}\left(\mathbf{Y}_{1}\right) b_{\mathbf{1}} .
\]

Подставляя представление (к) в выражение (л) и вновь используя уравнение (4.99), найдем выражение
\[
(\mathbf{Y})_{2}=\sum_{i=1}^{n} a_{i} \lambda_{i} \mathbf{A} \mathbf{X}_{M i} / b_{i}=\sum_{i=1}^{n} a_{i} \lambda_{i}^{2} \mathbf{X}_{M i} / b_{1} .
\]

На $k$-м шаге итерации получаем следующее выражение для вектора:
\[
(\mathbf{Y})_{k}=\sum_{i=1}^{n} a_{i} \lambda_{i}^{k} \mathbf{X}_{M i} / b_{1} b_{2} \ldots b_{k-1}
\]

в которое собственное значение $\lambda_{i}$ входит в $k$-й степени. Вынося за скобку общий множитель $\lambda_{1}^{k}$ и выделяя отдельно слагаемое, запишем содержащее $\mathbf{X}_{M 1}$ :
\[
(\mathbf{Y})_{k}=\lambda_{1}^{k}\left[a_{1} \mathbf{X}_{M 1}+\sum_{i=2}^{n} a_{i}\left(\lambda_{i} / \lambda_{1}\right)^{k} \mathbf{X}_{M i}\right] \mid b_{1} b_{2} \ldots b_{k-1} .
\]

В соответствии с выбранным порядком расположения собственных значений [см. выражение (и)] можно заключить, что величина $\left(\lambda_{i} / \lambda_{1}\right)^{k}$ стремится к нулю при увеличении номера $k$, поэтому можно записать
\[
(\mathbf{Y})_{k} \approx \lambda_{1}^{k} a_{1} \mathbf{X}_{M 1} / b_{1} b_{2} \ldots b_{k-1} .
\]

Из этого выражения видно, что $k$-е приближение вектора стремится к собственному вектору $\mathbf{X}_{M 1}$, поскольку остальные величины в этом выражении являются постоянными. Если в правой части выражения (о) выделить множитель $\lambda_{1}$, получим
\[
(\mathbf{Y})_{k} \approx \lambda_{1}\left(\lambda_{1}^{k-1} a_{1} \mathbf{X}_{M 1} / b_{1} b_{2} \ldots b_{k-1}\right)=\lambda_{1}(\mathbf{X})_{k} .
\]

Таким образом, $k$-е приближение для $\lambda_{i}$, задаваемое выражением (4.101), как это видно из выражения (o), принимает значения, близкие к $\lambda_{1}$. Благодаря указанному условию сходимости итерационного процесса численно наибольшее собственное значение $\lambda_{1}$ называют основной собственной частотой, а соответствующий этому значению вектор $\mathbf{X}_{\text {м1 }}$ – основным собственным вектором.

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

где
\[
\mathbf{F M} \mathbf{X}_{M 1}=\lambda_{1} \mathbf{X}_{M 1},
\]
\[
\lambda_{1}=\frac{1}{p_{1}^{2}} \text {. }
\]

Уравнение (4.103) представляет стандартную форму задачи на собственные значения [см. уравнение (4.99)]. В результате видим, что матрица
\[
\mathbf{A}=\mathbf{F M}
\]

является несимметричной, хотя здесь $\mathbf{M}$ – диагональная матрица с одинаковыми элементами на диагонали. Однако указанная потеря симметрии не существенна в том случае, когда основная форма колебаний определяется с помощью метода итераций.

В качестве числового примера рассмотрим применение уравнений движения в перемещения к трехмассовой системе (см. рис. 4.1,a), которая уже рассматривалась ранее с помощью уравнений движения в усилиях (см. пример 1 в п. 4.2). При $k_{1}=k_{2}=k_{3}=k$ и $m_{1}=m_{2}=m_{3}=m$ матрицы податливости масс для этой системы имеют вид
\[
\mathbf{F}=\delta\left[\begin{array}{lll}
1 & 1 & 1 \\
1 & 2 & 2 \\
1 & 2 & 3
\end{array}\right] ; \quad \boldsymbol{M}=m\left[\begin{array}{lll}
1 & 0 & 0 \\
0 & 1 & 0 \\
0 & 0 & 1
\end{array}\right],
\]

где, как и выше, имеем $\delta=1 / k$. В результате матрица А согласно выражению (4.105) принимает вид
\[
\mathbf{A}=m \delta\left[\begin{array}{lll}
1 & 1 & 1 \\
1 & 2 & 2 \\
1 & 2 & 3
\end{array}\right] .
\]

В качестве достаточно хорошего приближения для основной формы колебаний можно было бы взять суммы строк мартицы А. В результате получается вектор перемещений, обусловленных статически приложенными силами, которые, как и в методе Релея (см. п. 1.14), пропорциональны массам. Непрямой способ применения того же самого приема состоит в задании представления $(\mathbf{X})_{1}=\{1 ; 1 ; 1\}$ в качестве первого приближения для искомого вектора. Умножая вектор (X) $)_{1}$ слева на матрицу $\mathbf{A}$, согласно выражению (4.100) получим вектор $(\mathbf{Y})_{1}=m \delta\{3 ; 5 ; 6\}$. Первое приближение для $\lambda_{1}$, как следует из выражения (4.101), можно определить тремя различными путями. Для удобства проведения дальнейших вычислений разделим последнюю компоненту вектора (Y) $)_{1}$ на последнюю компоненту вектора (X) что дает $\left(\lambda_{1}\right)_{1}=\left(y_{n}\right)_{1} /\left(x_{n}\right)_{1}=6 m \delta$. Прежде чем перейти ко второму шагу итераций, пронормируем вектор (Y) $)_{1}$ путем деления каждой его компоненты на последнюю компоненту [см. выражение (4.102)], в результате получаем представление для второго приближения вектора ( $\mathbf{X})_{2}=\{0,500 ; 0,833 ; 1,000\}$. Когда используется нормирование подобного типа, делитель $b_{1}=6 \mathrm{~m} \delta$ приближенно равен собственному значению.

Второй шаг итерации состоит в использовании выражений (4.100)-(4.102) для получения соответственно $(\mathbf{Y})_{2},\left(\lambda_{1}\right)_{2}$ и $(\mathbf{X})_{3}$. Этот процесс итераций повторяется до тех пор, пока собственные векторы, определяемые на двух последующих шагах итераций, не совпадут с заранее заданной точностью. В табл. 4.1 приведены результаты расчетов, из которых видно, что в этом случае вектор, получаемый на пятом шаге итераций, совпадает с вектором, найденным на четвертом шаге, с точностью до третьего знака после запятой. Таким образом, окончательно можно принять
\[
\begin{array}{c}
\lambda_{1} \approx 5,049 m \delta, \quad p_{1}^{\mathrm{a}}=1 / \lambda_{1} \approx 0,198 \frac{k}{m}, \\
\mathbf{X}_{M 1} \approx\{0,445 ; 0,802 ; 1,000\},
\end{array}
\]

что полностью совпадает с результатами, полученными в примере 1 (см. п. 4.2), за исключением отличающихся способом нормирования векторов $\mathbf{X}_{\text {M1 }}$.

После того как найдена основная форма колебаний, можно продолжить этот процесс дальше и с помощью метода последователь-

ных приближений определить собственные значения и собственные векторы для высших форм колебаний. Если благодаря введению соответствующих подкреплений первая форма колебаний не реализуется, основной становится вторая форма колебаний. Если не реализуются ни первая, ни вторая формы, основной становится третья форма колебаний и т. д. Поскольку число собственных форм колебаний равно числу степеней свободы, введение дополнительных связей для соответствующих форм колебаний уменьшает число степеней свободы системы. Таким образом, можно ожидать, что матрица коэффициентов после второго шага итераций будет иметь порядок $n-1$, после третьего шага итерации $n-2$ и т. д. Однако при проведении числовых расчетов, а также при использовании ЭВМ удобнее сохранять порядок матрицы, равный $n$. Для этого применяется простой прием, суть которого объяснена ниже.

Ограничение на реализацию соответствующих форм колебаний в системе можно ввести, положив равным нулю соответствующие перемещения в главной форме. Из выражений (4.32) и (4.44а) следует соотношение
\[
\mathbf{X}_{\Gamma}=\mathbf{X}_{M}^{-1} \mathbf{X}=\mathbf{M}_{\Gamma}^{-1} \mathbf{X}_{M}^{\mathrm{T}} \mathbf{M} \mathbf{X},
\]

которое связывает главные координаты в векторе $\mathbf{X}_{\Gamma}$ с исходными координатами вектора. $\mathbf{X}$. Для того чтобы исключить появление первой формы колебаний, положим равным нулю компоненту $x_{\text {Г1 }}$ в выражении (т) для вектора $\mathbf{X}_{\Gamma}$ :
\[
x_{\Gamma 1}=\frac{1}{M_{\Gamma 1}} \mathbf{X}_{M 1}^{\mathrm{T}} \mathbf{M X}=0 .
\]

Если вектор $\mathbf{X}$ берется в качестве собственного вектора $\mathbf{X}_{\text {мі }}$ (где $i=2,3, \ldots, n$ ), то видно, что это условие положения дополнительной связи совпадает с условием ортогональности первой и высших форм колебаний по отношению к матрице $\boldsymbol{M}$. Для упрощения полагаем матрицу масс диагональной, тогда развернутая форма выражения (y) имеет вид
\[
M_{11} X_{M 11} x_{1}+M_{22} X_{M 21} x_{2}+M_{33} X_{M 31} x_{3}+\cdots+M_{n n} X_{M n 1} x_{n}=0 .
\]

Поскольку это равенство должно выполняться всегда, выразим перемещения $x_{1}$ через остальные перемещения следующим образом:
\[
x_{1}=-\frac{M_{22} X_{M 21}}{M_{11} X_{M 11}} x_{2}-\frac{M_{33} X_{M 31}}{M_{11} X_{M 11}} x_{3}-\cdots-\frac{M_{n n} X_{M n 1}}{M_{11} X_{M 11}} x_{n} .
\]

Подстановкой этого выражения для $x_{1}$ в исходные уравнения задачи на собственные значения [см. уравнение (4.99)] получаем $n$ уравнений для $n-1$ неизвестных. Первое из этих уравнений будет представлять собой линейную комбинацию остальных $n-1$ уравнений и поэтому может быть отброшено, в результате чего получается редуцированная система $n-1$ уравнений с $n-1$ неизвестными. Однако при определении второй формы колебаний можно проводить итерации для полной системы $n$ линейно зависимых уравнений, не понижая порядка матрицы, поскольку отсутствует первая форма колебаний уравнений. С этой целью представим выражение (ф), дополненное тривиальными соотношениями вида $x_{2}=x_{2}$, $x_{3}=x_{3}, \ldots$, в матричной форме
\[
\left[\begin{array}{c}
x_{1} \\
x_{2} \\
x_{3} \\
\ldots \\
x_{n_{-}}
\end{array}\right]^{-}\left[\begin{array}{ccccc}
0 & c_{12} & c_{13} & \ldots & c_{1 n}^{-} \\
0 & 1 & 0 & \ldots & 0 \\
0 & 0 & 1 & \ldots & 0 \\
\ldots & \ldots & \ldots & \ldots & \ldots \\
0 & 0 & 0 & \ldots & 1
\end{array}\right]\left[\begin{array}{c}
-x_{1}^{\prime} \\
x_{2} \\
x_{3} \\
\ldots \\
x_{n_{-}}
\end{array}\right],
\]

где
\[
c_{12}=-\frac{M_{22} X_{M 21}}{M_{11} X_{M 11}}, \quad c_{13}=\frac{M_{33} X_{M 31}}{M_{11} X_{M 11}}, \ldots, c_{1 n}=\frac{M_{n n} X_{M n 1}}{M_{11} X_{M 11}} .
\]

Первый элемент $x_{1}^{\prime}$ стоящего в правой части уравнения (х) векторастолбца является фиктивным перемещением, которое всегда умножается на ноль. В краткой форме это матричное уравнение имеет вид
\[
\mathbf{X}=\mathbf{T}_{\mathrm{B} 1} \mathbf{X}^{\prime},
\]

где матрица $\mathrm{T}_{\mathrm{B} 1}$ характеризует зависимость $x_{1}, x_{2}, x_{3}, \ldots, x_{n}$ от $x_{2}, x_{3}, \ldots, x_{n}$. Подставляя в соответствии с выражением (ч) представление для $\mathbf{X}_{M i}$ в правую часть уравнения (4.99), получаем
\[
\mathbf{A T}_{\mathrm{B} 1} \mathbf{X}_{M i}=\lambda_{i} \mathbf{X}_{M i} .
\]

Присутствие в этом соотношении матрицы $\mathbf{T}_{\mathrm{B} 1}$ обеспечивает наличие линейной зависимости, что необходимо для проведения итераций при нахождении второй формы колебаний. На каждом шаге итерации эта матрица умножается слева на выбранное представление для вектора перемещений и предполагается их ортогональность с первой формой колебаний. Однако для удобства эта матрица используется только один раз – на эту матрицу справа умножается матрица A. Матрица $\mathbf{T}_{\mathrm{B} 1}$ называется вычщщаюей, поскольку ее влияние сказывается в том, что она «вычищает» все относящееся

к первой форме колебаний и тем самым позволяет второй форме стать главной формой колебаний. Таким образом, можно записать

где
\[
\mathbf{A}_{1} \mathbf{X}_{M i}=\lambda_{i} \mathbf{X}_{M i},
\]
\[
\mathbf{A}_{1}=\mathbf{A T}_{\mathrm{B} 1},
\]

а матрица $\mathbf{T}_{\mathrm{B} 1}$ имеет вид $(\mathrm{x})$.
Для того чтобы показать, как применяется этот способ, определим собственное значение и собственный вектор, соответствующие второй форме колебаний системы, показанной на рис. 4.1, a. Зная первый собственный вектор $\mathbf{X}_{M 1}$, используем его компоненты для построения матрицы $\mathrm{T}_{\mathrm{B}_{1}}$ в следующей форме:
\[
\mathbf{T}_{\mathrm{B} 1}=\left[\begin{array}{ccc}
0 & -1,802 & -2,247 \\
0 & 1 & 0 \\
0 & 0 & 1
\end{array}\right] .
\]

Умножая матрицу А [см. выражение (p)] справа на матрицу $\mathbf{T}_{\mathrm{B} 1}$, получим
\[
\mathbf{A}_{\mathbf{1}}=m \delta\left[\begin{array}{rrr}
0 & -0,802 & -1,247 \\
0 & 0,198 & 0,247 \\
0 & 0,198 & 0,753
\end{array}\right] .
\]

Предположив, что нам ничего неизвестно относительно второй формы колебаний, возьмем первое приближение для вектора перемещений в виде $(\mathbf{X})_{1}=\{1 ; 1 ; 1\}$, что является плохим приближением для действительного значения собственного вектора. Более разумный выбор этого приближения позволил бы сойтись процессу итераций кочном тешению за меньшее число шагов.

Умножая вектор-столбец ( $\mathbf{X})_{i}$ на матрицу $\mathbf{A}_{1}$ (согласно выражению (4.100)), получим вектор (Y) $)_{1}=m \delta\{-2,049 ;-0,049 ; 0,951\}$. Разделив согласно выражению (4.101) последний элемент вектора (Y) $)_{1}$ на последний элемент вектора ( $\left.\mathbf{X}\right)_{1}$, найдем для $\lambda_{2}$ первое приближение $\left(\lambda_{2}\right)_{1}=0,951 \mathrm{~m} \delta$. Затем нормализуем вектор (Y) $)_{1}$ путем деления каждой его компоненты на значение последней компоненты [см. выражение (4.102) ], что в результате дает второе приближение для вектора перемещений $(\mathbf{X})_{2}=\{-2,155 ;-0,052 ; 1,000\}$. Эти результаты, а также результаты, получаемые на последующих шагах итераций, представлены в табл. 4.2, причем после пяти шагов итераций процесс еще не сходится. Продолжая этот процесс в течение одиннадцати шагов итераций, получим, что процесс итераций сходится при следующих значениях параметров:
\[
\begin{array}{c}
\lambda_{2} \approx 0,643 m \delta ; \quad p_{2}^{2}=\frac{1}{\lambda_{2}} \approx 1,555 \frac{k}{m} ; \\
\mathbf{X}_{M_{2}} \approx\{-1,247 ;-0,555 ; 1,000\} .
\end{array}
\]

За исключением способа, которым был пронормирован вектор $\mathbf{X}_{M_{2}}$, эти результаты совпадают с теми, которые были получены в примере 1 (см. п. 4.2).

После определения второй формы колебаний ее можно исключить из системы уравнений с помощью процедуры, аналогичной той, которая использовалась для исключения первой формы колебаний. Используя условия введения дополнительных связей $x_{\mathrm{r} 1}=0$ и $x_{\text {г } 2}=0$, выразим $x_{2}$ через компоненты $x_{3}, \ldots, x_{n}$ и представим соотношения между перемещениями в следующей матричной форме:
\[
\left[\begin{array}{c}
x_{1}^{\prime} \\
\ldots \ldots . \\
x_{2} \\
x_{3} \\
\ldots \\
x_{n}
\end{array}\right]=\left[\begin{array}{cc:ccc}
1 & 0 & 0 & \ldots & 0 \\
\ldots & 0 & d_{23} & \ldots & d_{n 2} \\
0 & 0 & 1 & \ldots & 0 \\
\ldots & \ldots & \ldots & \ldots & \ldots \\
0 & 0 & 0 & \ldots & 1
\end{array}\right]\left[\begin{array}{c}
x_{1}^{\prime} \\
x_{2}^{\prime} \\
\ldots \\
x_{3} \\
\ldots \\
x_{n}
\end{array}\right] \text {, }
\]

где
\[
\begin{array}{c}
d_{23}=-\frac{M_{33}\left(X_{M 11} X_{M 32}-X_{M 12} X_{M 31}\right)}{M_{22}\left(X_{M 11} X_{M 22}-X_{M 12} X_{M 21}\right)}, \ldots, \\
d_{n 3}=-\frac{M_{n n}\left(X_{M 11} X_{M n 2}-X_{M 12} X_{M n 1}\right)}{M_{22}\left(X_{M 11} X_{M 22}-X_{M 12} X_{M 21}\right)} .
\end{array}
\]

Стоящие в скобках выражения являются минорами первых двух столбцов матрицы $\mathbf{X}_{M}$, которые появляются в процессе получения решений для $x_{2}$, выраженных через $x_{3}, \ldots, x_{n}$. Дополнительное фиктивное перемещение $x_{2}^{\prime}$ фигурирует в качестве компонент вектора, стоящего в правой части уравнения ( $\left.\sigma^{\prime}\right)$. Это уравнение можно переписать в более компактной форме
\[
\mathbf{X}^{\prime}=\mathbf{T}_{\mathrm{B} 2} \mathbf{X}^{\prime \prime} .
\]

Матрица $\mathrm{T}_{\mathrm{B} 2}$ характеризует зависимость $x_{2}, x_{3}, \ldots, x_{n}$ от $x_{3}, \ldots, x_{n}$. Подстановкой выражения ( $\left.r^{\prime}\right)$ для $\mathbf{X}_{\text {мi }}$ в левую часть уравнения (4.106) получим
\[
\mathbf{A}_{2} \mathbf{X}_{M i}^{\prime \prime}=\lambda_{i} \mathbf{X}_{M i},
\]

где
\[
\mathbf{A}_{2}=\mathbf{A}_{1} \mathbf{T}_{\mathrm{B} 2} .
\]

Теперь уравнение (4.108), в котором отсутствуют члены, описывающие влияние первых двух форм колебаний, можно использовать в итерационной процедуре для определения третьей формы колебаний.

Применяятосанный подход к системе, изображенной на рис. 4.1, $a$, найдем что первое из выражений (в’)
\[
d_{23}=\frac{0,445 \cdot 1,000-(-1,247) 1000}{0,445(-0,555)-(-1,247) 0,802}=-\frac{1,692}{0,753}=-2,247,
\]

а матрица
\[
\mathbf{T}_{\mathrm{B} 2}=\left[\begin{array}{cc:c}
1 & 0 & 0 \\
\hdashline 0 & 0 & -2,247 \\
0 & 0 & 1
\end{array}\right] .
\]

Умножая матрицу А [см. выражение (э)] справа на матрицу, получим
\[
\mathbf{A}_{2}=m \delta\left[\begin{array}{ll:r}
0 & 0 & 0,555 \\
0 & 0 & -0,692 \\
0 & 0 & 0,308
\end{array}\right] .
\]

Третий столбец матрицы $\mathbf{A}$ пропорционален третьему собственному вектору, поэтому для определения последней искомой формы колебаний не требуется применять итерации. Более того, третья компонента этого вектора равна собственному значению. Таким образом, имеем
\[
\begin{array}{c}
\lambda_{3} \approx 0,308 m \delta, \quad p_{2}^{2}=1 / \lambda_{2} \approx 3,247 k / m, \\
\mathbf{X}_{M 3} \approx\{1,802,-2,247,1,000\} .
\end{array}
\]

Эти значения совпадают с результатами, полученными в примере 1 (см. п. 4.2), за исключением вектора $\mathbf{X}_{M 3}$, который здесь нормирован иначе.

Если система имеет два или более одинаковых собственных значения, они будут в равной степени главными, и собственный вектор, к которому будет сходиться итерационный процесс, в свою очередь зависит от выбранного первого приближенного значения для вектора перемещений. С помощью матрицы $\mathbf{T}_{\mathrm{B}}$ каждый последующий собственный вектор становится ортогональным к предыдущему собственному вектору, причем это имеет место и в том случае, когда имеются кратные собственные значения. Поскольку собственные векторы, соответствующие кратным собственным значениям, часто имеют нулевые компоненты, необходимо внимательно ознакомиться с проведением процесса ортогонализации, чтобы избежать деления на ноль.

В качестве простого примера системы с кратными корнями возьмем систему, рассмотренную в задаче 4.2 .9 (см. п. 4.2), и предположим, что три пружины, на которых закреплена масса, расположены вдоль осей $x, y$ и $z$. Положив $k_{1}=k_{2}=$ $=k_{3}=k$, найдем матрицу А для этой задачи:
\[
\mathbf{A}=\mathbf{F} \mathbf{M}=m \delta\left[\begin{array}{lll}
1 & 0 & 0 \\
0 & 1 & 0 \\
0 & 0 & 1
\end{array}\right] .
\]

Как можно убедиться простой проверкой, собственные значения для этой системы $\lambda_{1}=\lambda_{2}=\lambda_{3}=m \delta$. Однако собственные векторы здесь неизвестны и их надо определить, используя тот же подход, что и при решении предыдущего примера. Выбор приближенного выражения для вектора в виде $\left\{x_{1}, y_{1}, z_{1}\right\}=\{1,1,1\}$ позволяет удовлетворить уравнению (4.103), поэтому этот вектор становится первым собственным вектором системы. Первая матрица $\mathbf{T}_{\mathrm{B} 1}$ [см. выражения (х), (ц) и (ч)] имеет вид
\[
\mathbf{T}_{\mathrm{B} 1}=\left[\begin{array}{r:rr}
0 & -1 & -1 \\
0 & 1 & 0 \\
0 & 0 & 1
\end{array}\right],
\]

а матрицу А находим с помощью выражения (4.107):
\[
\mathbf{A}_{1}=\mathbf{A T}_{\mathrm{B} 1}=m \delta\left[\begin{array}{r:rr}
0 & -1 & -1 \\
0 & 1 & 0 \\
0 & 0 & 1
\end{array}\right] .
\]

Для определения второй формы колебаний опять возьмем первое приближение для вектора перемещений в форме $\left\{x_{1} ; y_{1} ; z_{1}\right\}=\{1 ; 1 ; 1\}$, тогда на втором шаге

итераций находим второй собственный вектор $\left\{x_{1}, y_{1}, z_{1}\right\}=\{-2,1,1\}$. Вторая матрица $\mathbf{T}_{\mathrm{B} 2}$ [см. выражения ( $\left.\sigma^{\prime}\right),\left(\mathrm{B}^{\prime}\right)$ и (г’)] принимает вид
\[
\mathbf{T}_{\mathrm{B} 2}=\left[\begin{array}{rr:r}
1 & 0 & 0 \\
\hdashline & 0 & -1 \\
0 & 0 & -1 \\
0 & 0 & 1
\end{array}\right] \text {. }
\]

C помощью этой матрицы и выражения (4.109) находим матрицу $\mathbf{A}_{2}$ :
\[
\mathbf{A}_{2}=\mathbf{A}_{1} \mathbf{T}_{\mathrm{B} 2}=m \delta\left[\begin{array}{ll:r}
0 & 0 & 0 \\
0 & 0 & -1 \\
0 & 0 & 1
\end{array}\right] .
\]

Таким образом, третий собственный вектор равен $\left\{x_{1}, y_{1}, z_{1}\right\}=\{0,-1,1\}$ и ортогонален двум остальным собственным векторам относительно матрицы $\mathbf{M}$. Подобная система собственных векторов не является единственной, но она удовлетворяет условиям ортогональности, выполнение которых необходимо для собственных векторов при использовании метода нормальных форм колебаний.

Итерационный процесс понижения числа степеней свободы системы, описанный выше, теоретически можно применять многократно до тех пор, пока не будут определены все частоты и формы колебаний системы со многими степенями свободы. Однако каждое собственное значение и собственный вектор, определяемые таким образом, являются только приближенными. Поэтому проводимая на каждом шаге ортогонализация будет неполной. Более того, каждое понижение числа степеней свободы сопровождается ошибками округления, которые накапливаются с каждым шагом. С вопросом о точности связано и то обстоятельство, что для получения большого числа частот и форм колебаний требуется выполнять необычно большое число арифметических операций. Следовательно, как об этом уже говорилось в начале данного параграфа, итерационный метод лучше всего использовать в том случае, когда требуется определить только несколько низших форм колебаний. Кроме того, необходимость выполнения большого числа арифметических операций в случае систем с очень большим числом степеней свободы требует применения ЭВМ, особенно тогда, когда трудно предугадать формы колебаний. Поэтому в приложении к книге дан текст программы на языке БЕЙСИК, под названием EIGIT3, которая позволяет вычислять три первые собственные значения и собственные векторы матрицы с помощью итерационного метода.

Хотя итерационный метод позволяет находить только несколько собственных значений и собственных векторов системы, это не мешает использовать метод нормальных форм колебаний при определении динамических перемещений в системе. Если найдены $n_{1}$ форм колебаний, где $n_{1} \leqslant n$, матрица форм (или $\mathbf{X}_{M}$ ) содержит вместо $n$ только $n_{1}$ столбцов. Такая прямоугольная матрица не имеет обратной матрицы, поэтому вместо выражений, содержащих обратные матрицы, следует использовать выражения (4.44а) и (4.44б), в которых имеются транспонированные матрицы $\mathbf{X}_{M}^{\mathrm{T}}$ и $\mathbf{X}_{\mathrm{H}}^{\mathrm{T}}$; при этом удается определить только $n_{1}$ первых нормальных форм колебаний, тогда как влиянием остальных форм колебаний на суммарное динамическое

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

Здесь не было рассмотрено применение итерационного метода к задачам на собственные значения, представленным в виде уравнений движения в усилиях [см. уравнение (4.17) ], поскольку главным при этом были бы наибольшие собственные значения $p_{n}^{2}$. В задаче, в которой проще определяются коэффициенты жесткости, а не податливости, можно всегда обратить неособенную матрицу жесткостей $\mathbf{S}$ и тем самым получить матрицу податливостей $\mathbf{F}$, которая имеется в уравнении (4.103). С другой стороны, для полуопределенной системы, матрица жесткостей которой является особенной, требуется проводить специальное исследование. В этом случае матрицы жесткостей и податливостей следует редуцировать путем перехода к новой системе координат, чтобы исключить формы движения как абсолютно жесткого тела, которые можно определить с помощью простого рассмотрения и составить процедуру для исключения этих форм.

Допустим известно, что колеблющаяся система имеет только одну форму движения как абсолютно жесткого тела, которая обозначена как первый собственный вектор $\mathbf{X}_{\text {м1 }}$. Эту форму исключаем из системы уравнений движения в усилиях, задав специальное условие дополнительного закрепления в виде $x_{\Gamma 1}=0$. Затем строим $y с е$ ченную матрицу преобразованил $\mathbf{T}_{\text {оп1 }}$, состоящую из подматрицы $n \times(n-1)$ матрицы $\mathbf{T}_{\mathrm{B} 1}$ справа от разделяющей линии в выражении (x). В данном случае воспользуемся соотношением
\[
\mathbf{X}=\mathbf{T}_{\text {оп1 }} \mathbf{X}_{\text {оп1 }},
\]

выражающим зависимость вектора $\mathbf{X}=\left\{x_{1}, x_{2}, x_{3}, \ldots, x_{n}\right\}$ от вектора $\mathbf{X}_{\text {оп1 }}=\left\{x_{2}, x_{3}, \ldots, x_{n}\right\}$. Последний вектор не содержит фиктивного перемещения $x_{1}^{\prime}$, которое использовалось в матрице $\mathbf{T}_{\mathrm{B} 1}$. Поскольку компоненты матрицы $\mathbf{T}_{\text {оп1 }}$ не зависят от времени, можно также записать
\[
\ddot{\mathbf{X}}=\mathbf{T}_{\text {on1 }} \ddot{\mathbf{X}}_{\text {on1 }},
\]

где $\ddot{\mathbf{x}}_{\text {оп1 }}=\left\{\ddot{x}_{2}, \ddot{x}_{3}, \ldots, \ddot{x}_{n}\right\}$. Подставляя выражения (4.110) и ( $\left(^{\prime}\right)$ в уравнения движения в усилиях, записанные для случая свободных колебаний [см. уравнение (4.30)], и умножая результат слева на матрицу $\mathbf{T}_{\text {oп1 }}^{\mathrm{T}}$, получим
\[
\mathbf{T}_{\text {on1 }}^{\mathrm{T}} \mathbf{M} \mathbf{T}_{\text {on1 }} \ddot{\mathbf{X}}_{\text {on1 }}+\mathbf{T}_{\text {on1 }}^{\mathrm{T}} \mathbf{S T}_{\text {on1 }} X_{\text {on1 }}=\mathbf{0}
\]

или
\[
\mathbf{M}_{\mathrm{ou1}} \ddot{\mathbf{X}}_{\mathrm{om} 1}+\mathrm{S}_{\mathrm{or} 1} X_{\mathrm{om1}}=\mathbf{0},
\]

где
\[
\mathbf{M}_{\mathrm{on1}}=\mathbf{T}_{\mathrm{on1}}^{\mathrm{T}} \mathbf{M T}_{\mathrm{on} 1}, \quad \mathrm{~S}_{\mathrm{on1}}=\mathbf{T}_{\mathrm{on1}}^{\mathrm{T}} \mathbf{S T}_{\mathrm{on1}} .
\]

Здесь $\boldsymbol{M}_{\text {оп1 }}$ и $\boldsymbol{S}_{\text {оп1 }}$ – квадратные симметричные матрицы порядка $n-1$, а уравнения (4.111) представляют собой уравнения движения, преобразованные к редуцированной системе координат, в кото-

рых исключены перемещения как абсолютно жесткого тела. В этой системе координат матрицу жесткостей $S_{0 п 1}$ можно обратить, в результате получим матрицу податливостей
\[
\mathbf{F}_{\text {oni }}=\mathbf{S}_{\text {onl }}^{-1} .
\]

Для подготовки применения ее в дальнейшем к решению задачи на собственные значения, соответствующей уравнениям (4.111), методом последовательных приближений запишем произведение
\[
\mathbf{A}_{\text {oп1 }}=\mathbf{F}_{\text {on1 }} \mathbf{M}_{\text {oп1 }} .
\]

Собственные векторы, определяемые с помощью подобной процедуры, преобразуются к исходным координатам согласно выражению (4.110).

Если в системе имеет место вторая форма движения как жесткого тела, определяемая собственным вектором $\mathbf{X}_{\text {м2 }}$, требуется ввести второе условие, ограничивающее перемещения, в виде $x_{\Gamma 2}=0$. Усеченная матрица преобразования $\mathbf{T}_{\text {оп } 2}$, соответствующая этому условию, представляет подматрицу порядка $(n-1) \times(n-2)$ матрицы $\mathbf{T}_{\mathrm{B} 2}$, расположенную ниже и правее разделительных штриховых линий в выражении ( $\left.\sigma^{\prime}\right)$. В этом случае имеем соотношение
\[
\mathbf{X}_{\text {oп2 } 2}=\mathbf{T}_{\text {оп2 }} \mathbf{X}_{\text {oч2 }},
\]

выражающее связь между векторами $\mathbf{X}_{\text {оп1 }}=\left\{x_{2}, x_{3}, \ldots, x_{n}\right\}$ и $\mathbf{X}_{\text {оп2 }}=\left\{x_{3}, \ldots, x_{n}\right\}$. Повторяя приведенные выше рассуждения, можем свести матрицы $\boldsymbol{M}_{\text {оп1 }}$ и $\boldsymbol{S}_{\text {оп1 }}$ к матрицам $\boldsymbol{M}_{\text {оп2 }}$ и $\boldsymbol{S}_{\text {оп } 2}$ с помощью следующих преобразований:
\[
\mathbf{M}_{\mathrm{on} 2}=\mathbf{T}_{\mathrm{on} 2}^{\mathrm{T}} \mathbf{M}_{\mathrm{on} 1} \mathbf{T}_{\mathrm{on} 2}, \quad \mathbf{S}_{\mathrm{on} 2}=\mathbf{T}_{\text {оп2 }}^{\mathrm{T}} \mathrm{S}_{\mathrm{on} 1} \mathbf{T}_{\text {оп2 }} .
\]

С другой стороны, подставив выражение (4.115) в соотношение (4.110), получим
\[
\mathbf{X}:=\mathbf{T}_{\text {оп2 }}^{*} \mathbf{X}_{\text {oп2 }},
\]

где
\[
\mathbf{T}_{\text {оп2 }}^{*}=\mathbf{T}_{\text {оп1 }} \mathbf{T}_{\text {оп2 }} .
\]

Комбинированная матрица $\mathbf{T}_{\text {оп2 }}^{*}$ порядка $n \times(n-2)$ выражает зависимость вектора $\mathbf{X}$ от вектора $\mathbf{X}_{\text {оп2 }}$ и равна произведению маттриц $\mathbf{T}_{0 n 1}$ и $\mathbf{T}_{\text {сп2 }}$ порядков соответственно $n \times(n-1)$ и $(n-1) \times$ $\times(n-2)$. С помощью этой комбинированной матрицы можно преобразовать матрицы $\boldsymbol{M}$ и $\boldsymbol{S}$ непосредственно в матрицы $\boldsymbol{M}_{\text {on2 }}$ и $\boldsymbol{S}_{\text {оп2 }}$ с помощью следующих преобразований:
\[
\mathbf{M}_{\mathrm{on} 2}=\left(\mathbf{T}_{\mathrm{on} 2}^{*}\right)^{\mathrm{T}} \mathbf{M T}_{\mathrm{on} 2}^{*}, \quad \mathrm{~S}_{\mathrm{on} 2}=\left(\mathbf{T}_{\mathrm{on} 2}^{*}\right)^{\mathrm{T}} \mathbf{S T}_{\mathrm{on} 2}^{*} .
\]

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

Для иллюстрации этого метода предположим, что первая пружина в системе, изображенной на рис. 4.1, a, имеет жесткость $k_{1}=0$. Если $k_{2}=k_{3}=k$ и $m_{1}=$ $=m_{2}=m_{3}=m$, то матрицы жесткостей и масс имеют вид
\[
\mathbf{S}=k\left[\begin{array}{rrr}
1 & -1 & 0 \\
-1 & 2 & -1 \\
0 & -1 & 1
\end{array}\right] ; \quad \mathbf{M}=m\left[\begin{array}{lll}
1 & 0 & 0 \\
0 & 1 & 0 \\
0 & 0 & 1
\end{array}\right] .
\]

Проверкой убеждаемся, что система является полуопределенной и что форма, coответствующая движению как абсолютно жесткого тела, определястся следующим вектором:
\[
\mathrm{X}_{M 1}=\{1 ; 1 ; 1\} .
\]

С учетом сказанного редуцированная матрица преобразования $\mathbf{T}_{\text {оп }}$ порядка $3 \times 2$ имеет вид
\[
\mathbf{T}_{\text {on1 }}=\left[\begin{array}{rr}
-1 & -1 \\
1 & 0 \\
0 & 1
\end{array}\right] .
\]

Используя эту матрицу, преобразуем матрицы масс и жесткостей с помощью выражения (4.112) к виду
\[
\begin{array}{c}
\boldsymbol{M}_{0 \Pi 1}=m\left[\begin{array}{lll}
-1 & 1 & 0 \\
-1 & 0 & 1
\end{array}\right]\left[\begin{array}{lll}
1 & 0 & 0 \\
0 & 1 & 0 \\
0 & 0 & 1
\end{array}\right]\left[\begin{array}{rr}
-1 & -1 \\
1 & 0 \\
0 & 1
\end{array}\right]=m\left[\begin{array}{ll}
2 & 1 \\
1 & 2
\end{array}\right] \\
\mathbf{S}_{0 \Pi 1}=k\left[\begin{array}{lll}
-1 & 1 & 0 \\
-1 & 0 & 1
\end{array}\right]\left[\begin{array}{rrr}
1 & -1 & 0 \\
-1 & 2 & -1 \\
0 & -1 & 1
\end{array}\right]\left[\begin{array}{rr}
-1 & -1 \\
1 & 0 \\
0 & 1
\end{array}\right]=k\left[\begin{array}{ll}
5 & 1 \\
1 & 2
\end{array}\right] .
\end{array}
\]

Обращая последнюю матрицу, получим матрицу податливостей
\[
\mathbf{F}_{\text {оп1 }}=\mathbf{S}_{\text {oп1 }}^{-1}=\frac{\delta}{9}\left[\begin{array}{rr}
2 & -1 \\
-1 & 5
\end{array}\right] .
\]

Далее, в соответствии с выражением (4.114) найдем произведение
\[
\mathbf{A}_{\mathrm{o \Pi 1}}=\mathbf{F}_{\mathrm{O \Pi 1}} \mathbf{M}_{\mathrm{oH1}}=\frac{m \delta}{3}\left[\begin{array}{ll}
1 & 0 \\
1 & 3
\end{array}\right] .
\]

Легко проверить либо прямым решением, либо методом последовательных приближений, что собственные значения этой матрицы суть $\lambda_{2}=m \delta$ и $\lambda_{3}=m \delta / 3$, а собственные векторы ( $\left.\mathbf{X}_{0 \text { п1 }}\right)_{M 2}=\{0,-1\},\left(\mathbf{X}_{\mathrm{OH1}}\right)_{N 3}=\{-2,1\}$. Преобразуя эти векторы к исходным координатам с помощью выражения (4.110), найдем
\[
\mathbf{X}_{M_{2}}=\mathbf{T}_{0 \Pi 1}\left[\begin{array}{r}
0 \\
-1
\end{array}\right]=\left[\begin{array}{r}
1 \\
0 \\
-1
\end{array}\right], \quad \mathbf{X}_{M 3}=\mathbf{T}_{0 \Pi 1}\left[\begin{array}{r}
-2 \\
1
\end{array}\right]=\left[\begin{array}{r}
1 \\
-2 \\
1
\end{array}\right],
\]

что совпадает с решением, полученным в последней части примера 1 в п. 4.2.

ЗАДАЧИ

4.7.1. Используя метод последовательных приближений, определить собственные значения и собственные векторы для системы, изображенной на рис. 4.2, $a$. Принять $m_{1}=m_{3}=m, m_{2}=2 m, l_{1}=l_{2}=l_{3}=l_{4}=l$.
Omвem: $\lambda_{1,2,3} \approx 1,707 \mathrm{ml} / \mathrm{T} ; 0,500 \mathrm{ml} / \mathrm{T} ; 0,293 \mathrm{ml} / \mathrm{T}$.
4.7.2. Решить задачу 4.2 .2 методом последовательных приближений, приняв $m_{1}=m_{3}=2 m, m_{2}=m$.
Omвem: $\lambda_{1,2,3} \approx 2,618 \mathrm{~m} \delta ; 1,000 \mathrm{~m} \delta ; 0,382 \mathrm{m \delta}$.
4.7.3. Решить задачу 4.2 .6 методом последовательных приближений, приняв $m_{1}=m_{2} \approx m, m_{3}=2 m$.
Omвem: $\lambda_{1,2,3} \approx 39,68 \alpha ; 2,815 \alpha ; 0,501 \alpha\left[\alpha=m l^{3} /(768 E I)\right]$.
4.7.4. Решить задачу 4.2 .7 методом последовательных приближений, приняв $m_{1}=m_{2}=3 m, m_{3}=m$.
Omвem: $p_{1}^{2}, 2,3 \approx 0,246 g / l ; 1,252 g / l ; 2,169 g / l$.
4.7.5. Решить задачу 4.2 .8 методом последовательных приближений, приняв $m_{1}=m_{2}=2 m, m_{3}=m$.
Ответ: $\lambda_{1,2,3} \approx 19,12 \alpha ; 4,000 \alpha ; 1,884 \alpha\left[\alpha=m h^{3} /(144 E I)\right]$.
4.7.6. Решить задачу 4.2.10 методом последовательных приближений, приняв $m_{1}=m$ и $m_{2}=3 m$.
Oтвет: $\lambda_{1,2,3} \approx 76,32 \alpha ; 8,978 \alpha ; 0,700 \alpha\left[\alpha=m l^{3} /(48 E I)\right]$.
4.7.7. Решить задачу 4.2.11, используя метод исключения форм движения как абсолютно жесткого тела, описанный в конце этого параграфа. Принять $m_{1}=m_{3}=$ $=m, m_{2}=2 m$.
Omвem: $p_{1,2,3}^{2}=0 ; 0 ; 6 E I /\left(m l^{3}\right)$.