О численном решении задач устойчивости прямых стержней с промежуточными упругими опорами Текст научной статьи по специальности «Физика»

УЧЕНЫЕ ЗАПИСКИ ЦАГИ

Том XIX

1988

№ 1

УДК 624.071.3

О ЧИСЛЕННОМ РЕШЕНИИ ЗАДАЧ УСТОЙЧИВОСТИ ПРЯМЫХ СТЕРЖНЕЙ С ПРОМЕЖУТОЧНЫМИ УПРУГИМИ ОПОРАМИ

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

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

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

В известных методах [2, 3, 5, 7] численного решения краевых задач механики для систем обыкновенных дифференциальных уравнений граничные условия заданы в начале и конце отрезка интегрироваия. В настоящей статье рассматривается алгоритм, позволяющий сводить краевую задачу для системы дифференциальных уравнений к системе обыкновенных уравнений для случаев, когда внутри интервала интегрирования в любом конечном числе точек заданы дополнительные граничные условия. Очевидно, при этом решения исходной системы уравнений будут разрывными в точках, где заданы дополнительные граничные условия, а число таких разрывов в каждой точке будет равно числу дополнительных условий. При этом величины скачков искомых функций в точках их разрыва являются неизвестными величинами, которые могут быть найдены при решении краевой задачи.

Такого рода задачи широко распространены в строительной механике, в частности, в настоящей статье рассматривается приложение-алгоритма к задачам устойчивости прямолинейных стержней на упругом основании с промежуточными упругими опорами.

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

Р. А. Мирошник

(1.1)

где X — параметр;

У =

У1

У2

Уп

/=

/л

1 и

\л /

(1.2)

Граничные условия на концах интервала интегрирования имеют

вид:

Здесь

?о[^(0)]=0 при х = 0; «¥>* [К (1)] = 0 при х = 1.

т <п.

(1.3)

(1.4)

Тт ' ^ ^п

Решения системы (1.1) У (х) могут иметь разрывы в р точках внутри интервала интегрирования с координатами Хц (0<Хг<1; 1=1,...,р) благодаря наличию дополнительных граничных условий в каждой из точек хг-

?Л^(*г). АК (л:^] == 0, (1.5)

