УЧЕНЫЕ ЗАПИСКИ ЦАГИ
Том XXVI 1995 М 1-2
УДК 629.735.33.015.4.027.24-567 629.7.015.4.023
ПОШАГОВО-ИТЕРАЦИОННЫЙ МЕТОД РАСЧЕТА НАПРЯЖЕННО-ДЕФОРМИРОВАННОГО СОСТОЯНИЯ БАЛКИ ПЕРЕМЕННЫХ КРИВИЗНЫ И ПОПЕРЕЧНОГО СЕЧЕНИЯ С УЧЕТОМ ГЕОМЕТРИЧЕСКОЙ НЕЛИНЕЙНОСТИ
А. И. Залюбовский, С. Р. Саркисов, П. А. Яковлев
Рассмотрен пошагово-итерационный метод расчета геометрически нелинейного напряженно-деформированного состояния (НДС) балки переменных кривизны и поперечного сечения с учетом действия инерционных нагрузок или собственного веса. Проведен сравнительный анализ результатов расчета с точными аналитическими решениями, с результатами экспериментов и расчетов конечно-элементных моделей.
Показаны устойчивая сходимость и высокая точность предложенного метода для расчетов жесткостиых характеристик рессорных стоек шасси легких самолетов.
В настоящее время в связи с интенсивным развитием малой авиации в России и широкою применения рессорных стоек шасси, представляющих собой плоскую криволинейную балку с переменными по длине радиусом кривизны, формой и площадью поперечного сечения, появилась необходимость определения НДС рессоры в геометрически нелинейной постановке. •
На этапе проектирования таких стоек шасси для расчетов посадочных нагрузок необходимо иметь достоверную информацию по ожидаемым жесткостным характеристикам, а именно диаграмму статического обжатия. На этой стадии, как правило, нет возможности определить ее экспериментально, и для некоторых типов рессор, встречающихся на практике, ошибка при определении нагрузок достигает 80%, которая обнаруживается на стендовых испытаниях, так как линейные модели расчета таких рессор по схеме криволинейной балки [1] приводят к большим погрешностям, обусловленным изменением формы конструкции в процессе нагружения. Аналитические методы решения [2] громоздки, сложны и пригодны лишь для очень простых
^«1 Ру^О
Рис. 1
моделей. Кроме того, аналитическое решение с изменением исходной геометрии рессоры в процессе проектирования, а также при переменных форме сечения и радиуса кривизны становится практически невозможным.
Из существующих программ расчета методом конечных элементов не все позволяют проводить нелинейный анализ, а те, которые позволяют, нуждаются в задании большого объема исходных данных, имеют низкую точность и требуют большого времени решения, которое превышает линейные решения в 100 раз [3], что подразумевает наличие вычислительной техники более высокого уровня. Например, в работе [4] рассмотрен криволинейный балочный конечный элемент, позволяющий учесть геометрически нелинейные деформации при решении задач статики и динамики конструкций. Однако значительное несовпадение результатов расчета с точным решением на примере простой модели, отсутствие сравнения с экспериментальными данными, неучет инерционных нагрузок, а также перечисленные выше недостатки не позволяют определить область применимости метода конечных элементов. .
Актуальность проблемы и привела к созданию метода, лишенного рассмотренных недостатков. Предлагаемый метод прост и нагляден, имеет высокую точность, устойчивую сходимость и требует малого времени решения на компьютерах низкого класса.
Рассматриваемая балка является балкой переменного сечения, ось которой имеет переменную кривизну в плоскости (рис. 1). Используя известное выражение для приращения кривизны балки:
1
записываем уравнение изогнутой оси балки переменной кривизны и поперечного сечения (в системе координат, связанной с неподвижным основанием) в виде
Предлагается следующий способ его решения. Разобьем ось балки на п участков, каждый из которых в общем случае остается криволинейным. Геометрия балки задается координатами узлов (у,-, Zi)- Для упрощения и повышения точности такого задания производится интерполяция координат узловых значений наиболее часто применяемыми кубическими сплайнами, которые менее подвержены осцилляциям и обеспечивают более хорошее качество аппроксимации по сравнению со степенными многочленами.
Кубический векторно-параметрический сплайн на каждом участке представляет собой многочлен с неизвестными коэффициентами. Краевые условия сплайна вместе с векторными алгебраическими уравнениями для нахождения коэффициентов сплайна образуют замкнутую систему уравнений, имеющую единственное решение.
Длину А/, участка выбираем при следующих допущениях: EJt, Mh rt = const. Углы отсчитываются от вертикали (ось Оу) против часовой стрелки. Такое же направление соответствует положительным изгибающим моментам и кривизне балки. Закрепление балки принимается консольным с податливостью заделки к3.
Известно, что прямолинейная балка постоянной жесткости EJ, длиной А/ и нагруженная постоянным изгибающим моментом М — = const (рис. 2), деформируется по дуге окружности радиусом р.
Кривизна оси балки в этом случае определится точным выражением
21^/2 г(.У) e(S) - r(S) ■ Е ■ F(S)
+
, 1 N(S)
r(S) F(S)
1 M
v p EJ
A(PM
z
Рис. 2
Угол поворота крайнего сечения балки (т. А)
Хорда дуги
Д h = 2р sin—
2
Высота сегмента, образованного изогнутой осью балки и хордой,
В случае когда балка в исходном состоянии (без нагрузки) уже имеет начальную кривизну с радиусом г = const (рис. 3), можно записать выражение для приращения кривизны балки [2]:
Второе слагаемое определяется наличием нормальных сил. Кривизна балки в этом случае определится как
где е — эксцентриситет нейтральной оси (рис. 4), возникающий вследствие нелинейного (гиперболического) закона распределения напряжений по высоте сечения балки.
Практически для любого типа сечений можно получить выражение для эксцентриситета. Для прямоугольного сечения [2]:
(1)
h
У
О
z
Рис. 3
ИЛИ
е »-
н1
12г
1 +
4 Н
15 4г
О выражение (1) стремится к обыч-
В [2] показано, что при 1 /г ному виду:
. 1 М
р " р " Ы'
Определим начальные параметры для каждого участка и узла балки. Величина хорды ДЛ/ участка и абсолютный угол ср/,/ наклона хорды
первого участка ДА,- = Дйо,- = , / = 1, я,
, ГлгИ
срАо! = аг&£ —±- +пт\ срЛо,- = фАо/_1 + Дфйо»;
ДфАо,- =агсг8
Д^_
Ау,-
■ arctg
Ду/+1
/ = 1, и,
где Ду, = у, - _ ь Дг, = £,■ - 2, - ь /и определяется с учетом положения
первого участка в квадранте осей координат; 0 — индекс, соответствующий параметрам недеформированной балки.
Если в качестве исходных данных имеются углы наклона касательных фо/ в точках (аппроксимация сплайнами дает такую возможность), то центральный угол дуги участка можно определить:
Дфо/ = дфм1 = Фо; - Ф0/-1 > * = 1»п-Кривизна участка определится как
д%
к ■ - — -кп - -
- . лФм0/ ’ 28111-----У-
/ = 1, и.
Если в качестве исходных данных имеются значения кривизны участков то центральный угол дуги участка будет равен
Аср Mt = Лфл/о,- = 2 arcsin
kri
ААо,-
, ' = 1, «>
где Аф^ш ^180°.
Угол наклона касательной в точке (у/, Zj)
А(Р Мої
Ф / = ФО / = Ф hi-----~---,
/ = 1, п.
Длина дуги участка будет
А/,- = Д/0/ = Дфм01 П . *' = 1, «•
Изменением длиной балки вследствие малых осевых деформаций можно пренебречь.
Параметры недеформированного состояния балки, приведенные выше, являются исходными данными для расчета ее поведения под нагрузкой.
Плечи сил Ру,Рг для расчета изгибающих моментов (см. рис. 1) на середине участка А// равны
т _ Zi +Zj-\ , _ Lyt ~ Zs---------X----+ Pi
1 - COS
г ... Уі-Уі-1..
~ys------5---+ Pi
1 - COS
АФMj 2
Дф Mt
с^Фа, ;
sin ФЛ,
Изгибающий момент на середине участка
Mi — Ру Lyi + Pi Lrf + Mx, і — 1, п,
а изгибающий момент в узле (у,-, Zi)
= Py(.Zs-Z,) + Pz &S ~ Уд + Мх, і = 1, и,
где (ys,Zs) — точка приложения сил Ру, Рг, Мх.
В общем случае каждая сила может иметь свою точку приложения. Нормальные и перерезывающие силы на середине участка выражаются как
Nt = Ру совфл. - Pz sin ф/j. ,
Qi = Ру sin Фhi + Pz COS<ph. , 1=1, n,
а нормальные и перерезывающие силы в узле (у„ z,)
Nj = Ру cos ф,- - Pz sin ф,-,
Qj = Ру sin ф,- + Pz cosф;-, / = 1, п.
Новая кривизна участка при действии изгибающего момента
kpi — к?/ + А кр(,
где
4 Mf kri , Nj .
Ak0j =—+ kri —, i = l, n.
pl eEFj EFj ’
Новый центральный угол дуги
= АФМ0, + Akpi А// * i = 1>п-
Новый абсолютный угол поворота сечения узла (yt, Zi)
ФО = АЧ>М0 = ^3 -Mo, i = о.
Ф/ = Ф»-1 + дФАо/ + А(РMi > * = !» я-
Новый абсолютный угол поворота хорды
АФ Mi . ,
ф*=ф/------у1-, i = l, л-
Новая величина хорды деформированного участка
. , - , . АФЛГ, . ,
А 1ц = 2kpi sin — —, i = l, п.
Новые координаты точек
У/ = У/ -1 + ДА, cos фаI,
Z, = 1 + ДА, sin фЛ/, /' = 1, и,
1 л „ при — = 0 имеем кусочно-линеиные участки, что позволяет значитель-
П
но снизить объем исходных данных, в случае упрощенных моделей, при очень незначительном снижении точности.
В ряде случаев, когда это требуется, несложно учесть влияние инерционной нагрузки (пу, nz, / = 1, и) или собственного веса балки и закрепленных на ней сосредоточенных грузов (пу = 1, пг = 0, i = 1, п) на НДС балки. Для этого к изгибающему моменту от внешних сил следует прибавить изгибающий момент от инерционной нагрузки или собственного веса, а к внешним силам прибавить инерционные силы
Mi = Mi + Mgj,
где изгибающий момент от инерционных сил на середине участка балки
71-1
= \ПУк+№ъ ~ LVc+0 + nZk+l(Lyi ~ Lyk+V>\
k=i
Изгибающий момент в узле
л-1
Mgi = ~^Gk+l • ~ Zi ~ LZk+1^ + nZ.k+1 (Уй ~ Уі ~ ^Ук+1 )]>
к=і
а вес участка балки
Gt = А1Г Ffj + m,g, i=l, n,
где mtg — вес сосредоточенного груза (считается закрепленным на середине участка); j — удельный вес материала балки; Ft — площадь сечения балки на м участке.
Нормальная сила на середине участка
Nt = Л} + Ngt, i=l, п.
Составляющая нормальной силы от инерционных нагрузок или веса на середине участка
П
Ngi = X • COS(PA, " nZk - Sin фАі ], і = 1,п, к=1
в узле балки
п
= ^Gk[nyk •cos<pi -nZk • sinсрЛ(], і = 0,и. fc=l
Перерезывающая сила на середине участка
Qt= Qi+ Qgi, і = 1, п.
Составляющая перерезывающей силы от инерционных нагрузок
П
Qgi = • sin9, + nzk ' c0ScPi ]» / = l>
k=1
в узле балки
П
Qgi = £Gk[«yk ■ sinф, + nlk ■ coscpl ], і = 0, и.
к=1
Нормальные напряжения определим с учетом гиперболического закона распределения напряжений по высоте сечения балки, поскольку
н
при существующей начальной кривизне (— < 10) линейный закон рас-
г
пределения дает существенную погрешность [1, 2]:
В случае наличия лобовой силы
му, ~ Рх[Ъ cos9a, + LZi sin Фа, ], ' = 1, я.
. Касательные напряжения определим по формуле Журавского, в частности, для прямоугольного сечения [1, 2]
30,
2В{ Hf
1-
4у
і = 1, п,
Зі>
1-
4дг
Д?
_ -^крг Т*У/ -
гккр
Фа,- - arctg
Угол закручивания рессоры
' М
где /х, , /кр/ — геометрические характеристики сечения; Я^кр —
момент сопротивления кручению; <7 — модуль сдвига материала рессоры.
Эквивалентные напряжения определим по третьей теории прочности:
<*экв, = >/ст? + 4т? ,
где X/ = т^. + ху1/ — для середины боковых сторон сечения; т,- = Хщ + Хщ — для середины верхней и нижней сторон сечения.
Сущность предлагаемого шагово-итерационного метода расчета заключается в том, что нагружение балки задается вектором внешних сил Р] в общем случае:
Р] = {^х] ’ Ру]> Рц ’ ^х] ’ {пУ } j ) •
Вектор внешних сил Р] может применяться с постоянным или переменным шагом всеми своими компонентами. На каждом шаге осуществляется процесс последовательных приближений, который ограничен требуемой точностью вычисления изгибающего момента в заделке балки е. Приращение компонентов вектора внешней нагрузки прекращается по критерию превышения напряжениями значения Стщ материала балки. Пошагово-итерационный процесс строится по следующему алгоритму:
— Определяются геометрические параметры недеформированной балки.
— Задается шаг по вектору внешней нагрузки.
— В случае превышения количества шагов по нагрузке пошаговоитерационный процесс завершается.
— Определяются изгибающие моменты и внутренние силы по геометрии балки на предыдущей итерации.
— Определяются геометрические параметры деформированной балки на текущей итерации.
— Определяются изгибающие моменты и внутренние силы на текущей итерации.
— Осуществляется сравнение изгибающих моментов в заделке на текущей и предыдущей итерации.
— В случае достижения заданной точности е определяются напряжения в сечениях балки, итерационный процесс завершается и вычисления продолжаются в п. 2 (при условии, что аэкв / < стш)-
— В случае недостижения заданной точности итерационный процесс возобновляется с п. 4.
— В процессе решения вся необходимая информация представляется в удобной для анализа форме.
/ 8> »* п2 из
1 2 3 1 2 3
0,15 -2,776 -2,775 -2,776 -0,441 -0,439 -0,439 6 2
0,3 -5,253 -5,247 -5,249 -1,710 -1,669 -1,700 6 2
0,45 -7,169 -7,160 -7,160 -3,649 -3,616 -3,616 6 2
0,6 -8,336 -8,333 -8,333 -6,020 -5,945 -5,946 6 2
0,75 -8,665 -8,694 -8,694 -8,532 -8,398 -8,400 6 2
0,9 -8,174 -8,281 -8,281 -10,88 -10,69 -10,69 6 2
1,05 -6,995 -7,231 -7,231 -12,80 -12,57 -12,57 6 2
1.2 -5,349 -5,758 -5,758 -14,07 -13,87 -13,87 6 2
1,35 -3,528 -4,114 -4,114 -14,58 -14,52 -14,52 6 2
1.5 -1,849 -2,546 -2,546 -14,34 -14,55 -14,55 6 2
1,65 -0,607 -1,264 -1,264 -13,50 -14,06 -14,06 6 2
1,8 -0,027 -0,405 -0,405 -12,33 -13,25 -13,25 7 2
1 — данные, полученные в работе [4]; 2 — данные по теории, приведенные в [4]; 3 — расчет по предлагаемому методу; п — количество итераций на каждом шаге; / — параметр,
„ М равный -==■
Ы
В некоторых случаях, когда итерационный процесс последовательных приближений затягивается, его можно существенно ускорить, если плечи действия сил Л» рг вычислять с учетом градиента деформаций, при У > 3:
^ (Л = Zs и -1) + ^ О ~ 1) О ~ 1) %; у я (Л = Уs 0' -1) + Д^ О -1) Дз'.у 0' -1)
А ^0-1)
А^(у -2)
ДгИ./ -1) = гяО) - -1) >
Ау3и-2)
Ау «г и -1) = Ху (Л ~ Уs и ~ 1) >
где у$, ^ — координаты точки приложения внешних сил; ^ < 1 — коэффициент гарантированного запаздывания деформаций, обеспечивающий апериодическую сходимость.
Для иллюстрации метода рассмотрено несколько примеров.
Пример 1 (взят из [2]). Определение силы Р, приложенной к концам вертикального гибкого стержня, приводящей к соприкосновению его концов. В силу симметрии рассмотрена верхняя половина стержня длиной /. Расчет проведен для балки прямоугольного сечения (Я = = 1,2 мм, В = 24 мм) длиной / = 750 мм, Е = 1,1 • 104 кГ/мм. Результаты расчета упругой линии без учета собственного веса балки в сравнении с [2] представлены на рис. 5. Относительная погрешность при определении силы Р составила 0,027%.
Проведен также расчет предлагаемым методом с учетом собственного веса пластины. Результаты расчета упругой линии на рис. 6 приведены в сравнении с результатами натурного эксперимента. Относительная погрешность определения силы Р составила 0,031%.
Пример 2 (взят из [4]). Рассмотрена задача изгиба консольной балки квадратного сечения (В = Я = 1 дюйм) длиной / = 12 дюймов, Е = 30 • 10б фунт/дюйм под действием изгибающего момента на конце. Балка была разбита на восемь участков.
8^, кГ/мм1
Результаты расчета в сравнении с результатами работы {4], полученными МКЭ в нелинейной постановке, приведенные в таблице, показывают высокую точность предлагаемого метода.
Пример 3. Рассмотрена задача определения диаграммы обжатия основной стойки шасси легкого самолета, представляющей собой арочную рессору переменной кривизны и поперечного сечения. Высота и ширина сечения по длине рессоры изменяются линейно в пределах: Н= (12 - 16) мм, В = (60 - 80) мм. Результаты расчета представлены на рис. 7 в сравнении с расчетом по линейной модели криволинейной балки [1] (рис. 8). Из анализа результатов следует, что для Р ошибка при применении линейной модели в случае определения обжатия составила 66,6%, а в случае определения напряжений 33%. Попытка расчета МКЭ с геометрически нелинейным анализом не привела к успеху в силу большого времени решения и низкой точности, что подтверждено в [3] на более простых моделях. Применение процедуры ускорения сходимости в данном примере дало снижение числа итераций на 30%.
Рассмотренные выше примеры наглядно демонстрируют высокую точность и хорошую сходимость для решения задач определения геометрически нелинейного НДС балки переменной кривизны и поперечного сечения. Кроме того, данный метод может быть успешно применен при определении НДС балок за пределами устойчивости (сжатоизогнутое состояние) практически для любых способов закрепления концов балки [5].
ЛИТЕРАТУРА
1. Тимошенко С. П. Сопротивление материалов. Т. 1.— М.: Наука.—1965.
2. Пономарев С. Д., Бидерман В. Л., Лихарев К. К. Расчеты на прочность в машиностроении. Т. 1.— М.: Машиностроение.— 1956.
3. Алберг Дж., Нильсон Э., Уолш Дж. Теория сплайнов и ее приложения / Пер. с англ.— М.: Мир,— 1972.
4. Зайцев С. Н. Криволинейный балочный конечный элемент, учитывающий геометрически нелинейные деформации // Ученые записки ЦАГИ.—
1991. Т. 22, № 5.
5. Ветчинкин В. П. Избранные труды. Т. 2.— М.: Изд-во АН СССР,- 1959.
Рукопись поступила 10/Ш1994 г.