Расчет коэффициентов интенсивности напряжений в вершине сквозной трещины с помощью специализированного гибридного конечного элемента Текст научной статьи по специальности «Физика»

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

Том XVIII 1 987 Мб

УДК 539.4.013.3 : 624.073

РАСЧЕТ КОЭФФИЦИЕНТОВ ИНТЕНСИВНОСТИ НАПРЯЖЕНИЙ В ВЕРШИНЕ СКВОЗНОЙ ТРЕЩИНЫ С ПОМОЩЬЮ СПЕЦИАЛИЗИРОВАННОГО ГИБРИДНОГО КОНЕЧНОГО ЭЛЕМЕНТА

А. Д. Дементьев

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

Оценена точность вычисления КИН с помощью предложенного специализированного элемента и влияние на нее параметров этого элемента.

Для расчета характеристик безопасной повреждаемости элементов конструкций (остаточной прочности и длительности роста трещин) широко используется линейная механика разрушения, в которой определяющую роль играет коэффициент интенсивности напряжений (КИН). В расчете КИН для элементов конструкции сложной геометрии и/или сложного нагружения применяются численные методы, в частности, метод конечных элементов (МКЭ). Применение обычных конечных элементов для расчета КИН наталкивается на трудности, связанные с непригодностью этих элементов для эффективного учета особенности в поле напряжений у вершины трещины типа 1 ¡Vг, где г — расстояние от вершины трещины.

Данное обстоятельство приводит к необходимости очень мелкой дискретизации области у вершины трещины на конечные элементы.

С целью уменьшения степени дискретизации используются .энергетические методы, например, метод виртуального роста трещины [1], либо специализированные конечные элементы, учитывающие особенность напряженно-деформированного состояния (НДС) в вершине трещины. В нашей стране и за рубежом для сквозных трещин в пластине, находящейся в условиях плоского напряженного состояния или плоской деформации, предложен ряд таких элементов [2—9], различающихся точностью аппроксимации НДС, типом функционала, используемого для получения матрицы жесткости, формой элемента, степенью выполнения условий совместности деформаций с окружающими конечными элементами.

Однако исследования предложенных в перечисленных работах элементов проведены, к сожалению, на достаточно узком классе тестовых задач, как правило, на полосе с центральной или краевой трещиной, находящейся под действием равномерного растяжения. Недостаточно изученным является также влияние на точность расчета КИН других факторов, например, размеров и формы элемента, положения кончика трещины внутри элемента, несовместности перемещений в промежутке между узлами по границе некоторых специализированных элементов [2—4]. Отсутствие данных более детального исследования этих элементов затрудняет оценку точности расчета КИН в реальных конструкциях.

В настоящей статье предложен совместный специализированный конечный элемент, позволяющий по сравнению с элементами [2—4, 6—7] увеличить точность решения задачи о трещине. При этом высокая точность расчета КИН достигается, как будет показано ниже, на грубой конечно-элементной сетке. Применение предложенного элемента требует меньшей вычислительной трудоемкости по сравнению с элементом работы [5]. В отличие от элементов работ [8, 9] данный элемент имеет произвольную форму и число узлов на каждой его стороне, что увеличивает эффективность применения элемента.

В данной статье для вычисления матрицы жесткости специализированного конечного элемента используется более простой гибридный функционал, а также приводятся результаты более широкого исследования этого элемента, необходимые для оценки точности расчета КИН в реальных конструкциях.

Матрица жесткости специализированного конечного элемента с трещиной. Для вывода матрицы жесткости элемента с трещиной используется гибридный функционал, основанный на выборе распределения

напряжений [10].

Применение этого функционала позволяет обеспечить совместность деформаций по границе между специализированным и обычными элементами. Функционал имеет следующий вид [10]:

п=£п„

п

