Решение систем линейных алгебраических уравнений с ленточными матрицами

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

Сильная разреженность матриц особенно характерна для математических моделей микроуровня при дискретизации краевой дифференциальной задачи на основе методов сеток. Замена частных производных системы дифференциальных уравнений отношениями конечных разностей на основе использования регулярной сетки приводит к ленточным матрицам, у которых ненулевые элементы располагаются только на k главных диагоналях (рис. 8.1). Такой вид имеют, например, матрицы (2.78) и (2.87), полученные в примерах 2.1 и 2.2 при моделировании теплопередачи.

Ленточнодиагональная матрица

Рис. 8.1. Ленточнодиагональная матрица

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

При этомaj о = 0Hfli>n+1 =0.

Метод прогонки включает два этапа: прямую прогонку и обратную прогонку.

Прямая прогонка заключается в вычислении прогоночных коэффициентов на основе рекуррентных выражений:

Значения коэффициентов и р принимают равными нулю.

При обратной прогонке последовательно вычисляют значения элементов вектора искомых переменных X, начиная с хп - рп+, по формуле

Если число диагоналей больше трех, применяют метод матричной прогонки [19].

Анализ статических состояний линейных технических систем

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

Анализ статического состояния проектируемой механической системы обычно сводится к определению геометрических координат и усилий или моментов в упругих элементах в состоянии покоя.

При построении математической модели механической системы координатный базис дифференциальных уравнений (8.2) могут составлять либо переменные типа потокам = (i^),i = 1 ,п, и типа потенциала^ = (Fyj), у = 1,2V, где и — вектор скоростей сосредоточенных масс; Fy — вектор потенциалов (усилий или моментов) упругих элементов; п — число степеней свободы; N — число упругих элементов; либо только переменные типа потока v=(i>,) и x = (xt), i-l,n, где xt — геометрические координаты сосредоточенных масс. В последнем случае для определения усилий в упругих элементах используют дополнительно компонентные уравнения, представляемые функциями вида Fy(x).

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

Математическая модель анализа статических состояний — система алгебраических уравнений. До решения уравнений необходимо оценить устойчивость и однозначность статического состояния объекта. Для этого используется спектр матрицы Якоби. Если вещественные части всех собственных значений отрицательны, объект обладает физической устойчивостью. Но при этом могут быть нулевые собственные значения. Это означает, что состояние покоя механической системы безразлично и достигается при любом произвольном положении объекта в фазовом пространстве. Система алгебраических уравнений в таком случае имеет неоднозначное решение. Такие системы называют полу определенными. Однозначное положение объекта в пространстве фазовых координат можно получить путем исключения избыточных степеней свободы, число которых равно количеству нулевых собственных значений. С этой целью накладывают ограничения на геометрические координаты для такого же количества сосредоточенных масс, зафиксировав их на постоянном уровне. Фиксировать следует координаты только концевых масс. В состоянии покоя принимают нулевые значения их координат, а при равномерном движении — постоянные значения скоростей этих масс.

Нулевые собственные значения характерны для механических систем с дифференциальным соединением сосредоточенных масс и систем без реактивных элементов.

Рассмотрим методику решения систем алгебраических уравнений полу- определейных объектов. В матричной форме математическая модель объекта имеет вид

где atj — элементы матрицы Якоби системы алгебраических уравнений; xt — фазовые координаты; bt — элементы вектора внешних воздействий.

Предположим, что из общего количества п собственных значений матрицы Якоби т собственных значений оказались равными нулю. Зафиксируем т концевых масс с номерами индексов kq, q - 1,пг. Математическая модель при этом преобразуется путем обнуления всех элементов столбцов и строк матрицы Якоби, содержащих диагональные элементы akk. Диагональные элементы akk при этом сохраняются. При анализе состояния покоя полагают xk =0 и обнуляют элементы вектора внешних воздействий этих строк: bk - 0. Если рассматривается равномерное движение, то принимают xk = const, a bk вычисляют из соотношения - akkxk.