где ДК(л^) = К+(л:() — У (хг)— векторы скачков разрывных функций У(х) в точках х,; У~(х1), У+(х1)-точке х1 до и после разрыва;

Ъ = '

значения функций У(х) в

уг</г; Д У (.*,•) =

У^(^)—УГ(хг)

(1.6)

У„ (*/) - У« (*/)

Количество граничных условий /', в каждой точке Х{ равно количеству функций У (х), имеющих разрывы в этой точке.

В граничных условиях (1.5) векторы <р(- могут зависеть как от функций У(х1), являющихся непрерывными в точке хь так и от величин скачков Д У(х{) разрывных функций У(х1).

При численном интегрировании системы (1.1) величины скачков ДК(.г,) являются неизвестными параметрами, которые могут быть найдены из рассмотрения граничных условий (1.3), (1.5).

Решение задачи проходит в следующем порядке. Вводится вектор неизвестных параметров А размерности я, состоящий из начальных значений функций У при х = 0 и величин скачков этих же функций в точке Хг'.

Ф!

ф-1 Ф-‘ 1=1 I- <1-8)

Ф,

В зависимости (1.7) размерность каждого из векторов ДУ(х1) (/ = = 1 будет равна /г<п, а в вектор У (0) включаются только не-

известные значения функций ¥ в начале интегрирования (л;=0).

Пусть решение системы (1.1), удовлетворяющее граничным условиям (1.3), (1.5), известно при значении параметра Х=Ао, тогда известен

и вектор .4 =-Д(А)(Х0). Предполагая, что А{к)(к0) является приближенным решением уравнений (1.1) при Х = Х0 + ДХ (ДХ —малая величина) по методу Ньютона с помощью итерационной зависимости находится уточненное значение Л(й+1) (Х0ДХ) в виде

А(к+1) (X + ДХ) = 4(й) (X) — (Л(й)) Ф (А(к)). (1.9)

Здесь XV — матрица Якоби:

д®! дФ,

да, ......... да3

^I- (1л°)

(ЗФ, <ЭФ9

да! ............ да3

Зависимость (1.9) позволяет осуществлять итерационный процесс нахождения вектора А, при условии, что матрица (1.10) не вырождена. Элементы матрицы (1.10), являющиеся частными производными функций Ф по неизвестным компонентам вектора .<4, находятся непосредственным дифференцированием для первых т координат вектора Ф и численно для остальных. Для формирования матрицы (1.10) необходимо численно интегрировать систему (1.1) (5 — т+ 1) раз.

Для нахождения решения системы (1.1) при произвольном значении, параметра К найденное значение вектора А (Хо+ЛЛ) предполагается приближенным решением при Х = Х0 + 2Д1, и с помощью (1.9) находится новый вектор А (Хо+2ДА).

Далее, в соответствии с изложенным, по шагам ДХ решение продолжается до произвольного значения К. Выбирая достаточно малые значения ДХ, всегда можно добиться сходимости итерационного процесса (1.9).

В большинстве задач строительной механики роль параметра X играет масштаб нагрузки; полагая Хо = 0, можно начинать решение с тривиального К = 0 или с близкого к нему.

В случае задания граничных условий только в точках л: = 0 и 1 (двухточечная краевая задача), когда уравнения (1.5) отсутствуют, описанный алгоритм полностью совпадает с [3].

Для случая, когда исходная система уравнений (1.1) является линейной, нахождение итераций по формуле (1.9) производится за один

шаг при задании любого начального приближения. При этом для двухточечной задачи объем вычислений как в описанном алгоритме, так и в известных методах [2] примерно одинаков, так как наиболее трудоемкая часть вычислений связана с интегрированием системы (5 — т+ 1) раз.

2. Для иллюстрации Метода рассмотрена распространенная в строительной механике самолета тонкостенная конструкция, показанная на рис. 1. Двухпоясная продольная силовая балка 1 состоит из пояса 2Г связанного с обшивкой 3, свободного пояса 4 и тонкой стенки 5. Балка 1 подкреплена промежуточными опорами 6, а по краям связана с поперечными балками 7. Исследуется устойчивость пояса 4.

Силовые факторы в балке 1 зависят от произвольной распределенной нагрузки q и сосредоточенных силовых факторов Р, при этом на-пряженно-деформированное состояние (НДС) балки предполагается известным. Задача определения НДС может быть решена либо традиционными методами строительной механики, либо для более сложных случаев методом конечных элементов (МКЭ).

В случае наличия в свободном поясе 4 сжимающих напряжений возможна потеря устойчивости рассматриваемой балки 1 путем ее опрокидывания.

Рассмотрим изгиб свободного пояса 4 в плоскости хОу, как балки на упругом винклеровском основании, при этом нахождение коэффициента постели упругого основания не представляет принципиальных трудностей (роль упругого основания здесь играет тонкая стенка 5). Промежуточные опоры 6 могут рассматриваться как упругие сосредоточенные опоры известной жесткости, а поперечные балки 7 — как жесткие опоры. Знание НДС позволяет считать, что нормальная сила в поясе 4 известна.

Различные подходы к расчету балок на упругом основании рассмотрены в [1, 4, 6], при этом в общем случае переменной нормальной силы и жесткости балки наиболее распространены энергетические методы. Эти методы эффективны в случае предсказания хотя бы приближенно формы потери устойчивости, что не всегда возможно.

В настоящей статье рассмотрен один из возможных алгоритмов решения задач устойчивости, связанный с численным интегрированием дифференциальных уравнений равновесия.

Линеаризованные дифференциальные уравнения изгиба балки в плоскости хОу [1] записываются в виде:

-гг = -ет + « ж = Ч+л'(»+Х>.<2-»

где х — координата элемента балки; у — прогиб балки; Ф—угол поворота; Ы, С?, М — нормальная, перерезывающая сила и изгибающий момент в балке; & — коэффициент упругости основания (постели); ЕІ — из-гибная жесткость балки; у', у"0 — производные от известного начального прогиба балки Уо(х).

Известная нормальная сила в балке представлена в виде

N — Мпах А^1 (л), (2.2)

где Л^тах — максимальное значение нормальной силы; N^(x)—известная функция, описывающая закон изменения по длине нормальной силы.

Уравнения (2.1) в безразмерном виде имеют вид:

аГуО

сіх0

с1<2 о

(1х<>

■ «У0;

сіх0

йШ°

= до + д'о(& + ^'),

(2.3)

EJ

Л/Гп ЛІІ АГ0 ^шах 1 „ £

ЛГ° = ; А'шах = —; а = -£у-; I - длина балки.

Введение начального прогиба балки у0 в уравнения равновесия

(2.3) позволяет получить неоднородную систему линейных дифференциальных уравнений, для которой численно может быть найдено решение краевой задачи, отличное от тривиального при любом значении Мтах в (2.2).При задании малых значений у°< 1 характер решения дифференциальных уравнений (2.3) показан на рис. 2. Как видно из рисунка, максимальные прогибы балки _у^ах неограниченно растут при приближении Л/тах к А^кр, где Л^р — критическое значение максимального нормального усилия в балке при заданном законе его изменения и заданных граничных условиях. При этом при ЛГ°ах > Л^р происходит смена знака у^ах.

Указанный характер поведения решения совпадает с результатами аналитического и приближенного исследования устойчивости с начальными несовершенствами для различных частных задач, приведенных в [1, 4, 6].

Как показывают расчеты, при малых значениях у°^1 характер кривой, приведенной на рис. 2, не зависит от закона изменения по длине _У° = У°(*), а разрыв этой кривой происходит при -^°ах = А^°р. Это позволяет находить численное значение уу£р в следующей последова-

тельности. Задается произвольный закон изменения начального прогиба y°(x)<3^1. При задании начальных прогибов намного меньше, чем необходимая точность удовлетворения граничных условий балки, функция у“(х) может и не удовлетворять заданным граничным условиям.

Для значения _/V°ax==0 с помощью изложенного в п. 1 алгоритма решается краевая задача для заданных граничных условий и вычисляется значение наибольшего прогиба у° ах. Исходя из необходимой точности определения дг° , выбирается шаг AN движения по 7V“ax и на каждом шаге после решения краевой задачи вычисляется у^ах. Смена знака .V0 на очередном шаге указывает на то, что в этом интервале находится /V® .

кр

Очевидно, целесообразно начальные шаги А N выбирать достаточно большими, а при смене знака ,у°ах вернуться назад, уменьшив AN, исходя из заданной точности нахождения д/°р. Таким образом, значение N°ax , при котором происходит потеря устойчивости балки, находится с любой наперед заданной точностью, а форма потери устойчивости находится в процессе интегрирования уравнений (2.3).

3. В качестве числовых примеров рассмотрено применение изложенного алгоритма к решению задач устойчивости двухопорных балок на упругом основании при различных законах изменения по длине № (£7 = const).

Пусть в плоскости xOz шарнирно закрепленной балки 1 (рис. 1) действует равномерно распределенная нагрузка q (схема 2, рис. 3), промежуточные опоры при этом отсутствуют. В поясе 4, рассматриваемом как балка на упругом основании, продольная сжимающая сила № будет изменяться по закону

iV0=4AC*(*2-x). (3.1)

Схема 2

^ пшпкшш—.—

/777 /777 /777/77? /777

^гпттп^гт°аШ\^^^Ш^Іах

/Яг ос =0 ^ ^ а =0 а =0

сппп "" - тип

а =5000

«=30000

<Х=1000рЯ’'* * = 1000

Л в л

ілл /777 /777л. 7ПППП "17 /777 оллл/

а. =30000

<х=гоооо

ж

йш ДШ7

Рис. З

лш?

С<

Для балки, шарнирно-опертой в плоскости хОг/, граничные

(1.3) имеют вид:

У>(0) = 0; М°(0) = 0; у°(1) = 0; Л1°(1) = 0,

а дополнительные граничные условия (1.5) отсутствуют.

Вектор неизвестных параметров Л имеет вид:

-СМ

0(0)

ста

а граничные условия (3.2) при х=\

?1=^0(<*1. а2) = 0;

(р2 = /И°(а1, а2) =0.

Вектор Ф (1.8) и матрица Якоби IV (1.10) имеют вид

Ф =

Фі

Фо

д?і д?і дії ^

даг да2

условия

(3.2)

(3.3)

(3.4)

(3.5)

Использование зависимости (1.9) позволяет за одну итерацию найти неизвестные значения 0(0) и <2°(0) при заданных значениях Д/“ах.

На кривой 2 рис. 3 приведена зависимость критической нагрузки №кр/№э от безразмерного параметра а (Л^=л;2 — значение эйлеровой критической силы для двухопорной балки без упругого основания, сжатой постоянной по длине силой №). Как показывает анализ результатов расчета, задание начального прогиба у° позволяет избежать трудностей, связанных с плохой обусловленностью задачи при больших значениях а [2]. Это связано с тем, что, когда убывающие решения (2.3) становятся у концов балки намного меньше функций у°0, у°', у°", определяющими для нахождения неизвестных начальных параметров становятся заданные начальные прогибы.

Кривая 1 на рис. 3 позволяет определить критическую силу при нагружении балки в плоскости хОг сосредоточенным изгибающим моментом, когда

Л°=-Л£ах*.

Кривая 3 рис. 3 построена для случая приложения в середине пролета в той же плоскости хОг сосредоточенной силы, при этом принято

ЛР = -2ЛС„*, 0<х<0,5; ) _

№ = 2Л/°ах(х-1), 0,5<х<1. }

Для всех указанных кривых на рис. 3 показаны также формы потери устойчивости при различных значениях а.

Кривая 4 на этом же рисунке построена для случая Л^°=А?шах = = соп51;, при этом для значений а>20 000 эта кривая хорошо согласуется с приведенным в [1, 4] аналитическим решением.

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

4. Рассмотрена балка на упругом основании с промежуточными упругими опорами. Концы балки приняты закрепленными шарнирно в обеих плоскостях хОг и хОу, а в плоскости хОг балка предполагается нагруженной равномерно распределенной нагрузкой. Сжимающая сила изменяется по зависимости (3.1). Пусть указанная балка имеет три промежуточные упругие опоры, расположенные равномерно в точках с

координатами л:? = 0,25; х\ = 0,5; х1 = 0,75. На концах балки граничные условия (3.2) и (3.4) сохраняются, а граничные условия (1.5) записываются в виде

?„ = Д (0,25)+ 5° У0 (0,25); |

Та1 = Д <2° (0,5) + 6°^° (0,5); - (4.1)

?«1 = А0° (0,75)+ 5°У (0,75), )

где Д<3°( ) —скачки перерезывающей силы 0й в точках, где име-

ш

ются промежуточные опоры; = —жесткость промежуточ-

ных опор, предполагаемая одинаковой.

Записываются векторы А и Ф

0(0)

0°(0)

ДО0 (0,25) А<2° (0,5)

А <2° (0,75)' у*(аи а2, а3, а4, а6у М°(аи а2, а3, а4, а5) = 1 а2) аз)

^2» ®з>

<р31(а,, а2, я3, а4,

(4.2)

Матрица Якоби IV имеет вид:

/

дФ) да!

дФ2

даг

д% да!

6>Ф4

да,

дФь

дал

аФ,

да2

ЭФ2

<ЭФ3

да2

дФ1

да2

дФ5

дао

_дФ1

д«з

<м>2

<?в3

1

дФ,

даг

<?Ф5

дао

(4.3)

Как и выше, один шаг итерации с использованием (1.9) позволяет найти значения вектора неизвестных параметров А. Аналогично изложенному составляются выражения для решения краевой задачи при любом количестве опор.

При наличии упругой опоры, препятствующей повороту в точке х9, граничное условие для скачка момента ДУИ°(л:9) имеет вид:

ДЛ4°(л°) + х° 9 (хУ} = 0,

(4.4)

где %°— безразмерная крутильная жесткость промежуточной опоры.

На рис. 4 приведены кривые 2, 3, 4 для нахождения критической силы при наличии соответственно одной, двух и трех промежуточных опор при описанных выше условиях. Кривые построены для абсолютно жестких опор (|°= 106).

Кривая 1, взятая с рис. 3, построена при отсутствии промежуточных опор, здесь же показаны формы потери устойчивости балок при различных значениях а.

Из рисунка видно, что при а>35000 все кривые практически совпадают, т. е. промежуточные опоры не оказывают влияния на критическую силу.

Рис. 4

Рис. 5

J______L

(x =1000; 5 =100

Л j

/777 7777 7777

a =10000 -,1=100

/Д7 ^77^77

ос=10000 ;£=/Ш a. =20000 ,(,=1000

л—JT JL% ^

дат ттгт

оІЧ0000 ]i=5000

7777 '/// 7777

10000

20000

30000

«

На рис. 5 для описанной выше балки построены зависимости для нахождения критической силы (кривые 1, 2, 3) при различных значениях жесткости промежуточных опор, соответственно равной 100, 1000, 5000, а также показаны формы потери устойчивости балки. Кривая 5, взятая с рис. 4, соответствует трем абсолютно жестким опорам, а кривая 1, взятая с рис. 3, соответствует отсутствию промежуточных опор. Как показывает анализ приведенных кривых, при большом количестве промежуточных упругих опор равной жесткости (начиная с трех опор) при определении критической силы жесткость промежуточных опор может быть распределена по всей длине балки.

Определяется

апР = а+Л°, (4.5)

где / — количество упругих опор.

Для определения критической силы в многоопорной балке можно пользоваться кривыми 4 и 5. Приняв а = аПр, по кривой 4 находится Л^р, соответствующее потере устойчивости с распределенными по длине упругими параметрами, а по кривой 5, построенной для заданного числа опор, для известного значения а находится соответствующее потере устойчивости с абсолютно жесткими опорами. Окончательное значение критической силы будет меньшее из двух найденных.

В заключение отметим, что изложенный алгоритм позволяет определять критическую силу в балках при произвольном законе изменения нормальной силы по длине, при произвольных граничных условиях, а также при любом количестве упругих либо жестких промежуточных опор. Последовательность расчета не меняется и для случаев, когда промежуточные опоры поворачиваются при приложении внешних нагрузок (см. рис. 1) по известному закону (этот закон находится при определении НДС).

Уравнение (4.1) при этом имеет вид:

^.= ДО»(^) + |О[уО{хо)_у0(л0)]==0> (4.6)

где УрО*”) =У°0(х°, Щ— известный из НДС балки закон смещения упругих опор.

ЛИТЕРАТУРА

1. А л футов Н. А. Основы расчета на устойчивость упругих систем.— М: Машиностроение, 1978.

2. Бидерман В. Л. Прикладная теория механических колебаний.—

М: Машиностроение, 1972.

3. Валишвили Н. В. Расчет оболочек вращения на ЭЦВМ. —

М: Машиностроение 1976.

4. В о л ь м и р А. С. Устойчивость деформируемых систем. — М.: Наука, 1967.

5. Давиденко Д. В. О приближенном решении систем нелинейных уравнений. — УМЖ, т. 5, № 2, 1953.

6. Тимошенко С. П. Устойчивость упругих систем. — М.: Наука, 1946.

7. Шаманский В. Е. Методы решения нелинейных краевых задач на ЭЦВМ. — Киев.: «Наукова думка», 1966.

Рукопись поступила 20/У 1986 г.