пл = § 9*Е(Т (IV + $ Т*ий8- | 7*и^5, (1)

2 *п Зп Зш

где я —статически возможный вектор напряжений; а* = [<зх <зу тлу], (* — означает операцию транспонирования);

Е—матрица упругих постоянных,

еи ^12 «18

Е = ^12 Є22 Є23

^18 е23 ^33

Уп — область элемента с номером п\

Т — вектор усилий, действующий по границе элемента,

Т* = [7х Т2\;

Ти Т2— проекции вектора Т на координатные оси х и у;

и — вектор перемещений на границе элемента 5Я, «* = \и1 и2];

и2 — проекции вектора перемещений на координатные оси, причем, и = и3 на 5л2;

— часть границы элемента, на которой задан вектор перемени й и$\

Т$ — вектор заданных усилий, действующих на части границы элемента 5п1;

^ = 5л1+5л! + 5„ з,

'„з — часть границы элемента, по которой он соединяется с соседними конечными элементами.

Условие стационарности функционала при варьировании векторов [ и приводит к следующим уравнениям Эйлера и естественным гра-1ным условиям:

Ез = Аи в области Уп\ и = и на границе

7=7, на 5„

7+ = 7- на 5я3;

(2)

д/дх О О д/ду д/ду д/дх

и — вектор перемещений в области , 7- — векторы усилий справа и слева от межэлементной границы.

Выберем вектор напряжений так, чтобы удовлетворить первому из авнений (2), тогда, используя формулу Остроградского—Гаусса, по-шм для случая отсутствия объемных сил:

П„= - — / (7* и + и* 7) с1Б + | 7* и йБ - ] 7* й йЗ .

(3)

Примем векторы 7, и, и в виде

Т = (4)

я=^р; (5)

и = (6)

Э Р — вектор неизвестных коэффициентов; # — вектор узловых ремещений элемента; Я, и, I* — аппроксимирующие матрицы. Подставив уравнения (4), (5), (6) в (3), будем иметь:

иъ и2 — проекции вектора перемещений на координатные оси, причем, и = Юу на 5л2;

5л2 — часть границы элемента, на которой задан вектор перемещений па\

Т5 — вектор заданных усилий, действующих на части границы элемента 5я1;

^л = ^/11 + ^п2 + з.

5Л з—часть границы элемента, по которой он соединяется с соседними конечными элементами.

Условие стационарности функционала при варьировании векторов а и и приводит к следующим уравнениям Эйлера и естественным граничным условиям:

Е<з = Аи в области 1/„;

где

и —и

Т=Т3

Т+ = Т-

А =

на границе 5Л; на 5л1; на 5л3;

д1дх О О д/ду д/ду д/дх

(2)

и

Т+, 7

вектор перемещений в области Уп\

-векторы усилий справа и слева от межэлементной границы.

Выберем вектор напряжений так, чтобы удовлетворить первому из уравнений (2), тогда, используя формулу Остроградского—Гаусса, получим для случая отсутствия объемных сил:

П „ =

/ (7* и + и* Т)(18+ ¡Т* и йЭ- ] Т*3 и йБ

(3)

Примем векторы 7, и, и в виде

7 = /?р; и=и$-, и = Ьд,

(4)

(5)

(6)

где Р — вектор неизвестных коэффициентов; (7 — вектор узловых перемещений элемента; /?, £/, Ь — аппроксимирующие матрицы. Подставив уравнения (4), (5), (6) в (3), будем иметь:

пР* ЯР + РЧ?«-/7*,

(?)

0=//?^й5; (9)

«л

/=■= I т1ь (¡Б . (10)

5Л1

Учитывая независимость вектора $ от номера элемента, с помощью условия стационарности функционала исключим его из уравнения (7):

НР-Од = 0; ] (И)

$ = Н-Хвд .

Подставляя уравнение (И) в (7), получим:

П п=±д*Кд-Рд,

где К — матрица жесткости элемента,

К=в* Н~10

Т7—вектор узловых сил элемента.

Для построения матриц /?, £/ используется решение задачи о полу-бесконечной трещине в анизотропной пластине. Это решение получено с помощью методов Лехницкого [11] в работе [10] и с некоторой модификацией используется в данной статье. Выражения для напряжений и перемещений согласно этому решению имеют вид:

а, = Ие £ к {С?"8 + (*2 (С + ( -1)* В) Й"2] а, * + к

+ * И С?-2 + $ (С - (-1)* В) $~2} а2к} ; а, = Ие ^ к { [СГ2 + (С+ ( -1)* В) С!“2] а, * +

*

+ ¿[С?-2 + (С-(-1)* Я) Й“8] а**} ;

Ие ^ к { [ 1^1 С?"2 + ц2 (С+ (-1)* В) С 2 2] а, к + к

+ ¿[¡.1 С?-2 + ца (С —( — 1)* В) С2Й-2] а2к\ ;

«. = 2^2 {[р^ + Р2(С + {-\)ЬВ)$\а1к + к

+ *[м?+л (С-(- 1)* В)й] а2й};

м2 = 2Не£{ [</,С? + £/2(С + (-1)*Д)С8*] я,* +

к

-И [<*,<:? + </, (С-(-1)*В) й] а2*}, где Р!, [а2, [а,, ц2 —корни характеристического уравнения е\\ I1,4 2е13 ¡I.3 -(- (2е12 4- 033) — 2е23 (I. + е2а = О

С* = **, ** = *£=1,2;

(12)

(13)

^■2 — М-1 .

_ У

{*•2 — 1*2

1*2 — И .

М-2 — ¡*2

Рь — еп Iх* 4“ е{2 ~~ е\ъ Рк;

е12 У-\ + «22 — ^23 Рк

£ = 1,2.

Приведенное решение справедливо для случая Решение для

изотропной пластины получается путем задания слабо анизотропного материала, что достигается изменением модуля сдвига на 0,1%. Вектор ¡3 в уравнениях (4) и (5) составляется из величин а^, а2и, удерживая равное их количество, с исключением коэффициента а22, который как и а12 соответствует постоянным напряжениям Ох для анизотропного материала и не дает вклада в напряженное состояние для ортотропного материала. Число элементов вектора р должно удовлетворять неравенству Щ>2ЫЯ— Л/т, где Л/д — число узлов в элементе, — число степеней свободы элемента как жесткого тела. Это неравенство обеспечивает дефект матрицы жесткости, равный ЛГТ [10].

Аппроксимирующая матрица Ь выбирается так, чтобы обеспечить совместность перемещений по границе между элементами.

КИН первого и второго рода Кг и Ки для сквозной трещины в анизотропной пластине определяются выражением [12]:

В расчетную систему МАРС [13] были введены два специализированных элемента с трещиной, матрицы жесткости которых вычисляются по изложенному выше алгоритму. Элемент для несимметричных задач имеет форму криволинейного восьмиугольника (рис. 1 , а), элемент для -симметричных относительно плоскости трещины задач соответствует

половине первого элемента (рис. 1,6). Число узлов на каждой стороне многоугольника может меняться от двух до пяти,, задавая таким образом закон аппроксимации перемещений и координат вдоль нее. При построении матриц Л, и для второго элемента используются только симметричные относительно плоскости трещины члены уравнений (12)

Параметрические исследования. В параметрических исследованиях использовались специализированные элементы с линейной и квадратичной аппроксимацией перемещений по контуру элемента. Проведены исследования изолированного специализированного элемента и исследования его в составе конечно-элементной модели пластины. Первые осу-

Рис. 1

и (13).

ществлялись путем нагружения этого элемента узловыми перемещениями, которые соответствовали точным перемещениям этих точек в задаче всестороннего растяжения и сдвигового нагружения бесконечной пластины с трещиной, длиной 21 (рис. 2,а), либо изгиба бесконечной пластины с двумя полубесконечными трещинами в своей плоскости (рис. 2,6). Точное решение этих двух задач получено с помощью функций напряжений Вестергаарда [14].

При исследовании изолированного элемента оценивалось влияние различных факторов на точность вычисления КИН. Как показали расчеты, погрешности расчета КИН, равные

Д,=

*і-*Г

хї

*н - *1,

К],

где К], К— теоретические значения КИН первого и второго рода соответственно, практически не меняются при числе коэффициентов Np >10 и числе точек интегрирования методом квадратур Гаусса по

Л1б

в)

м

21

М

Рис. 2

каждой стороне Ыи>>3. Следует отметить, что переход от интегрирования по области к интегрированию по контуру при вычислении матрицы жесткости элемента позволяет избежать влияния на точность интегрирования особенности НДС в вершине трещины. Зависимость погрешностей расчета КИН от относительных размеров элемента и положения кончика трещины внутри элемента для изотропной пластины приводится на рис. 3 и 4. Исследованные элементы показали достаточно высокую с практической точки зрения точность, особенно элемент с квадратичной аппроксимацией перемещений по контуру.

Конечно-элементные расчеты пластин с трещиной проведены с использованием МКЭ в форме метода перемещений, обычные конечные элементы имели линейную аппроксимацию перемещений по элементу, а специализированные — линейную аппроксимацию перемещений по контуру. В ходе этих расчетов решена задача о равномерном растяжении изотропной и ортотропной пластины с центральной трещиной, а также изотропной пластины с трещиной, расположенной под углом 45° к направлению растяжения.

Размеры изотропной пластины с центральной трещиной составляли 160X480 мм, ортотропной — 200x200 мм, размеры пластины с наклонной трещиной — 200X400. В силу симметрии первой задачи рассматри-

Рис. 3

а! и = 1, а — 20 мм, — 18 Изолированные испытания элемента: растяжение и сдвиговое нагружение плас-

О-----------О А»

тины,

О—О •--• а

а=1, Т-1, 7=90°, т= 1 (см. рис. 1, а) изгиб пластины, М = 1 (см. рис. 1, б) Исследования элемента в составе конечно-элементной модели пластины:

--------- Д1 растяжение пластины,

сг = 1 (см. рис. 5, а)

Рис. 4

а — 20 мм, I — 40 мм, - 18 Изолированные испытания элемента: растяжение и сдви-д, | говое нагружение пластины,

0=1, т = 1, 7=90°, т-1 (рис. 1, а)

Исследования элемента в составе конечно-элементной модели пластины:

___ . ---- дх растяжение пластины,

о = 1 (рис. 5» я)

___ .. --- д, растяжение пластины,

о = 1, сторона № 1 элемента (рис. 1, б) выходит на свободную границу пластины

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

Результаты расчета приводятся на рис. 3 и 4 и слабо зависят, как показали исследования, от абсолютных размеров элемента в диапазоне изменения параметра а, равном 3-^-30 мм. В большей части проведенных расчетов использовались специализированные элементы квадратной формы (либо половина квадрата для симметричных задач), однако как показали численные эксперименты, точность вычисления КИН слабо менялась даже при значительном отклонении формы элемента от квадратной.

Конечно-элементная модель ортотропной пластины подобна модели, приведенной на рис. 5, а. Комплексные параметры анизотропии в этом расчете выбраны равными = I (¿ — мнимая единица), ^2 = 1$,

б

\ ! І

є

( І \

б

\ \ і

а)

Ю

і)

\ \ \

Рис. 5. Конечно-элементные модели пластин с трещиной

1—специализированный элемент; 2—поверхность трещины

82=Е1/Е2, где Ей Е2-—модули упругости вдоль координатных осей х и у соответственно (ось х направлена вдоль трещины). Результаты расчета (поправка Я, учитывающая конечные размеры пластины), полученные при неизменной конечно-элементной сетке, приводятся в таблице вместе с поправкой Ни, вычисленной методом коллокаций [15]. Расчетная модель пластины с прямолинейной наклонной трещиной представлена на рис. 5, в. Длина трещины при одной и той же конечноэлементной сети менялась в диапазоне 20—60 мм.

£і/£2

1, мм 0, 5 0,9 1,5 4,5

н нк Н нк Н Н нк

40 1,29 1,30 1,23 1,23 1,19 1,18 1.14 1,12

50 1,39 1,46 1.31 1,35 1,25 1,28 1,19 1,20

60 1,56 1.64 1,44 1,50 1,36 1,41 1,27 1,31

Результаты расчета сравнивались с точным значением КИН для наклонной трещины в бесконечной пластине [14]. Отклонение расчетных значений КИН первого и второго рода для коротких относительно ширины пластины трещин (20—40 мм) не превысило 3%, а для ¡более длинных (50—60 мм), где сильнее сказывается влияние конечных размеров пластины, отклонение составило около 5%.

Проведенные исследования показали, что при относительном размере элемента 1/а, равном 1—6, и при относительной координате положения кончика трещины внутри элемента ^/а, равной 0,05—0,75, погрешность вычисления КИН первого и второго рода не превышает 1—7% и 1—8%. соответственно. Эффективность предложенного элемента позволяет рекомендовать его для решения практических задач.

I. Гришин В. И., Донченко В. Ю. Метод расчета коэффициентов интенсивности напряжений в элементах конструкции с криволинейными трещинами. — Ученые записки ЦАГИ, 1983, т. 14, № 2.

2. Д а ш е в с к и й Е. М. Анализ поля напряжений в окрестности трещин методом конечных элементов. — Строительная механика и расчет сооружений, 1973, № 4.

3. В у s к о v Е. The calculation of stress intensity factors using the finite element method with cracked elements. — International Journal of fracture Mechanics, 1970, vol. 6, N 2.

4. Aberson J. A. Anderson J. M. Cracked finite elements proposed for NiASTRAN. — NASA TMX-2893.

5. Партон В. 3., Борисковский В. Г. Динамическая механика разрушения. — М.: Машиностроение, 1985.

6. Н е n s h е 11 R. D., Show К. G-, Crack tip finite elements are unnecessary. — International Journal for numerical Methods in engineering, 1975, vol. 9, N 3.

7. Морозов E. М., Никишков Г. П. Метод конечных элементов в механике разрушения. — М.: Наука, 1980.

8. Tong P., Pi an Т. Н. Н., Lasry S. J. A hybrid-element approach to crack problems in plane elasticity. — International Journal for numerical Methods in engineering, 1973, vol. 7, N 3.

9. T о n g P. A hybrid crack element methods for solid continua. — International Journal for numerical Methods in engineering, 1969, vol. 1, N 2.

II. Лехницкий С. Г. Теория упругости анизотропного тела. — М.: Наука, 1977.

12. Парис П., Си Д ж. Анализ напряженного состояния около трещины. — В сб. «Прикладные вопросы вязкости разрушения. — М.: Мир, 1968.

13. Г а л к и н Д. С., Галкина Н. С., Гусак Ю. В., Зайцев С. Н., Иванов А. И., И в а н т е е в В. И., Кудряшов А. Б., Литвиненко А. А., Полищук В. А., Чу бань В. Д., Шевченко Ю. А. Многоцелевая автоматизированная расчетная система МАРС. Комплексы программ математической физики (материалы VIII всесоюзного семинара по комплексам программ математической физики). — Новосибирск, 1984.

14. Tad а Н„ Paris Р. С., I г win G. R. The stress analysis of cracks handbook. Del research corporation, 1973.

15. Bowie O. L., Freese С. E. Central crack in plane orithotropic rectangular sheet. — International Journal of fracture Mechanics, 1972, vol. 8, N 1.

Рукопись поступила 19/JX 1986