Пусть т -1, а одна из концевых масс имеет номер индекса 2. Примем Х2 - const; ai2 =0; a^j =0; i,j- 1 ,п; 022 ^0; =0 (состояние покоя). Уравнения

(8.69) после этих преобразований принимают вид

Рассмотрим типовые примеры анализа статических состояний линейных технических систем.

Пример 8.3. Определить начальные и конечные значения фазовых координат механической вращательной системы, необходимые для моделирования переходных характеристик и определения показателей качества переходных процессов. Динамическая модель системы приведена на рис. 8.2. Это простейшая модель трансмиссии машины, в которой энергия двигателя используется для ее перемещения и привода рабочих органов (например, трактор с приводом рабочих органов агрегатируемой сельхозмашины через вал отбора мощности, зерноуборочный комбайн, машина для уборки улиц).

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

Энергия, необходимая для функционирования рассматриваемой механической системы, подводится от внешнего источника (двигателя) с параметрами Мв1, сох и используется для совершения полезной работы потребителями энергии с параметрами Мв3, со3 и Мв4, со4. Часть энергии двигателя затрачивается на компенсацию внутренних потерь системы, характеризуемых моментом сопротивления Мв2. Ниже будет показано, что рассматриваемая система по л уопре деленная. Для обеспечения равновесного статического состояния такой системы значения моментов всех внешних воздействий должны быть согласованы между собой. Поэтому будем полагать заданными значения моментов Мв1, Мв2, Мв4, а Мв3 подлежит определению.

Динамические процессы системы описываются следующими уравнениями:

Для решения системы уравнений (8.71) с целью получения переходных характеристик необходимо задать начальные значения угловых скоростей юг0, i-l,n, начальные и конечные значения моментов в упругих элементах Му,о и MyjK, j = 1, N, где п — число степеней свободы системы, равное количеству сосредоточенных масс; N — количество упругих элементов.

При определении со,0 рассматривают либо состояние покоя системы, либо равномерное движение. В первом случае принимают сц0 =0, i - 1,п, во втором используется система уравнений (8.72). В статическом состоянии механической системы моменты диссипативных элементов Мду =0, j = 1,N, так как отсутствуют относительные скорости перемещений сосредоточенных масс. Тогда уравнения (8.72) получают вид

В системе уравнений (8.73) четыре неизвестных ю;, а уравнений всего три. Для их решения необходимо задать значение одной из угловых скоростей, например угловой скорости вала двигателя СЩ. Задавать угловые скорости рекомендуется концевым массам.

При определении моментов упругих элементов Муу и угловых координат сосредоточенных масс ф, необходимо, используя уравнения (8.71), получить математическую модель статического состояния системы. Для механических систем обычно рассматривается состояние покоя.

Решение задачи статики возможно двумя способами.

Первый способ заключается в составлении двух систем алгебраических уравнений. Одна из них — для определения Муу-, а вторая — для ф*. Первая получается на основе топологических уравнений, входящих в систему (8.71), полагая в них dcOi/dt =0, i = l,n, и учитывая, что в статическом состоянии Mpj = 0, у = IjV:

Так как значения Мв1, Мв2, Мв4 заданы, то неизвестными в уравнениях (8.74) являются Му1, Му2, Му3 и Мв3.

Для получения второй системы уравнений используются компонентные уравнения, входящие в систему (8.71), но только в интегральной форме. Интегрируя эти уравнения и учитывая, что Му1, Му2 и Му3 определены при решении уравнений (8.74), получаем

Решение этой системы уравнений осуществляется так же, как и системы (8.73). Так как число искомых неизвестных ф, превышает количество уравнений, то необходимо задать значение угла поворота одной из концевых масс (ф4, ф3 или ф4) и определить значения угловых координат остальных масс в начальном ф;0 и конечном ф статических состояниях системы.

Второй способ заключается в подстановке Му/ из уравнений (8.75) в топологические уравнения системы (8.71). В статическом состоянии системы топологические уравнения принимают вид (8.74). В результате получаем следующую систему уравнений:

Для проведения анализа состояния покоя объекта используем уравнения (8.76). Примем следующие значения его параметров: сх = 300; с2 = 1200; с3 = 500 (Cj в Н м/рад); u^i = 2; U22 = 2,5; И32 = 4; п43 = 1,5; Г|г = 0,95; i = 1,4. Внешние воздействия: Мв10 - 200; Мв1к - 300; Мв2 - 50; Мв4 - 375(Мв/ в Н • м).

При этих данных собственные значения матрицы Якоби: А,1=-18 462; Х2 = -2243; ^.3 = -746; ?i4 = 0. Наличие нулевого собственного значения означает, что система полуопределенная, поэтому для решения уравнений (8.76) необходимо лишить ее избыточной степени свободы. Зафиксируем положение сосредоточенной массы J3. Выполнив преобразования уравнений (8.76), предусмотренные изложенной выше методикой, получаем

Очевидно, что ф3 =0. Собственные значения матрицы Якоби этой системы уравнений: = -18 240; Х2 = -2405; Х3 = -771; А,4 = -34,9. Собственные значения матрицы Якоби системы алгебраических уравнений содержат только вещественные части. В данном случае все они отрицательны, следовательно, система устойчивая.

Решив систему уравнений (8.77), определим искомые значения угловых координат сосредоточенных масс системы срг, i - 1,п. Затем, используя уравнения (8.75), вычислим значения Муу, j = 1,N, а также неизвестное внешнее воздействие Мв3 = Му2Цз2Лз2-

