Analysis of shear flexible beam with both clamped ends by Ritz method

全文:

В настоящее время для производства авиационной и ракетно-космической техники широко применяются композиционные материалы. Композиты обладают высокой удельной прочностью и жесткостью, поэтому позволяют изготавливать конструкции с высокими показателями весового совершенства. Для композитов характерна низкая сдвиговая жесткость по отношению к трансверсальным напряжениям. Учет этой особенности приводит к усложнению математической модели за счет повышения порядка разрешающих уравнений из-за введения в рассмотрение углов трансверсального сдвига. Тотальное господство численных методов, и прежде всего метода конечных элементов (МКЭ), в расчетной и инженерной практике не уменьшает интереса исследователей к приближенным методам, которые в ряде случаев в силу более простого алгоритма реализации позволяют получить приемлемые решения. В настоящей работе рассматривается расчетная модель балки, базирующаяся на математической теории пластин Рейсснера-Мидлина. В расчете, выполняемом с помощью метода Ритца, учитывается транс-версальный сдвиг. Данный подход требует вариационной постановки задачи, при которой в качестве исходного функционала фигурирует выражение потенциальной энергии деформации конструкции. Для сдвиговой модели балки это выражение имеет вид E = ~f f N-+M-- + Qyl dx-fqwdx, (1) 2^ ^ dx dx ) 0 где 0 - угол наклона сечения, измеряемый в плоскости изгиба; u - перемещения точек начальной плоскости вдоль продольной оси балки; w - прогибы точек начальной плоскости; у - угол сдвига или осреднен-ная по высоте сечения деформация трансверсального сдвига; q - погонная балочная нагрузка. Угол наклона сечения 0 складывается из угла сдвига у и угла изгиба dw/dx [1]: I = y dw dx (2) Для внутренних силовых факторов справедливы следующие физические соотношения: _du ~d 0 ~du ^d 9 ^ N = B- + C-, M = C- + D-, Q = Ky , (3) dx dx dx dx где В - жесткость балки при растяжении-сжатии вдоль продольной оси; С - смешанная жесткость; D изгибная жесткость; К - жесткость балки при трансверсальном сдвиге. 54 Вестник СибГАУ. № 1(53). 2014 Рассмотрим балку, закрепленную по обоим концам и нагруженную равномерно распределенной балочной нагрузкой. В случае симметрично армированного слоистого пакета смешанные жесткости C будут равны нулю. Принимая это обстоятельство во внимание, подставим физические соотношения (3) в функционал (1). В результате получим ( < . -2 / . л\2 Л B 2 f + D d y dx d 2 w ^ dx2 + K y2 dx - j qw dx . (4) . | 2nx y = F sin Iu = U sin I-I. (6) X = {U W F} R = {0 r 0}T где j |(sinh (kl) - sin (kl))osh (kx) - cos (kx))) ^ J [-(cosh (kl) - cos (kl) )(sinh (kx) - sin (kx))J dw dx У1 Щ w2 dw dx У 2 u2 ние размером 1 x 1 см. Модуль упругости материала E = 180 ГПа, модуль сдвига G = 5 ГПа. Коэффициенты жесткости вычисляются по формулам B =- Ebh Ebh Граничные условия на защемленных краях выражаются следующими соотношениями: при x = 0 и x = l: w = 0, = 0, y = 0 . (5) dx С учетом данных граничных условий примем следующие распределения для кинематических параметров: {(sinh (kl) - sin (kl))(cosh (kx) - cos (kx)) - | > , -(cosh (kl) - cos (kl))(sinh (kx) - sin (kx)) 2nx' l У 'к ~ Реализация процедуры метода Ритца сводит вариационную задачу к системе алгебраических уравнений вида AX = R, (7) где A - матрица третьего порядка, компоненты которой вычисляются программно; X - вектор неизвестных, состоящий из амплитудных значений в разложениях (6); R - вектор нагрузки: D =--, K = Gbh . (8) 1 - д2 12(1 - д2) Зададим погонную балочную нагрузку q = 100 Н/м и выполним статический расчет двумя вышеописанными методами. Результаты решения задачи с помощью МКЭ и в балочных функциях (методом Ритца) имеют хорошее совпадение по прогибам. Расхождение по центральному прогибу не превышает 1 %. Выполним серию расчетов для той же самой композитной балки, увеличивая высоту ее сечения (при сохранении неизменной ширины 1 см). Результаты расчетов представим таблично (табл. 1). Прогибы, вычисленные с помощью метода Ритца и МКЭ Таблица 1 Размеры сечения, b x h , см x см Центральный прогиб, мм Разница, % Нагрузка, q, Н/м МКЭ Метод Ритца 1 х 1 1,605 1,619 0,87 100 1 х 2 2,098 2,099 0,05 1000 1 х 3 1,336 1,318 1,41 2000 1 х 4 1,235 1,198 2,94 4000 1 х 5 1,406 1,340 4,73 8000 r = q| Решая систему уравнений (7), определим компоненты вектора X, а следовательно, и распределения кинематических параметров (6). Решения, полученные с помощью приближенного метода Ритца, будем сравнивать с результатами численных расчетов, выполняемых методом конечных элементов. Основными узловыми неизвестными являются: прогиб w, изгибной угол наклона сечения, угол трансверсального сдвига и продольное перемещение u. Таким образом, вектор узловых параметров в двухузловом элементе имеет вид Как видно по этим данным, результаты решений двумя методами имеют хорошее совпадение по центральным прогибам. Тем не менее нельзя не заметить, что по мере увеличения высоты сечения балки имеет место тенденция к увеличению разницы прогибов, вычисленных разными методами. Этому обстоятельству можно дать следующее объяснение. При увеличении толщины (высоты сечения) балки возрастает роль сдвигового фактора. В решении, полученном методом Ритца, также учитывается трансверсальный сдвиг, но для точного представления функции прогибов в балке большой толщины одной базисной балочной функции оказывается недостаточно, что приводит к некоторому расхождению с более точным (МКЭ) вариантом расчета. Рассмотрим теперь задачу об изгибе трехслойной балки. Будем полагать, что жесткость несущих слоев существенно выше жесткости промежуточного слоя. В этом случае можно допустить, что параметры жесткости балки B и D обеспечены несущими слоями: 2E„ bt „ E„b(H3 - h3) „ „, / 2t . h (9) B = ■ D = V“/2 Разрешающие уравнения теории МКЭ, соответствующие нашей задаче, представлены в [2], а вывод основополагающих матриц - в [3]. Предположим, что балка, изготовленная из однонаправленного углепластика, имеет квадратное сече- 1 -цн 12(1 -цн) ' / Gн Gз где H - полная толщина пакета; h - толщина слоя заполнителя; t - толщина каждого из несущих слоев; b -ширина сечения; EK - приведенный модуль упругости материала несущих слоев; GK и Gs - модули сдвига материалов несущих слоев и заполнителя соответственно. Выполним с помощью методов Ритца и конечных элементов расчет трехслойной балки, защемленной по 55 Математика, механика, информатика обоим краям и нагруженной равномерным погонным балочным усилием q. Несущие слои толщиной в 1 мм -композитные (Ен = 180 ГПа, дн = 0,3). Материал заполнителя - изотропный с модулем упругости Ез = 5 ГПа и коэффициентом Пуассона цз = 0,3. Модули сдвига материалов несущих слоев и заполнителя определяются по классической формуле G = Е/(1+д). Результаты расчетов для моделей с различной толщиной слоя заполнителя представлены в табл. 2. Они свидетельствуют о том, что приближенное решение методом Ритца с использованием в качестве базисной только одной балочной функции дает достаточно точные результаты по центральному прогибу для модели с относительно жестким слоем заполнителя. Причем на всем интервале изменения толщины слоя заполнителя (от 1 до 5 см) расхождение приближенного решения (методом Ритца) и численного (МКЭ) не превышает 1 %. Для расчета с помощью МКЭ задействована специально разработанная авторская программа [4]. Таблица 2 Прогибы трехслойной балки, вычисленные с помощью метода Ритца и МКЭ Толщина слоя заполнителя h , см Центральный прогиб, мм Разница, % Нагрузка, Н/м МКЭ Метод Ритца 1 2,215 2,233 0,81 100 2 3,110 3,127 0,55 500 3 2,928 2,925 0,10 1000 4 3,424 3,406 0,53 2000 5 4,525 4,483 0,93 4000 Ез, ГПа Центральный прогиб, мм Разница, % МКЭ Метод Ритца 5 4,119 4,120 0,007 1 5,003 4,826 3,54 0,5 6,109 5,696 6,76 0,1 14,983 12,143 18,96 0,05 26,150 19,107 26,93 представления прогибов одной балочной функцией. Для получения более точного решения методом Ритца необходимо увеличить число базисных балочных функций в представлении прогиба. Однако при этом значительно возрастают вычислительные затраты и исчезает то преимущество приближенного решения (методом Ритца) перед численным (МКЭ), которое имело место при одночленной аппроксимации функции прогибов, а именно, простота реализации. Рассмотрим теперь задачу об устойчивости защемленной балки при осевом сжатии. Функционал устойчивости для сдвиговой модели балки имеет вид l Г / , \2 ( - -о \2 d у dx d 2 w Л dx2 Еще один численный эксперимент проведен с целью определения точности приближенного решения методом Ритца при различных значениях жесткости заполнителя. Для этого проведена серия расчетов двумя методами (Ритца и МКЭ) той же самой трехслойной балки с композитными несущими слоями толщиной в 1 мм и слоем заполнителя толщиной 2,5 см. Анализ результатов (табл. 3) позволяет сделать вывод о том, что при малой жесткости слоя заполнителя приближенное решение методом Ритца (с одной базисной функцией для прогиба) дает значение центрального прогиба, существенно отличающееся от прогиба, рассчитанного с помощью МКЭ. Таблица 3 Прогибы трехслойной балки при различной жесткости слоя заполнителя где N0 - усилие осевого сжатия. В функционале устойчивости отсутствует продольное перемещение и, поэтому из распределений (6) здесь актуальны выражения для прогиба w и угла трансверсального сдвига у. Реализация процедуры метода Ритца сводит вариационную задачу к обобщенной задаче на собственные значения AX + N0BX = 0, (10) где A и B - матрицы второго порядка, компоненты которых вычисляются программно; X - вектор неизвестных, состоящий из амплитудных значений в разложениях (6): X = {W F}T . Решая систему уравнений (10), определим критическое значение усилия N0. Форма потери устойчивости задается балочной функцией, фигурирующей в представлении прогиба (6). Полученные с помощью приближенного метода Ритца критические величины будем сравнивать с результатами численных расчетов, выполняемых методом конечных элементов. Основными узловыми неизвестными являются: прогиб w, изгибной угол наклона сечения и угол трансверсального сдвига. Вектор узловых параметров элемента имеет вид s- -Ь I.f а dw У1 w2 I -J У 2 dx Неточность приближенного решения методом Ритца в данном случае обусловлена неадекватностью Вывод системы уравнений устойчивости произвольного балочного элемента выполняется программно [3], а соответствующие данной задаче матрицы жесткости представлены в [5]. Предположим, что балка, изготовленная из однонаправленного углепластика, имеет квадратное сечение размером 1 х 1 см. Модуль упругости материала Е = 180 ГПа, модуль сдвига G = 5 ГПа. Коэффициенты изгибной и сдвиговой жесткости вычисляются по формулам (8). Критическое усилие сжатия, определенное при решении системы (10), соответствующей методу Ритца, оказалось равным 6623,83 Н. Конечноэлементный расчет для 50-узловой модели дал величину в 6444,44 Н. Относительная разница этих значений составляет 2,78 %, что свидетельствует о высокой Т 56 Вестник СибГАУ. № 1(53). 2014 точности приближенного решения задачи устойчивости методом Ритца. В табл. 4 представлены критические нагрузки, определенные в серии расчетов композитной защемленной балки при варьировании высотой сечения. Здесь кроме результатов численного (МКЭ) и приближенного (метод Ритца) решений приводится аналитическое решение, вычисленное по формуле N кр = 4П Ци не учитывающее трансверсальный сдвиг. Значения критических усилий для сплошной композитной балки Таблица 4 h, мм МКЭ, Н Метод Ритца, Н Относи тельная разница, % Аналитическое решение, Н 5 811,43 835,73 2,99 813,43 10 6444,44 6623,83 2,78 6507,43 20 50100,62 51099,08 1,99 52059,45 30 161479,00 162811,12 0,82 175700,65 40 360150,83 358043,35 0,59 416475,61 50 653690,92 640293,63 2,05 813428,93 60 1,040-106 1,004 106 3,46 1,406 106 70 1,509-10° 1,438-106 4,71 2,599-106 80 2,050-106 1,932-106 5,76 3,332-106 90 2,648-106 2,475-106 6,53 4,744-106 для модели с относительно жестким слоем заполнителя. Причем на всем интервале изменения толщины слоя заполнителя (от 1 до 10 см) расхождение приближенного решения (методом Ритца) и численного (МКЭ) не превышает 3 %. Таблица 5 Критические усилия для трехслойной балки с различной толщиной слоя заполнителя h3, мм Метод Ритца, Н МКЭ, Н Отн. разница, % 10 4803,09 4677,02 2,70 20 17147,78 16762,82 2,30 30 36667,77 35986,80 1,89 40 62974,52 62035,97 1,51 50 95708,23 94620,84 1,15 60 134535,24 133474,27 0,79 70 179145,68 178321,76 0,46 80 229251,44 228910,52 0,15 90 284584,22 285045,82 0,16 100 344893,93 346478,19 0,46 Таблица 6 Критические усилия для трехслойной балки с различной жесткостью заполнителя Как видно по представленным результатам (табл. 4), на большей части диапазона изменения высоты сечения балки расхождения критических величин, найденных двумя методами (МКЭ и метод Ритца), не превышают 3 %, что подчеркивает точность метода Ритца при решении задачи устойчивости балок с относительно высоким значением трансверсальной сдвиговой жесткости. Вместе с тем оценивая результаты аналитического решения (четвертый столбец, табл. 4), можно заключить, что уже при данных значениях трансверсальной сдвиговой жесткости нельзя не учитывать в постановке трансверсальную податливость. Рассмотрим теперь задачу об устойчивости трехслойной балки. Коэффициенты изгибной и сдвиговой жесткости вычисляются по формулам (9). Определим с помощью методов Ритца и конечных элементов критическую осевую нагрузку для трехслойной балки, защемленной по обоим краям. Несущие слои толщиной в 1 мм - композитные (Ен = 180 ГПа, цн = 0,3). Материал заполнителя - изотропный с модулем упругости Ез = 5 ГПа и коэффициентом Пуассона ^ = 0,3. Модули сдвига материалов несущих слоев и заполнителя определяются по классической формуле G = Е/(1+д). Ширина сечения балки 10 мм. Результаты расчетов, выполненных с помощью оригинальной программы, для моделей с различной толщиной слоя заполнителя представлены в табл. 5. Они свидетельствуют о том, что приближенное решение методом Ритца с использованием в качестве базисной только одной балочной функции дает достаточно точные результаты по критическому усилию Ез, ГПа Метод Ритца, Н МКЭ, Н Относительная разница, % 5 26035,86 25504,87 2,08 1 22223,03 22439,30 0,96 0,5 18830,16 19508,33 3,48 0,1 8833,13 9540,25 7,41 0,05 5613,39 5822,22 3,59 0,04 4829,16 4872,84 0,90 0,03 3985,85 3831,66 4,02 0,02 3076,50 2684,66 14,60 0,01 2093,06 1414,85 47,94 Еще один численный эксперимент проведен с целью определения точности приближенного решения задачи устойчивости с помощью метода Ритца при различных значениях жесткости заполнителя. Для этого проведена серия расчетов двумя методами (метод Ритца и МКЭ) той же самой трехслойной балки с композитными несущими слоями толщиной в 1 мм и слоем заполнителя толщиной 2,5 см. Анализ результатов (табл. 6) позволяет сделать вывод о том, что при малой жесткости слоя заполнителя приближенное решение методом Ритца (с одной базисной функцией для прогиба) дает значение критической нагрузки, существенно отличающееся от величины, рассчитанной с помощью МКЭ. Неточность приближенного решения методом Ритца в данном случае также обусловлена неадекватностью представления прогибов одной балочной функцией. Для решения задачи о собственных колебаниях получен следующий функционал: 1J 2 J0 dx--ю2 J" и2 + w2 ) + D, dx = 0, 57 Математика, механика, информатика где Bp, Dp - параметры, характеризующие инерционные свойства; ю - круговая частота колебаний. Реализация процедуры метода Ритца при выбранной аппроксимации кинематических параметров (6) сводит вариационную задачу к обобщенной задаче на собственные значения AX + ю2 BX = 0, (11) где A и B - матрицы третьего порядка, компоненты которых вычисляются программно; X - вектор неизвестных, такой же как в системе (7). Решая систему уравнений (11), определим первую частоту собственных колебаний балки. Соответствующая форма колебаний задается балочной функцией, фигурирующей в представлении прогиба (6). Полученные с помощью приближенного метода Ритца собственные значения будем сравнивать с результатами численных расчетов, выполняемых методом конечных элементов. Вывод системы уравнений состояния произвольного балочного элемента представлен в [6]. При этом компоненты матриц жесткости и инерции вычисляются программно [3]. Предположим, что балка, изготовленная из однонаправленного углепластика, имеет квадратное сечение размером 1 х 1 см. Модуль упругости материала Е = 180 ГПа, модуль сдвига G = 5 ГПа. Коэффициенты жесткости вычисляются по формулам (8). Плотность материала р = 1500 кг/м3. Параметры инерции определяются как Bp = pbh, Dp = р bh 12 (12) с точностью до 20 знака после запятой. Отметим, что аналитическое решение (13) не учитывает влияния трансверсального сдвига. Значения критических усилий для сплошной композитной балки Таблица 7 h, мм Метод Ритца, Гц МКЭ, Гц Относительная разница, % Аналитическое решение, Гц 10 117,30 117,10 0,17 118,04 20 230,35 228,83 0,66 236,08 30 335,66 330,93 1,43 354,114 40 430,99 420,97 2,38 472,15 50 515,42 498,17 3,46 590,19 60 589,07 563,14 8,47 708,23 70 652,76 617,09 5,78 826,27 80 707,66 661,70 6,94 944,30 90 755,05 698,49 8,10 1062,34 100 796,14 728,90 9,22 1180,38 Первая частота собственных колебаний, определенная при решении системы (11), соответствующей методу Ритца, оказалась равной 117,3 Гц. Конечноэлементный расчет для 50-узловой модели дал величину в 117,1 Гц. Относительная разница этих значений составляет 0,17 %, что свидетельствует о высокой точности приближенного решения задачи о собственных колебаниях балки методом Ритца. Напомним, что циклическая f и круговая ю частоты колебаний связаны формулой ю = 2 п f. В табл. 7 представлены частоты, определенные в серии расчетов композитной защемленной балки при варьировании высотой сечения. Здесь кроме результатов численного (МКЭ) и приближенного (метод Ритца) решений приводится аналитическое решение, вычисляемое по формуле (13) где i - номер моды; F - площадь сечения балки. Фигурирующий в (13) параметр k - это корень частотного уравнения, которое для защемленной по краям балки имеет вид cosh (kl)cos (kl) -1 = 0 . Для первых четырех мод этот параметр имеет значения 4,730; 7,853; 10,966; 14,137. В расчетной программе в действительности эти корни определены Как видно по представленным результатам (табл. 7), на половине диапазона изменения высоты сечения балки (от 10 до 50 мм) расхождения в значениях первой частоты, найденных двумя методами (МКЭ и метод Ритца), не превышают 3,5 %, что подчеркивает точность метода Ритца при решении задачи о собственных колебаниях балок с относительно высоким значением трансверсальной сдвиговой жесткости. Однако нельзя не отметить тенденцию увеличения расхождения при дальнейшем увеличении высоты сечения. Рассмотрим теперь модальный расчет трехслойной балки. Определим с помощью методов Ритца и конечных элементов критическую осевую нагрузку для трехслойной балки, защемленной по обоим краям. Несущие слои толщиной в 1 мм - композитные (Ен = 180 ГПа, цн = 0,3, рн = 1500 кг/м3). Материал заполнителя - изотропный с модулем упругости Ез = 5 ГПа, коэффициентом Пуассона дз = 0,3 и плотностью рз = 1000 кг/м3. Модули сдвига материалов несущих слоев и заполнителя определяются по классической формуле G = Е/(1+д). Ширина сечения балки 10 мм. Параметры инерции трехслойки определяются следующим образом: BP= 2Рн bt + Рз bh , DP .Рн b (H 3 - h3) , Рз bh 12 12 Результаты расчетов, выполненных с помощью оригинальной программы [4], для моделей с различной толщиной слоя заполнителя представлены в табл. 8. Они свидетельствуют о том, что приближенное решение методом Ритца с использованием в качестве базисной только одной балочной функции дает достаточно точные результаты по первой частоте для модели с относительно жестким слоем заполнителя. Причем на всем интервале изменения толщины слоя заполнителя (от 1 до 10 см) расхождение приближенного решения (методом Ритца) и численного (МКЭ) не превышает 2,5 %. Еще один численный эксперимент проведен с целью определения точности приближенного решения 58 Вестник СибГАУ. № 1(53). 2014 задачи о собственных колебаниях трехслойной балки с помощью метода Ритца при различных значениях жесткости заполнителя. Для этого проведена серия расчетов двумя методами (метод Ритца и МКЭ) той же самой трехслойной балки с композитными несущими слоями толщиной в 1 мм и слоем заполнителя толщиной 2,5 см. Анализ результатов (табл. 9) позволяет сделать вывод о том, что при малой жесткости слоя заполнителя приближенное решение методом Ритца (с одной базисной функцией для прогиба) дает значение первой частоты, существенно отличающееся от величины, рассчитанной с помощью МКЭ. Таблица 8 Значения первой частоты собственных колебаний для трехслойной балки с различной толщиной слоя заполнителя h, мм Метод Ритца, Гц МКЭ, Гц Относительная разница, % 10 107,29 107,05 0,22 30 185,99 184,65 0,73 50 236,94 234,13 1,20 70 275,96 271,50 1,64 90 307,82 301,59 2,07 100 321,81 314,68 2,27 Таблица 9 Значения первой частоты собственных колебаний для трехслойной балки с различной жесткостью слоя заполнителя Ез, ГПа Метод Ритца, Гц МКЭ, Гц Относительная разница, % 5 383,83 381,55 0,60 1 354,68 345,54 2,65 0,5 326,53 312,07 4,63 0,1 223,70 197,93 13,02 0,05 178,34 149,65 19,17 0,01 108,90 71,15 53,06 Таким образом, мы видим, что метод Ритца обеспечивает приемлемые результаты решения задачи о собственных колебаниях податливых при трансвер-сальном сдвиге балок и для сплошных композитных и трехслойных балок с относительно высокими значениями трансверсальной сдвиговой жесткости. На примере задач об изгибе, устойчивости и собственных колебаниях защемленной с двух сторон балки, податливой при трансверсальном сдвиге, исследована точность метода Ритца с удержанием одного члена базисных балочных функций при аппроксимации прогибов. Выявлены границы применимости этого способа решения для расчета сплошных композитных балок и трехслойных балок с податливым заполнителем.

作者简介

Vladimir Nesterov

Siberian State Aerospace University named after academician M. F. Reshetnev

Email: misternester@gmail.com
Candidate of Engineering Sciences, associate professor, associate professor of the Department of Computer modeling

参考

补充文件