УЧЕНЫЕ ЗАПИСКИ Ц А Г И Т о м XV 19 8 4
№ 6
УДК 629.7.015.624.073
ИССЛЕДОВАНИЕ НАПРЯЖЕННО-ДЕФОРМИРОВАННОГО СОСТОЯНИЯ ОБОЛОЧКИ С ОТВЕРСТИЕМ МЕТОДАМИ КОНЕЧНЫХ ЭЛЕМЕНТОВ И ГОЛОГРАФИЧЕСКОЙ ИНТЕРФЕРОМЕТРИИ
Т. К■ Бегеев, В. И. Гришин, В. С. Писарев
Рассматривается расчетно-экспериментальный метод исследования напряженно-деформированного состояния элементов конструкций, основанный на совместном применении методов голографической интерферометрии и конечных элементов. Приводятся результаты определения напряжений на контуре кругового отверстия в тонкостенной цилиндрической оболочке, которые в пределах погрешности метода согласуются с данными известных теоретических и экспериментальных исследований.
При исследовании местной прочности элементов авиационных конструкций, связанном с оценкой концентрации напряжений в окрестностях отверстий, подрезов, соединений, широкое применение получил метод конечных элементов (МКЭ) [1]. Однако в ряде случаев эффективность анализа напряжений с помощью МКЭ существенно ограничивается сложностью расчетных схем в подобных задачах, обусловленной необходимостью учета часто весьма неопределенных статических и кинематических граничных условий вдали от концентраторов. По этой причине значительный интерес вызывает возможность создания расчет-но-экспериментальных методов исследования напряженно-деформированного состояния (НДС). Такие методы могут быть основаны, например, на комбинации поляризационно-оптического метода и различных способов аналогового или численного интегрирования дифференциальных уравнений теории упругости [2, 3].
Одним из перспективных сочетаний представляется совместное использование МКЭ и метода голографического интерферометрии (МГИ) [4], позволяющего бесконтактно измерять трехмерные поля перемещений поверхности реальных объектов сложной формы. В данной работе рассматриваются принципиальные вопросы, связанные с применением подобного подхода для определения НДС в районах концентраторов типа отверстий в тонкостенных элементах конструкций.
1. Напряженно-деформированное состояние в упругом теле V под действием внешних р и массовых / сил должно удовлетворять дифференциальным уравнениям равновесия в области V:
статическим граничным условиям на границе
М<з^=р,
кинематическим граничным условиям на границе 52:
и,=*и0,
0)
(2)
(3)
уравнениям связи между деформациями и перемещениями в области V:
где Вг и £>0 — дифференциальные операторы [5]; М — матрица, элементы которой зависят от ориентации единичного вектора нормали к границе ¿>2 области V; С — матрица закона Гука: и и 8— напряжения и деформации соответственно; ¿/2, 1/0, I/—векторы перемещений.
Подставляя выражение (4) в (5), а полученное соотношение в (1) и (2), можно получить:
Решение системы уравнений равновесия (6) относительно неизвестных перемещений и в области V находится из условия стационарности функционала полной потенциальной энергии конструкции [6]:
где Т — индекс транспонирования.
Минимизация функционала (7) ведется по перемещениям в узлах конечноэлементной сетки. Точность решения системы уравнений (6) для определенной аппроксимации поля перемещений в элементе непосредственно связана с количеством узлов, выбираемых для построения дискретной модели, и точностью задания нагрузок р либо перемещений и. По этой причине при расчете реальных конструкций с многоочаговыми концентраторами приходится, если это вообще возможно, значительно расширять рассматриваемую область до размеров, где проявление концентрации сглаживается, что приводит к значительным трудозатратам по составлению исходной информации для решения задачи и накладывает ограничения на класс применяемых ЭВМ.
Использование данных, полученных с помощью голографической интерферометрии, для задания граничных условий вида (3) в окрестности концентратора напряжений позволяет существенно уменьшить рассматриваемую область решения. При этом возможны различные
А, и=ь,
физическим соотношениям в области V: .
С г = а,
(4)
(5)
Д, С />0 £/ = / — в области Т' М С О0 и — /— на границе [/2 = ио — на границе
(6)
щи)=\<ю0 и (IV— Л ит/йУ- [ и^рйБ,
(7)
V
V
варианты подхода к проведению эксперимента. Однако в любом случае целесообразно на натурный объект или его модель наносить сетку, в узлах которой определяются перемещения, адекватную конечноэлементному разбиению. Компоненты вектора узловых перемещений И0 на границе 52 желательно измерять во всех контурных узлах. Если же это сделать не удается (часто в силу технических причин), то недостающие компоненты вектора £/0 приходится, получать с помощью интерполяции.
2. Точность определения НДС, а следовательно, практическая целесообразность применения комбинированного расчетно-экспериментального метода, зависит от точности измерения перемещений на границе 32. Среди известных подходов к интерпретации картин интерференционных полос наиболее универсальным является способ расшифровки по абсолютным порядкам [7]. В рамках этого подхода неизвестный вектор перемещений X произвольной точки исследуемой поверхности определяется из решения системы линейных алгебраических уравнений вида:
КХ = 1Ы, (8)
где К — матрица чувствительности размерности 3X3, элементы которой равны проекциям разности единичных векторов наблюдения и освещения на координатные оси х;(г'= 1, 2, 3); Х=\Хи Х2 Хд}; ./V = {Л^, Л^, ЛМ— вектор абсолютных порядков полос; X —длина волны источника освещения.
Наиболее эффективным практическим инструментом для оценки погрешностей компонентов вектора перемещений АХ, (¿=1, 2, 3) и их минимизации является подход, основанный на принципах нестатистического планирования эксперимента ¡[8, 9]. Основная идея этого подхода заключается в том, что для измерения компонентов вектора X с наибольшей точностью необходимо минимизировать число обусловленности матрицы К, которое в пределе равно единице. При расшифровке интерферограмм по абсолютным порядкам полос чувствительность к компоненте перемещений, нормальной к поверхности исследуемого объекта, превосходит чувствительность к касательным компонентам. В этом случае матрица чувствительности интерферометра неравновесна, а оценки погрешностей по числу обусловленности становятся весьма консервативными. Поэтому для анализа точности измерений матрицу чувствительности целесообразно представить в виде:
где />! и — диагональные матрицы размерности 3X3, которые выбираются так. чтобы элементы масштабированной матрицы чувствительности К* имели одинаковый порядок величин.
Погрешности определения компонентов Х1 будут минимальны для голографического интерферометра с матрицей чувствительности /С, для которой существуют матрицы и такие, что матрица К* становится ортогональной. Схема подобного интерферометра представлена на рис. 1. На этом рисунке ось хи вдоль
которой направлен единичный ‘вектор освещения £$, совпадает с нормалью к исследуемой поверхности в точке О. Для обеспечения ортогональности матрицы К* необходимо, чтобы углы а;. и
определяющие направления единичных векторов наблюдения e¡ (/'= 1, 2, 3), были связаны соотношениями:
<!>! = <1)2; <Xj + «2 + 2а3 =
№ ■
S 2 . » ё 2 ’ sin —
2
(0 = О < 6 < 90°);
0 < < 90°;
0 < <pr < 2arctg
COS i
(9)
т = 1, <7 = 3, г—1, если 0 < 6 <60°, и т —— 1, <7 = 1, г — 3, если 60< 9 < 90°.
Элементы масштабирующих матриц Ох и /)2 выглядят следующим образом:
А(1) = Д(2)=1; А(3) =
cos2 — 4-
2 sin2
cos2 *
sin2
DA 1)=—
cos 0 l^cosI
]/ 2 cos — (1 + cos tpj)
D2 (2) = D2 (3)
(10)
Основной вклад в погрешности АХ1 {¿= 1, 2, 3) вносят ошибки в измерении компонентов вектора правой части системы уравнений (8). Вклад погрешностей определения элементов матрицы К при-
мерно на порядок меньше и на практике его обычно не учитывают. Тогда, если вектор чувствительности (е} — е8) постоянен по всей исследуемой поверхности (что обеспечивается при коллимированном освещении объекта и телецентрическом наблюдении по всем трем направлениям), величины погрешностей компонентов вектора перемещений можно оценить по формулам:
|ДЛГ,|<£2(»)ХД п, 1 = 1, 2, 3, (11)
где Д> (г) — элементы матрицы — (10); Д п — погрешность измерения абсолютных порядков полос, которая предполагается равной
для всех направлений е; (у = 1, 2, 3).
Если же вектор чувствительности меняется по поверхности, например в распространенных на практике случаях, когда расстояния от источника освещения и точек наблюдения до объекта конечны, то оценки (11) справедливы только для начала координат О (см. рис. 1). Однако и в этом случае в окрестности точки О будут обеспечены условия оптимальности измерений в рамках ограничений, накладываемых на параметры конкретного интерферометра. Для величин погрешностей АХг справедливы следующие оценки:
| Д Л^С/З Ог(1) сопс! К*—, Ал, ¿=1,2,3, (12)
I «I ^ 2 V / [| О* || ЦЛИ1
где || • || — нормы векторов и матриц, согласованные друг с другом; соп<1 К* = ||/С*|||1(/С*)_1|| — число обусловленности матрицы К*\ /)2(г)— элементы матрицы /)2, которые определяются по формулам (10); Х* = О^Х, N* = 0^ — масштабированные вектора перемещений и абсолютных порядков полос соответственно; О*—вектор с компонентами 6,:
2’3
(к}/ — элементы матрицы К)-
Рассмотренные оптимальные схемы интерферометров, параметры которых удовлетворяют соотношениям (9), представляет собой метрологическую базу, необходимую для практического использования расчетно-экспериментального метода анализа НДС конструкций.
3. Принципиальным моментом при анализе эффективности расчетно-экспериментального подхода для определения НДС является возможность проведения тестовых экспериментов. Это, связано с тем, что для получения точностных характеристик метода и рекомендаций по его применению необходимо сравнение полученных в подобных экспериментах результатов с известными теоретическими решениями или с данными других экспериментальных методов. В рассматриваемом подходе проблема проведения корректных тестов занимает особое место вследствие высокой (порядка длины волны видимого света) чувствительности голографических интерферометров к деформационным перемещениям и особенно к паразитным смещениям исследуемого объекта, как целого. Все это предъявляет жесткие требования к обеспечению адекватности условий деформирования и схемы нагружения образца в эксперименте,тем условиям, при которых получены теоретические и экспериментальные результаты, используемые для сравнения.
В данной работе объектом исследования служит тонкостенная ци-диндрическая оболочка с круговым вырезом. Выбор этого распространенного в летательных аппаратах конструктивного элемента, во-первых, дает возможность воспроизвести наиболее неблагоприятную с метрологической точки зрения ситуацию, которая может встречаться на практике, — а именно, наличие больших градиентов деформаций на малых по площади участках криволинейной поверхности — и, во-вторых, позволяет использовать для сравнения результаты как теоретических [10, 11], так и экспериментальных [12] исследований.
Реальный объект из материала Д16Т с цилиндрической рабочей частью длиной 100 мм (внешний радиус ^? = 30 мм, толщина стенки /2=1,5 мм), где расположено одно отверстие (радиус г0=\2 мм), выполнен в виде специального образца (рис. 1), что дает возможность для его нагружения растягивающими усилиями применять прецизионное устройство, результаты экспериментальной проверки которого приводятся в работе ;[8]. Эти результаты свидетельствуют, что устройство обеспечивает одноосное растяжение рабочей части образца без ее смещений, как целого, с интерференционной точностью. При нагружении массивная часть образца, обозначенная на рис. 1 буквами аа, жестко закреплена и не смещается, а растягивающие усилия равномерно распределены по окружности на фланцевой части (сечение ЬЬ).
4. При исследовании тонкостенных элементов конструкций с помощью МКЭ объемные интегралы в выражении (7) заменяются поверхностными, а поверхностные — контурными. Соответственно граница области V—52, где задаются перемещения при решении системы уравнений (6), представляет собой контур на поверхности объекта.
Контур, который используется для определения НДС в окрестности отверстия расчетно-экспериментальным методом, нанесен на поверхность оболочки и обозначен на рис. 1 буквами АВС^И. Этот контур отвечает полугеодезической системе координат (ки ср, х^ с началом в центре отверстия О (ось х1 направлена от плоскости чертежа — рис. 1). Шаг сетки в направлении ср равен я/8 в направлении — 5 мм.
Для измерения перемещений применяется голографический интерферометр, схема которого и ее привязка к поверхности исследуемого объекта показаны на рис. 1 (начало координат О, совпадает с центром отверстия на поверхности оболочки). Параметры оптимальной схемы интерферометра, согласно соотношениям (9), выбираются следующим образом:
-]>1 = <1»2 = 25°; ф3 = 58,6°; а1 = 26,6°; а2 = 90° — а^, а3 = 45°,
что обеспечивает одновременное наблюдение контура АВС^И по всем направлениям е,(/=1, 2, 3). Точечный источник освещения расположен на оси Х\ на расстоянии 0,83 м от исследуемого объекта. Расстояния от точек наблюдения, определяемых векторами и е2, до начала координат — 0,4 м, от точки наблюдения с вектором еъ — 0,3 м.
Величины компонентов вектора перемещений ^{(1=1, 2, 3) в узлах сетки на контуре АВСопределяются при приращении нагрузки (2=1185±4Н (первая экспозиция проводится для ненагруженного образца). Соответствующая интерферогрямма оболочки, полученная по направлению наблюдения е3, показана на рис. 2. Справа на фотографии виден упругий резиновый мостик, который служит для иденти-
Рис. 2
фикации порядкового номера и знака интерференционных полос. Абсолютные порядки в узлах координатной сетки определяются с помощью линейной аппроксимации.
Распределение компонента перемещения Х2 на контуре ABCi и распределение компонентов Xi и Х3 по контуру BCi представлено на рис. 3 кривыми 2, 1 и 3 соответственно. Эти значения используются в дальнейшем в качестве граничных условий при анализе НДС на контуре отверстия AD с помощью расчетно-экспериментального подхода. Погрешности определения компонентов, оцененные по неравенству (12) для Art = 0,15, не превосходят 0,1 мкм для Х{ и 0,3 мкм для Х2 и Х3.
5. Применение комбинированного подхода к анализу НДС элементов конструкций предполагает наличие соответствующего расчетного инструмента. В этом качестве используется специализированный комплекс программ расчета местной прочности СКП ФИТИНГ [1]. Комплекс представляет собой совокупность сервисных программ и рабочих модулей, предназначенных для решения задач теории упругости, пластичности и механики разрушения. Обмен между модулями осуществляется посредством банка данных, который содержит всю необходимую для работы каждого модуля информацию. Там же после окончания работы модуля хранится выходная информация. Модульное построение СКП ФИТИНГ позволяет реализовать автономную работу отдельных блоков и предусматривает дальнейшее развитие комплекса посредством включения новых модулей.
Математическое обеспечение СКП ФИТИНГ написано на языке ФОРТРАН-IV и адаптировано на ЭВМ БЭСМ-6 и ЭВМ серии ЕС. В комплексе реализован МКЭ в сочетании с методом подконструкций, что позволяет рассматривать широкий класс задач местной прочности, включающий элементы конструкций сложной формы, разнообразные кинематические и статические граничные условия. Аппроксимация конструкций осуществляется на основе одномерных, двумерных и объемных конечных элементов, позволяющих решать задачи в рамках теории упругости или строительной механики.
В настоящем исследовании для анализа НДС тонкостенных конструкций применяются треугольные несовместные элементы, воспринимающие как мембранные, так и изгибные нагрузки, которые были впервые предложены в работе [13] и получили развитие для расчета оболочек в работе [6].
Для того, чтобы проверить работоспособность подобных элементов в окрестности концентратора напряжений и оценить степень дискретизации, необходимую для получения приемлемой точности, на первом этапе проводится конечноэлементный расчет исследуемой оболочки. При расчете, учитывая геометрическую симметрию рабочей части образца, рассматривается четвертая часть оболочки, разбиение .которой на конечные элементы показано на рис. 4.
Решение задачи проводится с помощью двух подконструкций, которые образованы разделением четверти оболочки на две равные части по линии МХН. Первая подконструкция (АМ^РИ) включает 289 узлов — контур отверстия АО разбит на 16 равных частей; также на 16 отрезков, шаг которых увеличивается при удалении от отверстия, разбиты линии АМХ и /Ж Каждый четвертый узел на линиях АИ и ВСх соответствует узлам координатной сетки, в которых определяются перемещения методом голографической интерферометрии (см. рис. 2). Вторая подконструкция (ШЬМ^ состоит нз 112 элементов, которые получаются при равномерном разбиеции линий ЬМ1 и Мна 7 и 8 частей.
Рис. 4
Нагрузка в виде единичного растягивающего напряжения вдоль оси х2 задается на контуре /-У. На этом же контуре, исходя из реальных условий нагружения, приравниваются нулю компоненты вектора перемещений Х1 и Х3, а также углы поворота сечений вокруг осей х1 и х3. Контур АЬ закрепляется в направлении оси х2, т. е. принимается гипотеза Х% = 0 (что не соответствует условиям деформирования реальной детали в окрестности отверстия). Углы поворота сечений на этом контуре вокруг осей х1 и х3 считаются равными нулю в силу симметрии задачи. По соображениям симметрии принимаются равными нулю также компонент вектора перемещений Х3 и углы поворота сечений вокруг оси х2 на линиях Ы и /Ж
Полученное в результате расчета распределение мембранной о” и изгибной а» составляющих нормальных напряжений на контуре отверстия АО для единичного модуля упругости Е и коэффициента Пуассона у = 0,3 представлено на рис. 5 кривыми 3 к 4 соответственно. На этом же рисунке показаны мембранная и изгибная составляющие
нормальных напряжений для аналогичной ^ параметр р = 1,17, р2 =
= Г°^—щ---------— ) бесконечной оболочки, приведенные в работе
[10] (кривая 2 и значки А). Значками • и О обозначены данные по мембранным и изгибным напряжениям соответственно из работы
[11] для бесконечной оболочки с параметром г0/1//?й = 1,81 при > = 0,3.
Максимальная величина мембранных напряжений из МКЭ-расчета хорошо согласуется с экспериментальной, полученной с помощью метода оптически чувствительных покрытий [12]. Сравнение перечисленных результатов свидетельствует, что выбранный тип конечного эле-
Рис. 5
мента и степень дискретизации исследуемого объекта обеспечивают достаточную точность определения концентрации напряжений с помощью МКЭ. Наблюдаемые различия в величинах напряжений незначительны и, по-видимому, связаны с влиянием конечных размеров образца.
Для анализа НДС на контуре отверстия расчетно-эксперименталь-ным методом используется фрагмент конечноэлементной сетки из 119 узлов, граница которого S2 обозначена на рис. 4 буквами ABCJ3. При этом с целью минимизации количества необходимой экспериментальной информации три компоненты вектора перемещений узлов задаются только на части контура ВС1. На линии ЬС\ сохраняются кинематические граничные условия, которые применяются при МКЭ-расче-те четверти оболочки. Перемещения на контуре AD не описываются. На линии АВ, лежащей в плоскости симметрии оболочки, задается компонент перемещения Х2, так как при растяжении реального образца, закрепленного в сечении аа (см. рис. 1), ее точки перемещаются в направлении оси х2 на различные расстояния, т. е. условие Х2 = const не выполняется.
Компоненты вектора перемещений в узлах конечноэлементной сетки, не совпадающих с узлами координатной сетки, где определяются перемещения методом голографической интерферометрии, необходимые для формирования на границе S2 вектора узловых перемещений
и0 из системы уравнений (6), устанавливаются посредством интерполяции. Интерполирующие кривые, показанные на рис. 3, проводятся через экспериментальные точки с помощью сглаживания от руки. Определенные таким образом компоненты вектора перемещений X в
узлах сетки после нормировки на величину-----------------------(Е= 7Х
ХЮ4 МПа) формируют используемый при расчете в качестве граничных условий вектор узловых перемещений 1/0.
6. Мембранная о“ и изгибная а" составляющие нормальных напряжений на контуре отверстия, полученные с помощью комбинированного подхода, показаны на рис. 5 кривыми 1 и 5 соответственно. Эти результаты, как показывают оценки, в пределах погрешностей, обусловленных ошибками в задании граничных условий, геометрических размеров образца (Я, /г, г0), механических свойств материала (Е, у) и нагрузки (2, перекрываются с результатами МКЭ-расчета четверти оболочки практически на всем контуре АО. Следует, однако, отметить полученное при комбинированном подходе смещение точки, где достигаются максимальные растягивающие мембранные напряжения, внутрь контура АО. Такой факт может быть связан с реальными условиями закрепления образца, но, с другой стороны, может иметь место и влияние погрешностей, с которыми задаются компоненты вектора узловых перемещений ио в направлении оси х2 на линии АВ.
Вообще говоря, вопрос о влиянии погрешностей задания вектора [/0 на точность определения НДС связан, кроме точности экспериментального определения компонентов Хь с проблемой рационального выбора контура 52. Эта проблема состоит в том, что с точки зрения трудозатрат желательно уменьшать область V в окрестности концентратора, ограниченную контуром 32. Однако очевидно, что существует критический размер области V, за которым погрешности задания вектора 1/0 могут оказать значительное влияние на результаты определения деформаций и напряжений. Анализ этого вопроса, так же как и вопроса выбора наилучшего способа интерполяции экспериментальных данных в узлах конечноэлементной сетки, требует проведения широкого круга параметрических исследований, которые выходят за рамки данной работы.
Результаты исследования НДС оболочки с отверстием показывают, что расчетно-экспериментальный метод является эффективным средством определения концентрации напряжений в реальных тонкостенных элементах конструкций. Использование данных голографических измерений, полученных с помощью оптимальных схем интерферометров совместно со специализированным комплексом программ расчета местной прочности ФИТИНГ, обеспечивает точность определения напряжений, сопоставимую с точностью известных аналитических, расчетных и экспериментальных методов. Применение рассмотренного комбинированного подхода наиболее перспективно при анализе НДС элементов конструкций с многоочаговыми концентраторами напряжений.
ЛИТЕРАТУРА
1. Барышников В. И., Гришин В. И., Донченко В. Ю., Тихонов Ю. В. Применение метода конечных элементов к исследованию местной прочности элементов авиационных конструкций.—Ученые записки ЦАГИ, 1983, т. XIV, № 1.
2. П р е й с А. К. Определение напряжений в объеме детали по данным измерений на поверхности. — М.: Наука, 1979.
3. Barishpolsky В. M. A combined experimental and numerical method for the solution of generalized elasticity problems. — „Experimental Mechanics“, 1980, vol. 20, N 10.
4. Pa с с о x а А. А. Применение голографической и спекл-интерферо-метрии для определения напряженно-деформированного состояния элементов конструкций, дефектоскопии и дефектометрии. — В кн.: Физические, основы голографии. — Л., 1978.
5. Галкина Н. С., Г ришин В. И., Донченко В. Ю. Исследование напряженно-деформированного состояния элементов авиационных конструкций и их соединений. — Труды ЦАГИ, 1979, вып. 2012
6. 3 е н к е в и ч О. С. Метод конечных элементов в технике. — М.: мир, 1975.
7. В г i е г s J. D. The interpretation of holographic interferograms (review). — „Optical and Quantum Electronics“, vol. 8, N 6, 1976.
8. П и с a p e в В. С., Яковлев В. В., Щепинов В. П. Оценка точности определения компонентов вектора перемещений в методе голо-графической интерферометрии. — В кн.: Физика и механика деформации и разрушения. — М.: Энергоиздат, вып. 9, 1981.
9. П и с а р е в В. С., Яковлев В. В., Индисов В. О., Щепинов В. П. Планирование эксперимента по определению деформаций методом голографической интерферометрии. — Журнал технической физики, т. 53, вып. 2, 1983.
10. V a n Dyke P. Stresses about a ascular hole in cylindrical shell. — AIAA Journal, vol. 3, N 9, 1965.
11. Г узь A. H., Чернышенко И. С., Чехов Вал. Н. и др. Цилиндрические оболочки, ослабленные отверстиями. — Киев, Наукова думка, 1974.
12. Александров А. Я., Ахметзянов М. X., Ракин А. С. Исследование упруго-пластического деформирования оболочек с вырезами и усилениями методом фотоупругих покрытий. — «Прикладная механика», т. 2, вып. 3, 1966.
13. В a z е 1 е у G. P., Cheung J. К., Ivons В. М., Zienkie-w i с z О. С. Triangular elements in bending-conforming and nonconforming solutions. — Proc. Conf. Matrix Methods in Struct. Mech., Air Fovce Inst, of Tech., Wright Patterson A. F. Base, Ohio, October, 1965.
Рукопись поступила 14/IV 1983