Получены следующие значения ф/} Муу и Мв3:

  • ? в исходном состоянии равновесия (при Мв10 =200 Н м) фх = 1,32812; Ф2 =0,33073; ф4 =—0,13039; Му 200; Л7у2 — 158,75; ^уз —263,158; Мв3 = 603,25 (фг в рад; Myj и Мв3 в Н • м);
  • ? в конечном состоянии равновесия (при Мв1к =300 Н • м) (pj =3,54167; Ф2 — 1,27083; ф4 — 0,49635; Л7у4 =300; Л7У2 =610; Л7у3 — 263,158; Л7в3 — 2318.

Пример 8.4. Определить значения фазовых координат в состоянии покоя технической системы, динамическая модель которой дана на рис. 6.6, а, принимая два различных уровня внешнего воздействия на сосредоточенную массу т2 (FB20 и FB2к).

Отличительной особенностью этой системы является возможность возникновения виртуальной связи на массу т2.

Математическая модель технической системы приведена в параграфе 6.5. Она описывается уравнениями (б.ЗО)-(б.ЗЗ). Полагая dvi/dt-0, dXi/dt = 0 и Н =0 и подставляя значения Fyj из выражений (6.32) в уравнения (6.31), получаем следующую систему линейных алгебраических уравнений относительно неизвестных геометрических координат Xi".

Если при решении системы уравнений (8.78) при заданных значениях FBi и Fb2 удовлетворяется условие отсутствия виртуальной связи, т.е. О < х2 < L2b, то полученные результаты принимаются в качестве решения. Найденные значения х10, x20 и х, хзатем подставляются в выражения (6.32) и определяются Fyl0, Fy20 и Еу1к, Fy2K.

Если же при решении уравнений (8.78) условие (6.29) окажется нарушенным, то это означает, что на массу наложена виртуальная связь. Тогда при условии х2 <0 принимают х2 =0, а при х2 х2 = L2b, исключают из системы (8.78) второе уравнение, описывающее равновесное состояние массы т^, на которую наложена виртуальная связь, и решают оставшиеся уравнения с принятым ограничением на координату х2. Затем вычисляют значения Fyj по формулам (6.32).

Примем xi=x2=0 при FBl -Fb2 -0. Параметры системы: =200 Н/м; с2 =100 Н/м; L2b= 0,055 м.

Используя уравнения (8.78), определим координаты масс ;c10, х20 и х, х при Ев10 =2Н,Ев1к =6НиЕв2 = 2 Н. В результате получаем: х10 =0,02м;х20 =0,04м; х =0,04м; х = 0,06м.Таккак х >Ь, то следует принять значение х - L2b и по первому уравнению системы найти значение хы. В данном случае х =0,0383 м. Тогда получаем следующие значения усилий в упругих элементах: Еу10 =4 Н; Fy20 =2 Н; Еу1к =7,67 Н; Еу2к=1,67Н.

Пример 8.5. Определить распределение температуры по толщине стенки некоторого технического объекта, левая граничная поверхность которого имеет постоянную температуру Тв1 (граничное условие I рода), а правая охлаждается водой, имеющей постоянную температуру Тв2 =20 °С (граничное условие III рода). Стенка двухслойная: материал первого слоя — металлокерамика типа МК-5, а второго — сталь. Толщина первого слоя 1 мм, а второго 2 мм. Используя равномерную сетку, выделим 6 дискретных элементов. Динамическая модель теплового объекта аналогична представленной на рис. 5.17, а. Параметры дискретных элементов приведены в табл. 8.1.

Таблица 8.1

№ слоя

Геометрические размеры слоев

Теплофизические коэффициенты слоев

Плотность материала р, кг/м

Толщина 1, м

Площадь А, м2

С, Дж/(кг • К)

X, Дж/(с • м • К)

1

0,5 10“3

0,04

638

15,5

6200

2

Тоже

Тоже

638

15,5

6200

3

»

»

510

43,6

7800

4

»

»

510

43,6

7800

5

»

»

510

43,6

7800

6

»

»

510

43,6

7800

Математическую модель динамических процессов в тепловой системе получаем на основе уравнений (5.97) и (5.100), составленных для объектов с граничными условиями I и III рода соответственно:

где T;, i- 2,6— температура внутренних слоев; Тв1 — температура первого (левого) граничного слоя: Тх = Тв; Тв2 — температура внешней среды, соприкасающейся с граничной поверхностью шестого (правого) слоя стенки; стh ? = 2,6— теплоемкость дискретных слоев; цт/, ?=1,5— коэффициент теплопроводности при индуктивном теплообмене между слоями стенки; цкт2 — коэффициент теплопроводности при конвективном теплообмене шестого (граничного) слоя стенки с внешней средой.

Значения параметров cTi и цт, вычисляются по формулам (5.88) и (5.89), а параметр Цкт2 — по формуле

где а — коэффициент теплообмена внешней среды со стенкой объекта через конвекцию; Ав2 — площадь граничной поверхности теплообмена.

Графики распределения температуры Tпо толщине стенки

Рис. 8.3. Графики распределения температуры Tt по толщине стенки:

1 — при T*i =100 °С; 2 — при T*i =120 °С; п — номер дискретного элемента стенки

При охлаждении потоком воды а = 500...20 000 Вт/(м2 • К). Примем а = 5000. В результате параметры тепловой системы имеют следующие значения: рт1 = рт2 = 1240 Вт/К; Мтз» •••» Ит5= 3488 Вт/К; цкт2 = 200 Вт/К; ст1 = ст2 = 79,1 Дж/К; ст3, ...,ст6 = 79,6 Дж/К.

Используя систему уравнений (8.79), построим математическую модель анализа статических состояний теплового объекта, полагая dTJdt = 0, i -2,6 и опуская общие делители уравнений cTi, получаем:

Решение этой системы уравнений осуществлялось методом Гаусса. На рис. 8.3 приведены графики, иллюстрирующие распределение температуры по толщине стенки при двух постоянных значениях температуры Тв1 на левой граничной поверхности: 100 и 120 °С. Собственные значения матрицы Якоби системы уравнений (8.81): -12 044; -7476; -3372; 1775; -181. Число обусловленности матрицы Якоби ц = 66,541. Отрицательные вещественные части всех собственных значений свидетельствуют об устойчивости статического состояния системы.

 
Посмотреть оригинал
< Пред   СОДЕРЖАНИЕ   ОРИГИНАЛ     След >