УДК 624.04:517.926.7
ВРЕМЕННОЙ АНАЛИЗ УПРУГОПЛАСТИЧЕСКИХ КОНЕЧНОМЕРНЫХ ДИССИПАТИВНЫХ СИСТЕМ ПРИ НЕСТАЦИОНАРНЫХ ВОЗДЕЙСТВИЯХ
А.Н. Потапов
Изложен аппарат технической реализации метода временного анализа диссипативных конструкций из идеально упругопластического материала при нестационарном процессе. Рассмотрен пример неупругих колебаний плоского трехэтажного каркаса под действием кратковременной нагрузки большой интенсивности.
Введение
Временной анализ реакции конечномерных систем с внутренним трением представляет исключительно сложную и трудоемкую задачу, особенно в тех случаях, когда требуется проводить учет пластических деформаций в элементах конструкций при нестационарных воздействиях. Поэтому аналитические решения в подобных задачах крайне редки и связаны с упрощающими предпосылками. Упрощения относятся либо к выбору расчетных динамических моделей с числом степеней свободы не более одного-двух [1, 2], либо к выбору модели упруговязкого сопротивления, обычно соответствующей классическому пропорциональному (однородному) демпфированию [3].
В [4] изложены основные положения теории временного анализа, основанной на разработанных методах анализа матричных линейных и нелинейных уравнений. Эта теория обладает достаточно широкой областью приложений. Объектами исследования являются конечномерные (дискретные) диссипативные системы (ДДС). Учет внутреннего трения проводится на базе линейной модели вязкого сопротивления. Разрешающие уравнения как упругой, так и неупругой реакции строятся в замкнутом виде через интеграл Дюамеля, причем в задачах колебаний ДДС с нелинейной восстанавливающей силой закон деформирования подчинятся идеально упругопластической диаграмме Прандтля. Важно подчеркнуть, что эффективность подхода не зависит ни от типа демпфирования (однородное, неоднородное), ни от характера внешнего воздействия и соответствующего ему режима колебаний.
Теоретические основы
Общая идея подхода состоит в использовании теории промежуточных состояний (в рамках теории пластического течения). Весь процесс не-упругого анализа реакции разбит на такие интервалы, внутри которых динамические параметры системы неизменны. Переход из одного интервала в другой происходит при скачкообразном изменении параметров жесткости в критических временных точках соответствующих открытию или закрытию шарниров текучести. Поэтому на каж-
дом отдельном интервале используется схема упругого решения, а при переходе через критические точки решения на смежных интервалах стыкуются посредством начальных условий, обеспечивающих непрерывность движения упругопластической диссипативной системы. В результате весь неупругий процесс представляется рядом последовательно изменяющихся квазиупругих решений.
Условие равновесия динамической квазиуп-ругой системы на произвольном интервале времени t е [4 /w], где tb tn-i - время перехода в промежуточные состояния, выражается матричным уравнением движения
мп о + с т + K(mt) = Pit), (i)
где
P(t) = -R0(t,) + R\td + P(t).
Здесь: М= diag (ти ..., т„), С = Ст, Щ) = Щ)Т -матрицы инерции, демпфирования и жесткости соответственно; Л0(^)> R (*д - векторы составляющих восстанавливающих сил (предельных значений и остаточных усилий); Y(t), P(t) - векторы искомых перемещений и заданной внешней нагрузки (все обозначения приняты в соответствии с работой [4]).
Уравнения реакции системы на интервале t е [t„ ti+j] (г = 0, 1,...) имеют вид:
Y(t) = fy?Zk (0, Щ = &кЩ% (0, (2)
/Ы1 *=1
г*(0 = ф£(*-О^[-ЗД+й]+
(3)
где Uk--U\ = M(S*-Sj) (/= 1,2;/**),Ф*(0= eSk‘;
Y0, 70 - векторы начальных условий (перемещений и скоростей);
= (4)
‘i
Yk{>0 = |ф1 ~ *)?(*) dr ■ (5>
h
Уравнения реакции (2), (3) обобщают интеграл Дюамеля на область физически нелинейных
систем, позволяя на квазиупругом интервале определять значения перемещений и скоростей узлов упругопластической конструкции при произвольной нагрузке. Матрицы Бк (к = 1,2) суть корни характеристического уравнения движения собственных форм (матричного квадратного уравнения): МБк2 + СБк + К(ф = 0 (/' = 0, 1,...), которое вытекает из однородного дифференциального уравнения, соответствующего (1).
Технический анализ проблемы
Проблема реализации разрешающих уравнений (2), (3) требует построения аналитических зависимостей для интегралов (4), (5). Это наиболее предпочтительный вариант, хотя не исключена возможность использования численных процедур. В любом случае сложности временного анализа переносятся на область технической реализации задачи. Как показали исследования, формализация интегралов (4), (5) на базе строгого математического аппарата оказывается чрезвычайно плодотворным шагом, что позволяет задействовать весь технический арсенал методов анализа сложных ДДС, ориентируя их на современную вычислительную технику и делая доступными для практического применения.
Вычислительная формула интеграла (4) зависит от условий состояния квазиупругой системы [4]. При условии невырожденного состояния (скг 5* Ф 0) интеграл (4) представляется конечной формулой ./*(/) = 05*Т)"’1[Ф*Т(/ -*/)-£]. В случае вырожденного состояния при полной диссипации системы (с1е15* = 0, сЫ 5; ф 0), интеграл (4) имеет вид функционального ряда
„0+1_/р+1\
Вычислительная формула (6) может быть преобразована к более простому виду, не содержащему матричных функций-сомножителей, что позволит эффективнее осуществлять практический анализ уравнений (2), (3). Покажем это.
Поскольку для фундаментальной матрицы Ф*(0 справедливо разложение в бесконечный ряд по степеням матричного аргумента %
00 №)*
р=о Р-
то проведем почленное перемножение элементов ряда в (6). Выполнив группировку слагаемых перед одинаковыми степенями матричных корней 5*, получим скалярные функции ар(г) (р = 0, 1, ...) в виде степенных рядов:
р {-\)я (1ч+х
^ = Ъ (4+1)!(р-9)! ' (7>
в результате чего формула (6) запишется так
\Р
•4(0= Е(5*) ар({)-
р=о
(8)
Проведем более детальное исследование функции (7).
Вынесем величину 1/(р + 1)! за знак суммы. Замечая, что факториальный множитель, входящий в <зг-й член ряда, является биномиальным коэффициентом многочлена (і - і,) р+1, то есть
с9+1-
'~'Р+1 —
-?+1
где Ср+\ - число сочетаний изр + 1 по д + 1, будем
иметь
«,(')- ^ і (-о 'с;:; г) .
В развернутой форме записи получим многочлен со степенями 1Р~Ч^+Х
•••-(-■Г1 с£«*/’ -(->)'.
где
Ср°+1 =4+1 -С2р+1 + ...+(-1)/’“ СрР+1 +(-1)Р С$ .
Далее учтем некоторые свойства биномиальных коэффициентов:
р+1
сД, = с;+7, Ср+\ = Ср+1 =1, X (-1) =0.
9=0
согласно которым
1 р+1
В соответствии с формулой бинома Ньютона имеем
р+1
9=0
Сопоставив правые части двух последних выражений, приходим к очевидному соотношению яДО = (? - ?() р+1/(р +1)!. Отсюда, учтя (8), найдем окончательно
ч’МГ1
(9)
Таким образом, получено аналитическое выражение интеграла (4), обеспечивающее разложение У*(0 в одинарный ряд с более простой (по сравнению с (6)) структурой входящих в него членов. С технической точки зрения организация вычислений ряда (9) намного проще, нежели для ряда (6), поскольку при программировании задачи для многочленов вида (/ - ?,)/Я’1 гораздо удобнее использовать схему вложенных умножений.
Вычисление интеграла (5) связано с характером внешнего воздействия. Рассмотрим действие кратковременной нагрузки по синусоидальному закону Р(г) = втб /-Р0, где 0 = Еп!1а, Е - единичная матрица, - время действия нагрузки; Р0 - вектор амплитуд внешних сил. Тогда из (5)
П (0 = 4 (Оро> 4 (О = (' -T)sm(0T)dr . (10)
Техника вычисления интеграла 1$) изложена в [4] и сводится к исследованию непрерывного уравнения Сильвестра (t е [г„ ?,+1] (/' = 0, 1,...))
+ (11) правая часть которого определяется выражением Fiit) = Ф*т(*-#/)[3*т sm0/, + 0 cos0/,] -- sin0r + 0 COS0?].
Здесь матрицы Sk (k = 1, 2) являются постоянными величинами на интервале времени t е [4 tM\. Для скалярной матрицы 0 из (11) вытекает обычное матричное линейное уравнение
[(S'/)2 + е2ш=+ е2] = Ш
из которого непосредственно следует
т = [№т)2 + QVFdt) = Fit) mTf + 02]-1.
Пример расчета
Рассмотрим пример упругопластических колебаний 3-этажного каркасного здания с металлическими колоннами и железобетонными перекрытиями, находящегося под действием кратковременной нагрузки большой интенсивности синусоидальной формы P(t) - sin (Ert/ta)P0, где ta = 0,8 с -продолжительность импульсной нагрузки. Вектор амплитуд динамических сил равен Р0 = [250, 30, 170]т (кН). Расчетная схема здания имеет вид сдвиговой модели (рис. 1), согласно которой ригели считаются абсолютно твердыми (заключающими в себе массы этажей) и жестко соединенными с невесомыми колоннами. Поперечные сечения колонн (сталь С255) нижнего этажа - сквозные, верхних этажей - сплошные.
На основе инерционных и жесткостных характеристик конструкции сформированы следующие матрицы: М= йіщ (2,97,2,97,2,1) (кНс2/см),
К =
289,96 -99,62 0
-99,62 297,97 -198,36
0 -198,36 198,36
93,88 -32,04 0
С = -32,04 95,17 -64,16 -10"2 (кНс/см).
0 -64,16 64,97
Демпфирование является неоднородным и учитывается по модели [4]: С = (КТ + ТК)12, где Т= <йщ (Гь /2, *3); Ь = (й'лгХ»?;/ /-й)0’5, т,, гй - диагональные элементы матриц Т, М, К; 3= 0,1 - логарифмический декремент колебаний.
Значения предельно упругих перемещений концов колонн по этажам здания составили У0= [1,57, 1,72, 1,09]т см. Вследствие высокой скорости деформирования (ё = 3-10“3 > 10"4 с-1), реакция системы на интервале * € [0, /а] вычислялась с учетом изменения прочностных свойств материала. Превышение динамического предела текучести [2] над статическим обусловило возрастание предельно упругих перемещений на активном участке нагружения на 5 %.
Численная реализация упругопластической задачи проводилась в универсальной математической системе МаЙ_,аЬ. Шаг временного анализа был принят А? = 0,01 с. При поиске критических временных точек минимальный шаг доходил до 10-7 с. Основные этапы временного анализа по вычислению параметров неупругой реакции приведены в табл. 1-3.
Табл. 1 содержит информацию о времени образования и закрытия пластических шарниров в несущих элементах этажей здания, включая повторную текучесть и разгрузку в колоннах второго этажа (см. второй полуцикл диаграммы деформирования на рис. 2).
Таблица 1
№№ этажей Начало текучести ti, С Начало разгрузки h, с
здания 1-йп/ц 2-йп/ц 1-йп/ц 2-йп/ц
1 0,4235 0 0,8548 0
2 0,6039 2,1113 1,1570 2,1174
3 0 0 0 0
(кН/см),
Рис. 2. Диаграммы деформирования «восстанавливающая сила - относительное перемещение» конструктивных элементов этажей здания
В табл. 2 представлены значения максимальных относительных перемещений каждого этажа здания (вектор Y^), остаточных относительных и абсолютных перемещений (векторы У* =Ymsx-YQ, У*), накопленных в процессе реагирования системы.
В табл. 3 дана некоторая информация о промежуточных состояниях квазиупругой системы на интервалах t е [/,, ?/+i] (переходы из одного промежуточного состояния в другое происходят в моменты времени г,).
На диаграмме «восстанавливающая сила -относительное перемещение» (рис. 2) изображены жесткости этажей конструкции в упругопластической стадии. Причем диаграмма деформирования 2-го этажа состоит из двух полуциклов, каждый из которых содержит площадку текучести. Первый полуцикл диаграммы построен при динамическом предельно упругом перемещении, равном Уо2 = 1,05-1,72 = 1,806 см, так как начало текучести (t2 = 0,6039 с) по времени оказалось в активной зоне нагружения. Высота площадки текучести соответствует предельному значению восстанавливающей силы Rm - 179,9 кН. Второй полуцикл построен при Rq2 - -171,3 кН и статическом пре-
дельно упругом перемещении у02 = 1,72 см, поскольку начало текучести по времени (/5 = 2,1113 с) относится к этапу свободных колебаний здания. Вследствие необратимых деформаций колебания этажей происходят со смещением относительно первоначального недеформированного положения здания.
На рис. 3 приведены осциллограммы параметров неупругой реакции здания при нестационарном процессе от импульсного воздействия. Колебания второго этажа совершаются с учетом гистерезисного характера диаграммы деформирования (рис. 3 а). Движение этажа на отрезке времени ? е [0, 1,492 с] соответствует первому полу-циклу диаграммы. Асимптоты аъ Ь\ на графиках фиксируют уровень остаточных смещений (абсолютный и относительный), устанавливающий новое положение равновесия этажа при его колебаниях. Начиная с момента времени 4, происходит понижение этого уровня на величину 0,00653 см (асимптоты а2, Ь2) вследствие разгрузки и смены знака пластических деформаций при реализации второго полуцикла диаграммы деформирования (? > 1,492 с). При I = г3 абсолютные перемещения второго (и третьего) этажа возрастают на величину
Таблица 2
этажей здания Наибольшие перемещения ^пих > СМ Относительные и абсолюта в к-х полуциклах диаг] ле остаточные перемещения рамм деформирования
Y*, см Y*, см
1-й п/ц 2-й п/ц за оба п/ц 1-й п/ц 2-й п/ц за оба п/ц
1 2 3 3,5662 6,2685 < 1,1445 1,9177 4,4626 0 0 -0,0065 0 1,9177 4,4560 0 1,9177 6.3803 6.3803 0 -0,0065 -0,0065 1,9177 6,mi 6,3131
Таблица 3
Параметры реакции j Время tj (с) начала открытия или закрытия пластических шарниров
t\ =0,4235 h=0,6039 ta=0,8000 h =0,8548 «*=1,1570 *5=2,1113 4=2,1174
см 1 2 3 1,6485 0,0778 0,4853 2,4953 1,8060 0,5805 3,4223 3,5129 0,6604 ЗМІЇ 4,1001 0,4130 2,4418 6,1826 0,5814 1,3901 2,7425 -0,3579 1,4108 2,7360 -0,3580
Ї‘(Ь)- см/с 1 2 3 5,2486 7,3361 -0,5868 5,1660 9,1227 2,3307 2,3943 10,1631 -3,6806 0 11,1249 -4,9809 -3,2352 0 3,6847 3,9375 -2,1533 0,4419 2,8848 0,0036 -0,4571
Го=М> см 1 2 3 1,6485 1,7263 2,2116 2,4953 4,3013 4,8818 3,4223 6,9352 7,5955 3,4877 7,5878 8,0008 2,4418 8,6244 9,2058 1,3901 4,1326 3,7747 1,4108 4,1468 3,7889
Y0=(y0j), см/с 1 2 3 5,2486 12,5847 11,9979 5,1660 14,2887 16,6194 2,3943 12,5574 8,8768 0 11,1249 6,1440 -3,2352 -3,2352 0,4495 3,9375 1,7842 2,2261 2,8848 2,8883 2,4312
1 313,778 133,869 133,869 -179,909 0 171,341 0
Ro, кН 2 0 179,909 179,909 179,909 0 -171,341 0
3 0 0 0 0 0 0 0
1 0 0 0 365,012 -79,536 -79,536 -78,885
Д\кН 2 0 0 0 0 444,548 444,548 443,897
3 0 0 0 0 0 0 0
Уі{і), см
10
8
6
4
2
0
-2
-4
-6
-8
г.д.гй1 1
/ * Йч /4-л ’ч\ /А*
.** *Х 1 А* а ;Л/' '! і \ *2 , V / »
и/\ \ 1 1 і і і 1 1 1
1 1 і і 1 » 1 1 1
■ л 1 г 1 \ і
% \ /
н V и І6
О
і, с
УА 0. см 10 8 6 4 2 О -2 -4 -6 -8
О
/\ г.ічу --Ду-
1 \
-ІС *
V >
к 1 1 1
і \ у
V
г, с
а)
б)
м/с
В)
*,С
ЩкН 8
-4
-6
А
А Л і \
м й А Й&Л л Л Л- гТ\ л
тжтмттш И1НШШ
2М т V *4
\ 1ч \1
і, с
О
2 3
г)
Рис. 3. Осциллограммы параметров упругопластической реакции каркасного здания (цифры на графиках - номера этажей; пунктир - упругая реакция здания): а - перемещений 2-го этажа (сплошная линия - абсолютные; штрих-пунктир - относительные; точки - статические перемещения); б - перемещений всех этажей; в - скоростей; г - диссипативных сил
1,9177 см. Это вызвано развитием и последующим накоплением пластических деформаций в концевых сечениях колонн первого этажа (асимптота сх).
Более общую картину упругопластической работы этажей конструкции иллюстрируют осциллограммы перемещений (рис. 36), скоростей (рис. Зв) и сил внутреннего трения (рис. Зг).
Экстремальные пики на осциллограмме диссипативных сил образуются в моменты времени, соответствующие появлению текучести. При этом наибольшие значения силы (положительная 1*1 = 6,65 кН при t = 0,86 с и отрицательная Рх = -5,26 кН при * = 0,98 с) приходятся на интервал текучести t е [?2, когда происходит рост пластических деформаций в колоннах второго этажа.
Таким образом, полученные результаты свидетельствуют о высокой эффективности метода временного анализа в задаче колебаний конструк-
ций с учетом вязких и пластических свойств материала при нестационарных воздействиях.
Литература
1. Гольденблат И.И., Николаенко Н.А. Расчет конструкций на действие сейсмических и импульсивных сил. — М.: Госстройиздат, 1961. - 320 с.
2. Попов Н.Н., Расторгуев Б.С. Динамический расчет железобетонных конструкций - М.: Стройиздат, 1974. —208 с.
3. Чернов Ю.Т. Исследование нелинейных систем при кратковременных динамических воздействиях // Строительная механика и расчет сооружений. - 1982, ИяЗ.—С. 35-40.
4. Потапов А.Н. Динамический анализ дискретных диссипативных систем при нестационарных воздействиях. - Челябинск: Изд-во ЮУрГУ, 2003. -167 с.