УЧЕНЫЕ ЗАПИСКИ Ц А Г И
Том IV 197 3 № 3
УДК 629.735.33.015.4:539.4.012
МЕТОД РАСЧЕТА УПРУГИХ ВОЛН В НЕОДНОРОДНЫХ
БАЛКАХ
И. А. Ляховенко
Для задач о поперечных волнах в неоднородных по длине балках получены уравнения типа уравнений Тимошенко и описан алгоритм их численного решения на основе метода характеристик при различных комбинациях граничных условий. В качестве иллюстрации даны примеры расчета шарнирно опертой на двух опорах конической балки и однородной балки, опертой на одном конце и свободной на другом, при действии распределенной динамической нагрузки.
Поперечные колебания балок исследуются чаще всего на основе инженерной теории изгиба, причем при выводе основного дифференциального уравнения колебаний (уравнения Бернулли— Эйлера) в поперечных сечениях балки исключаются податливость на сдвиг и продольные силы инерции. При достаточно длинных изгибных волнах указанные факторы не играют существенной роли, но при уменьшении длины волны их влияние становится решающим [ 1 ] — [3], [9].
Недостатки теории Бернулли —Эйлера были устранены Тимошенко путем учета деформаций сдвига и инерции вращения элементов балки.
Область применимости уравнений Тимошенко обследована теоретически Е. С. Сорокиным и А. С. Архиповым [4] с помощью точных уравнений теории упругости (рассматривалась задача о свободных поперечных колебаниях балки как плоская задача теории упругости).
Надежные эксперименты с высокочастотными колебаниями были проведены методом голографической интерферометрии [5]. Показано, что найденные формы и частоты колебаний хорошо согласуются с теорией балки Тимошенко вплоть до длин полуволн порядка высоты балки. Поэтому можно ожидать, что и для балок переменного сечения, рассматриваемых в настоящей статье, теория Тимошенко обеспечивает приемлемую для практических целей точность определения внутренних усилий и перемещений.
ОСНОВНЫЕ УРАВНЕНИЯ ДЛЯ НЕОДНОРОДНЫХ БАЛОК*
Обозначим полный прогиб балки в сечении с координатой х в момент времени / через ии(х,і), составляющую прогиба от изгиба — через 1юп(х,(), составляющую прогиба от сдвига — через {х, і), тогда можем записать: XV (х, і) = ж>и (х, і) + хгіс (х, і).
0) & > г)
Угол р между касательной к изогнутой оси балки и осью х без учета сдвига определяется как дхюп/дх. Соотношения для изгибающего момента и перерезывающей силы в сечении х имеют
Здесь Е1(х) — жесткость балки на изгиб; (х) — жесткость на сдвиг (С —модуль 2-го рода); £ — коэффициент, зависящий от формы сечения балки.
В работе [8] рассмотрена задача определения & для разных форм поперечного сечения. Оказалось, что значения £ лежат в сравнительно узком диапазоне. Неопределенности в выборе /г, заключают авторы работы [8], вообще не так серьезны при той степени приближения, которую имеет модель сдвига балки.
Здесь будет принято, что к является некоторой слабоменяю-щейся функцией координаты х.
За положительные направления х, чю, УИ, (3, д* примем указанные на фиг. 1 (в левом верхнем углу). Уравнения равновесия имеют в этом случае вид
Здесь т(х) — рЕ(х) — погонная масса балки; /р (х) ^ рі(х) — массовый момент инерции слоя балки единичной длины.
* Под неоднородными балками, как и в работе (6], будем подразумевать балки с переменными по их длине жесткостными и массовыми характеристиками.
Фиг. I
ВИД [1, 7]:
(1)
(2)
Для расчета неоднородных балок удобно выражения площади поперечного сечения и момента инерции представить в виде
р{х) = р0/(х) и І(х) = І0і(х), где Г0 и /,
значения соответствующих величин мерные функции координаты х.
Преобразуем соотношения для момента и силы, а также уравнения (2) после подстановки ний (1) к следующему безразмерному виду:
0 —некоторые амплитудные (константы), а / и і — безраз-
_. _ дЭ
М = і (х)
дх
<3 = — **(*)/Ц-^
перерезывающей в них соотноше-
(3)
д
дх
д
— — I дт о
-/(•*)
д2 ии дР
— <7* (х, і) ;
Здесь
дх
т = — ; л : '"о
(4)
£/о/л> ’
£/о/г02 ’
ЕГо/4
•!й=
к (х) (7
Опуская черту над безразмерными символами, представим систему уравнений (4) для случая непрерывного изменения / и I в следующем виде:
дг т 1 д2 т л , _ч / дда п\ ' п/ л.
* (.-*■> Ч~1
•х.2 (х) дР дх \ дх
д2р дх2
Здесь обозначено
л (*) = —и-/«
ж--см
дт
дх
дт. (х)
т-Цх) -й«^-
(5)
дх
у. (х) дх
В{х)
ді (х)
С(х) = **(х)^- , Р(х,і) =
і (х)
і (х) дх
Я* (х, О
/(■*)
(6)
Из уравнений (4) как частный случай получаются известные уравнения Тимошенко для балок постоянного сечения.
ЧИСЛЕННЫЙ АЛГОРИТМ РЕШЕНИЯ ОСНОВНЫХ УРАВНЕНИЙ НА ОСНОВЕ МЕТОДА ХАРАКТЕРИСТИК
В общем случае для неоднородных балок получить аналитическое решение уравнений (5) невозможно, поэтому возникает задача разработки алгоритма численного решения. Как известно [1], уравнения Тимошенко относятся к гиперболическому типу, поэтому для построения алгоритма численного решения уравнений можно воспользоваться методом характеристик.
Естественно считать неизвестные величины ту и р непрерывными функциями координат л, і, производные которых могут
претерпевать разрывы. Для тех областей плоскости х1, в коточастные производные от ад и р непрерывны, можно
рых первые написать
сі
дт
дх
дт
Ж
д2 т , =='д^(ІХ +
д - т = дхді
сіх
а -£■
д2 да
_д=р_
дхд(
сН\
<Н:
а
ді
д2$
дхд(
СІХ + ^йі-
(7)
Шесть уравнений (5) и (7) могут быть использованы для определения шести вторых производных, если сами функции да, р и их первые производные заданы вдоль некоторой кривой. Однако при задании функций да, 8, дт/дх, дт/дЬ, дфідх, д$/ді вдоль некоторых направлений решения уравнений (5) и (7) будут неопределенными. Такие направления называются характеристическими, а линии, идущие вдоль таких направлений, —характеристиками. На характеристиках вторые производные могут претерпевать разрывы.
На характеристиках определитель Л системы уравнений (5) и (7) равен нулю:
і
ііх2 = 0.
(8)
Приравнивая нулю выражение в первых скобках (8), получим уравнения прямых линий на плоскости ■ с1х
сН
1,
(9)
которые назовем с{ - и сГ -характеристиками. Вдоль этих направлений можно получить характеристические соотношения, если приравнять нулю определитель матрицы коэффициентов системы (5), (7) со столбцом коэффициентов при неизвестном д2$/дх'2, замененным на столбец правых частей. Раскрывая получающийся таким путем определитель и учитывая (9), найдем
дт
Л
д(
- 3
дх) 1
Я1йх — 0,
(10)
Верхние знаки соответствуют характеристике с?, нижние—характеристике сГ-
Приравнивая нулю выражение во вторых скобках (8), получим уравнения кривых йх/сН = + * (х), которые назовем С2- и с-Г-характе-ристиками.
Определяя далее д^т/дх2 из системы (5), (7), получаем, что вдоль характеристик с} и с7 справедливы соответственно соотношения
дт
. . д~х
где
г. <^3 я / ..ч ( дт ^ 1
/<2 :
іІ
дт
~дГ
дх
■ ХСІ
Л(х)
СІХ = 0,
(П)
дх
Р{х, і).
Равенство нулю знаменателей и числителей в решениях для д2ги!дхд(, д2,№1д£2, д2р/сЬсдЛ д2р/д/2 приводит к уравнениям, идентичным (10) и (11). Так как рассматриваются только непрерывные т и {$, то соотношения
, дю , , ди>
ат — —г- ах -4—а1
дх 1 д/
ар
(12)
справедливы вдоль любого направления плоскости хі. В тех областях плоскости хі, в которых первые производные ОТ 10 и р непрерывны, соотношения (Ю)—(12) представляют собой систему уравнений и могут быть использованы для нахождения шести неизвестных її), р, дъи/дх, дъо/дЬ, д$/дх, др/'ді, если заданы соответствующие граничные и начальные условия.
Для численного решения этой системы построим на плоскости х/ сетку характеристик. Так как мы имеем дело с семейством четырех однопараметрических семейств линий—характеристик положительного и отрицательного наклона, то, вообще говоря, можно сетки двух семейств характеристик строить независимо одну от другой. Тогда значения искомых функций в узлах сетки определяются с помощью сложных интерполяций, что приводит, в конечном счете, к потере точности результатов. Целесообразно поэтому на плоскости хі сначала построить сетку характеристик сіх = + сІЇ с шагом /г, а затем характеристики <іх = + ч-сИ вписать в основную сетку,, полагая при этом в пределах одной ячейки х постоянной*. Тогда интерполяционные формулы будут едины для всех ячеек.
Расчетные формулы с учетом пересчета найдем для произвольной ячейки N (см. фиг. 1,а). Отрезки прямых и би-
характеристики первого семейства, N8 и — характеристики второго семейства. Через точки N и Б проведем пунктирными линиями характеристики третьего (ЛуУ и 5£>) и четвертого (Сб и N8) семейств.
Введем обозначения дт)дх = р, дт/д( = д, д$/дх = и, д$/д1=:У и перепишем в следующем виде характеристики системы и соответствующие им соотношения вдоль характеристик:
= 1, йю — йи = [С (х) (р — р) + В (х) и\ сИ , йф = (V и) сИ\
йх (іі
~ = — 1, йюсій = \С (х) (р— р) -)- В (х) и\ сІЇ, с1$ = (г>~и)сИ-,
йх
СІІ
-х, сід — мір
сітю = (д + *р) Ш]
Лх
Ж
•—х, сід + -/.йр = — -х2
и — А (х) (р — Р)-------— Р{х, і)
СІІ,
(13)
сі,ш = (д — тир) сіі.
і
* Так как на участках непрерывного изменения сечения балки величина к, а следовательно, и % мало изменяются, то принятое допущение не приводит к существенным погрешностям.
Дифференциальные соотношения (13) вдоль соответствующих характеристик следует заменить приближенными разностными выражениями. Из полученной таким образом системы алгебраических уравнений определяются неизвестные а, р, V, (3, q и w в точке N через значения этих функций (с использованием линейной интерполяции) в узлах W, S и Е ячейки основной сетки характеристик. Отметим при этом, что в процессе составления разностных выражений, заменяющих дифференциальные (13), значения величин, стоящие в скобках в правых частях соотношений (13), на соответствующих отрезках характеристик заменяются среднеарифметическими из крайних значений на этих отрезках. Такая процедура получения расчетных формул эквивалентна пересчету, который, как известно [10], повышает точность на один порядок. Опуская громоздкие алгебраические выкладки, выпишем окончательные расчетные формулы для произвольной внутренней (полевой) точки области переменных (х, t)\
uN = us + vE-vw-Р/у— Р$ +' ЯЕ Qw '
1 \ СсЛ2 '
vtt =---C7h* Г*-4 (vs + w + ve) +
»+—'
Co h. fr
4—2—(Pn 4* Ps — ^Ps) 4“11E — wu/+ -у \Bs{uN 4- us) +
4~ СE (p£ (3£) -|- Cv,(P\)y Рк7) + Bp uE 4- В wuw ]},
Р.У + ~2~ 0VN + Vw -j- Vs 4- VF)\
Qn = Qs + x2 \Pe ~~ Pf ~ (p W— + IT (^Л- + Ps +
4“ P\v ~r Pg ~Ь [qg qw ~ СЭл/ Ps)l "4 ^ e (Pe Pc) 4“ ^w^P'X' Pup)>
™,v 4“ (Яц 4" Яs 4~ Я w ~1’ Яe)'
Выпишем формулы для точек границы области xt. Рассмотрим вначале полуячейку WNE, примыкающую к границе t = 0 (см. фиг. 1,6). Пусть в начальный момент балка покоилась и напряжения в ней равнялись нулю. В этом случае
то = -ЁЕ-=-^.==8 = -^- = -^- = 0 (14)
дх dt ■ дх dt { )
Заменив соотношения (13) вдоль соответствующих характеристик разностными выражениями и решив с учетом (14) получающуюся таким путем алгебраическую систему уравнений, найдем:
un=Pn'=vn= P/v = °;
Яы~ 2 [ N ' 2
wn==~Y(Jn-
Выпишем формулы, справедливые на левой и правой границах области переменных xt. На фиг. 1, в, г изображены полуячейки SEN и SWN, примыкающие к соответствующим границам. Рассмотрим случай шарнирного опирания балки по обоим концам и
случай, когда один конец балки оперт, другой свободен. На шарнирных опорах смещение, скорость и изгибающий момент равны нулю, следовательно, ,
Представляя в разностном виде соотношения (13) вдоль соответствующих характеристик и решая с учетом (15) получающуюся таким путем алгебраическую систему уравнений, найдем для точек левой границы
Выпишем, наконец, формулы для случая, когда правый конец балки (х — \) свободен. На свободном конце момент и перерезывающая сила равны нулю, следовательно,
Используя равенство (16) и соотношения (13), записанные в разностном виде вдоль соответствующих характеристик, найдем:
Аналогичным образом могут быть получены формулы и для других видов граничных условий.
Выведенные формулы позволяют шаг за шагом находить в каждой точке /V сетки характеристик на плоскости переменных х1 шесть неизвестных функций и, р, V, [3, д, то. При этом безразмерные момент и перерезывающая сила в точке N определяются по формулам (3).
На основе изложенного для ЭЦВМ БЭСМ-6 составлена программа расчета деформаций и напряжений в неоднородных балках
= Я3 = = WN^ Яы = чы= О-
(15)
С с Л2 ч
+ С5 /г яЕ — Р5)--------— (ъЕ + ъ3) | ;
Рл/ = ^5 + \ + ^5 + 2 юе)-
Аналогичные формулы для точек правой границы:
Р N = Рэ 2 Я]хп
Vд, V5 4" - | 2 иц, 4- И. [CW (Рур Р^) 4"
Со Л2 )
-|- с5 и. (р— с/^ р^) 2 “I- ^5) | >
= ^5 ~2~
иЫ=и5'~Рм — Р5 = °-
(16)
+ хЗ н ^ (2 - X) (Р„ + Р5);
2____^
+ -2" ^ + 2 <?», + —(д„ + я5) •
под действием динамических нагрузок. Исходные данные задачи, такие как нагрузка, жесткостные и массовые характеристики, задаются в виде таблиц для определенного числа сечений. Промежуточные значения соответствующих величин вычисляются с помощью интерполяции. Шаг интегрирования выбирается по усмотрению программиста, что дает возможность путем изменения шага и повторного расчета оценивать точность численного решения. Такая оценка может быть выполнена по формулам, указанным в работе [10].
Обозначим через Ал удлинение полуволны (если >. — удлинение всей балки, п— число полуволн, на которые разбивается балка при
поперечных колебаниях с п полуволнами, то 'кп=~^)- Шаг интег~
рирования должен быть достаточно малым по сравнению с длиной полуволны наивысшей из гармоник, вносящей в суммарный эффект волн существенный вклад. Можно показать, что для правильного учета гармоники с длиной полуволны, равной Хп, заданный шаг 'п должен быть меньше Хл примерно на один порядок. Так как минимальное значение ~кп равно примерно высоте балки, то можно ограничить снизу начальный шаг величиной порядка десятых долей г.
ПРИМЕРЫ ЧИСЛЕННОГО РАСЧЕТА
1. Рассмотрим вынужденные колебания шарнирно опертого усеченного конуса с углом раствора 2 а под действием распределенной нагрузки. Начало координат принято в центре большего основания, и ось л: направлена по оси конуса. Площадь и момент инерции (безразмерные) поперечного сечения с координатой к выражаются при этом следующими формулами:
f(x) = (1 — лМ^а)2, г(х) = (1 — лПа*)4.
Далее, используя (8), получим
л<*) = -п^тв>
В расчете принято Х = 20, а = 0,04, что соответствует значению -^- = 0,4 (/— высота усеченного конуса, /. — высота полного
конуса); нагрузка по времени нарастает линейно до значения £ = 20, далее мгновенно снимается (время действия нагрузки составляет примерно 0,06 периода первого тона колебаний). Амплитуда нагрузки (безразмерное значение) принята равной единице, параметр х — равным 0,7.
Расчет выполнен для трех различных значений шага /г (0,2; 0,1; 0,05) и по полученным данным с помощью метода Рунге [10] сделана экстраполяция на нулевой шаг. Результаты расчета представлены в виде графиков.
На фиг. 2 для различных моментов времени I показаны кривые прогибов, Для максимальных амплитуд формы кривых прогибов весьма близки к формам собственных колебаний усеченного конуса по первому тону, а разница в значениях времени, отвечающая максимальным отклонениям балки в противоположные стороны, практически совпадает с полупериодом собственных колебаний по первому тону. Максимальная погрешность расчета за период не превышает 1%.
На фиг. 3 представлены эпюры изгибающих моментов для различных моментов времени. Из поведения кривых следует, что
и при действии равномерно распределенной нагрузки в балке возбуждаются высокие тона колебаний со значительной амплитудой, если время изменения нагрузки достаточно мало.
2. Рассмотрим вынужденные колебания опертой на одном и свободной на другом конце однородной балки с удлинением Х = 30 при действии распределенной динамической нагрузки. Нагрузка по времени задавалась по форме кривой, изображенной на фиг. 4 (в верхнем левом углу). Время нарастания нагрузки до максимума принималось равным 0,3 от общей продолжительности ее действия. Продолжительность действия нагрузки принималась в расчете равной 45, что составляет примерно 0,12 периода первого тона собственных колебаний.
На фиг. 4 изображены эпюры изгибающих моментов для различных значений времени. Для рассматриваемого примера по теории Бернулли — Эйлера существует точное решение в виде ряда по формам собственных колебаний. Это решение в виде графиков изгибающего момента представлено на фиг. 4 (пунктирные линии) для двух характерных значений времени (t = 63 и 99). Существенное количественное различие решений данной задачи по теориям Бернулли — Эйлера и Тимошенко вполне объяснимо на основе изложенного выше.
ЛИТЕРАТУРА
1. Пономарев С. Д., Бидерман В. Л., Л и х о р е в К. К., Макушин В. М., Малинин Н. Н., Ф е о д о с ь е в В. И. Расчеты на прочность в машиностроении. Т. III. М., Машгиз, 1959.
2. Duwez Р. Е., СI а г k D. S. and Bohnenblust N. F. The behavior of long beams under Impact loading. J. of Appl. Mech., 17,
No 1, 1950. (См. также сб. „Механика". М., Изд. иностр. лит., 1950).
3. С л е п я н Л. И. Нестационарные упругие волны. Л., „Судостроение", 1972.
4. Архипов А. С.. Сорокин Е. С. Исследование свободных поперечных колебаний балки как плоской задачи теории упругости.
В сб. „Строительная механика". М., Стройиздат, 1966.
5. Aprahamian R., Evensen D. A. Applications of holography to dynamics: high-frequency vibrations of beams. J. of Appl. Mech., vol. 37, series E, No 2, 1970. (См. также „Прикладная механика", 1970, № 2).
6. Б а л а б у x Л. И., Колесников К. С., Я а р у б и н В. С.,
А л ф у т о в Н. А., У с ю к и н В. И., Ч и ж о в В. Ф. Основы строительной механики ракет. М., „Высшая школа", 1969.
7. Л я х о в е н к о И. А. Поперечные волны в составных балках. Труды ЦАГИ, № 1085, 1967.
8. МГ-Л dl in R. D. and Deresiewicz H. Timoschenko’s shear coefficient for fluxural vibrations of beams. Proceedings of the Second U. S. National Congress of Applied Mechanics, ASME, 1954.
9. G a r r e I i с k J, М., В e n v e n i s I e J. E. Comparison of beam impact models. AIAA. J., No 4, 1970. (См. также „Ракетная техника и космонавтика", т, 8, № 4, 1970).
10. Жуков А. И. Применение метода характеристик к численному решению одномерных задач газовой динамики. Труды ма-тем. ин-та им. В. А. Стеклова, № 58. М., 19R0.
Рукопись поступила 7/ VIII 1972 г.