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

Известия Тульского государственного университета Естественные науки. 2010. Вып. 2. С. 86-96

Механика :

УДК 539.3

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

В.А. Мазин, В.Л. Михайлова, Л.Г. Сухомлинов

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

Ключевые слова: плоская задача теории упругости, вариационно-разностная процедура, напряжения около отверстий.

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

Как известно [1], плоская задача теории упругости формулируется для области (будем обозначать ее как 5), лежащей в плоскости Оху прямоугольной декартовой системы координат Охуг. Упомянутая область 5 представляет собой либо срединную поверхность тонкой пластинки (случай плоского напряженного состояния), либо поперечное сечение протяженного вдоль оси Ог цилиндра (случай плоского деформированного состояния). Указанные пластинка и цилиндр находятся под действием поверхностных и объемных сил, интенсивности которых будем обозначать как / и /, а их проекции на оси Оху — соответственно как дх, ду и fxх, /у. Для компонент вектора перемещений вводим обозначения их и иу, а для компонент тензоров напряжений и деформаций, входящих в формулировку плоской задачи, — обозначения

^xx, @уу, @ху и ^xx, £уу, £ху.

В случае изотропного материала вытекающая из закона Гука связь между обозначенными напряжениями и деформациями может быть представлена в виде

Охх — А1 £хх + А2^уу,

Оуу — А2^хх + А1^уу, (1)

@ху — 2О^ху ■

Коэффициенты О, А1, Л2 линейных зависимостей (1) выражаются через модуль Юнга Е и коэффициент Пуассона V материала исследуемого тела согласно следующей схеме:

Е .

1) О 2(1+v) ’

2) в случае плоского напряженного состояния

Е

А1 — -----2, А2 — vАl;

1 — V2

3) в случае плоского деформированного состояния

А1 — 2О + А, А2 — А,

где

vE А —

(1 + V)(1 - 2v)

Деформации £xx, eyy, exy, входящие в выражения (1), связаны с перемещениями их и uy в каждой точке исследуемого тела соотношениями Коши

вида / ч

dux duy 1 ( dux , duy \

£xx = -gx- eyy = ly' exy = 2 lay + -ax)- (2)

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

Ux\ru = Ux' иУ lru = Uy ■ (3)

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

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

(^хх^хх + Оуу¡еуу + 2&хуЗ^ху^$5 ---

— II(^хЗих + ¡уЗиу)d5 + J (дхбих + Яу¡иу)dГ■ (4)

5 Гд

Здесь с!5 — элемент площади области 5, а dГ — элемент длины границы

Г. Интеграл по Г берется только по части Гд границы Г, учитывая, что на участке Гад в соответствии с граничными условиями (3) должно быть

5их\ги — 0, 5иу Г — 0. (5)

Входящие в уравнение (4) вариации деформаций с использованием равенств (2) могут быть выражены через вариации перемещений так, что

д5их , д5иу 1 / ддих , ддиу\

хх — -х, ¡еуу — беху («)

Далее принимаем в качестве объекта исследования заданную в первой четверти плоскости Оху прямоугольную область 5 с размерами Ьх, Ьу. На рис. 1 представлена схема разбиения данной области на М х N прямоугольных элементов 5(^ (г — 1, 2,...,М; ] — 1, 2). Для определенности считаем, что участок Гад границы Г рассматриваемой области 5 включает левую и нижнюю стороны прямоугольника 5, а участок Гд — правую и верхнюю его стороны.

С учетом этого вариационное уравнение (4) для рассматриваемой прямоугольной области может быть представлено в виде М N гг

£ £// (ах^:„ + ат^ ),18 —

г=1 7=5(г,з)

М N пп М с х^+1

х х у у х

ЕЕ (¡хЗих + ¡у ¡иу )d5 + ^ / [дх(х)5их + Яу (х)6щ ]dx + (7)

г=1 ]=1^гЛ г=1-> х

8(г,з)

N Г уз+1

+ ^ [Ях(у)Зих + Яу (y)¡Uy]dy■

._1 ^у1

О XI Х2 ■ ■ • -V, Л',-1 . . . Д.'л/ Л.'д/+1 X

Рис. 1. Прямоугольная область с сеткой прямоугольных элементов

Предполагая, что величины ^ — х%+1 — хг и 1у^ — у7+1 — у7 (г — 1, 2, ■ ■ ■

..., М; ] — 1, 2, ■ ■ ■, N) достаточно малы, вычисление интегралов в уравнении (7) приближенно выполним по значениям подынтегральных функций в серединах соответствующих элементарных участков интегрирования. В результате получим

М N

Х'' \^[^(г,з)¡Ф’7) + п(г’Ляе(г’3) + 2Тг(г’3)¡?(г’3)]1(г)1(з) —

хх хх уу уу ху ху х у

І=1 3 = 1

М N

М

ЕЕ® ’3) биХ’3 + /У ’3 ' ’3) ] 1х^1у3 +ЕІ* (и(і'))6и.і') + »у (и(і) )Ц; ^} +

=1 3=1

=1

N

(8)

+ Х^[»х(и(з))5и.з) + »У (У(з))6иу)]1Уз)-

3=1

Здесь «волной» отмечены величины, определяемые в серединах соответствующих участков интегрирования. При этом

ии

(і)

= и + ]$/2, и(з) = Уз + 43)/2.

(9)

Вводим далее обозначения их7, и^7 для значений перемещений в узловых точках рассматриваемой сетки прямоугольных элементов с координатами хг, у7 (г — 1, 2, ■ ■ ■, М + 1; ] — 1, 2, ■ ■ ■ N + 1). Входящие в вариационное уравнение (8) деформации, относящиеся к середине элемента 5(г,7), определяем на

основе соотношений (2), вычисляя при этом соответствующие частные производные с применением центрально-разностных схем. В результате имеем

^ — о, 5(их+17 + их+17+1 — их — их^/ф, еуу) — о, 5(иу7+1 + иу+1,7+1 — иу7 — иу+1,7)/1(7),

— о, 25(их7+1 + <+17+1 — их7 — их+17 )/1(7) + (10)

+ 0, 25(иу+17 + иу+17+1 — и7 — иУ^^/^^

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

их"1 — 0, 5(их^+1 + uX+1’N+1), и^ — 0, 5(uУ’N+1 + и-+1^+1), и(х7) — 0, 5(иМ+1,7 + иМ+1,7+1), и7) — 0, 5(иМ+1,7 + иМ+1,7+1), (11)

а для середин прямоугольных элементов —

ихг’7) — 0,25(их7 + их+1’7 + игх7+1 + их+1’7+1),

и(г’7) — 0,25(иг7+иг+1’7+игу7+1+и^1^1) (12)

Для напряжений в серединах элементов в соответствии с соотношениями упругости (1) записываем

и(у,з) — \ иг,7) + \

Vхх — А1ехх "Т” А2еуу ,

Оуу7 — А2ихх7) + А^уу^, (13)

и(у,7) — 2^иг,з)

V ху — ¿^еху ■

С использованием связей (10)—(13) окончательно приходим к формулировке вариационного уравнения (8) в терминах узловых перемещений. С учетом того, что часть узловых перемещений задается граничными условиями (3), приравнивая коэффициенты при вариациях неизвестных узловых перемещений в левой и правой части указанного вариационного уравнения, получаем разрешающую систему линейных алгебраических уравнений относительно неизвестных узловых перемещений. Решаем эту систему методом Гаусса, а затем определяем параметры напряженно-деформированного состояния в средних точках элементов расчетной модели с использованием связей (10), (12), (13).

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

Отметим, что для сформулированной модели (как и в случае метода конечных элементов [3]) важным является вопрос о реакции элемента модели

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

их — ф^Д иу7 — —фф/2, их+1,7 — Ф7/2, иУ+1,7 — ф1хУ)/2, (14)

их7+1 — —ф^/2, и^1 — —фх/2,

их+1,7+1 — —ф^/2, иу+1,7+1 —

Подставляя выражения (14) в соотношения (10), убеждаемся, что жесткий поворот элемента не приводит к деформациям в его середине.

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

подставляя выражения (14) в формулы (12). Таким образом, в отличие от ме-

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

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

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

Проблему применения данной процедуры к анализу напряжений вокруг отверстий рассмотрим на примере области с центральным круговым отверстием (см. рис. 2), растягиваемой вдоль оси Ох равномерно распределенной

по ее левой и правой сторонам нагрузкой (с интенсивностью дх = ^д). При этом будем рассматривать следующие случаи: а) случай протяженной вдоль осей Ох и Оу области (а = Ь = то); б) случай бесконечно протяженной вдоль оси Ох полосы (а = то), имеющей ширину Ь, сопоставимую с радиусом Я отверстия; в) случай квадратной области (а = Ь), размер Ь которой сопоставим с радиусом Я отверстия.

У*

2 Ь г л Я.

ч

) J

>

• 2а

Рис. 2. Схема прямоугольной области с центральным круговым отверстием, равномерно растягиваемой вдоль оси Ох

С учетом симметрии схемы, представленной на рис. 2, численное моделирование проводим лишь для правой верхней четверти заданной области. При этом на участках ее границы, совпадающих с осями симметрии, формулируем условия: ux =0, qy = 0 (для точек, лежащих на оси Oy); uy = 0, Qx =0 (для точек, лежащих на оси Ox).

Вводим в рассмотрение полярные координаты г, в и соответствующие компоненты тензора напряжений аг, ад, агд. Поскольку контур отверстия свободен от нагрузок, для точек области, лежащих на этом контуре (r = R), должно быть аг = агд = 0. Другими словами для точек области при r = R напряжения аг и ад являются главными напряжениями (при этом аг = 0). Исходя из этого процесс вычисления интересующих нас напряжений аг и ад на контуре отверстия строим следующим образом.

По найденным с помощью вычислительной модели напряжениям ахх, ауу, аху в точках контура отверстия определяем соответствующие главные напряжения аі и а2, используя известные формулы (см. [1])

2аху \

ахх ауу /

а1 = ахх cos2 а + ауу sin2 а + аху sin 2а, (15)

22 а2 = ахх sin а + ауу cos а — аху sin 2а.

Здесь а — угол между осью Охи направлением действия напряжения аі. При корректной работе вычислительной процедуры определяемые таким обра-

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

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

Результаты численного моделирования с применением описанной процедуры представлены на рис. 3,4 в виде зависимостей окружного напряжения а$ на контуре отверстия от угла в (см. рис. 2). Здесь же для сравнения представлены экспериментальные данные и результаты аналитических решений, взятые из монографии [2]. При этом цифрами на рис. 3 обозначены случаи: 1 — а = Ь = то (вариационно-разностное решение); 2 — а = Ь = то (точное аналитическое решение); 3 — а = то Я/Ь = 0.5 (вариационно-разностное решение); 4 — а = то, Я/Ь = 0.5 (приближенное аналитическое решение); 5 — а = то, Я/Ь = 0.5 (эксперимент). Аналогичные цифры на рис. 4 обозначают случаи:

1 — Я/Ь = 0.7 (вариационно-разностное решение); 2 — Я/Ь = 0.7; (приближенное аналитическое решение); 3 — Я/Ь = 0.8 (вариационно-разностное решение); 4 — Я/Ь = 0.8; (приближенное аналитическое решение).

Сравнивая кривые 1 и 2 на рис. 3 можно видеть, что результаты численного моделирования для случая бесконечной области практически совпадают с известным аналитическим решением /д = 1 — 2еов2в (см. [1]). Небольшие отклонения (порядка 6 %) наблюдаются лишь в малых окрестностях углов в = 0° и в = 90°. Здесь следует пояснить формирование расчетной схемы при моделировании бесконечной области.

Оценки показывают, что в точках области, отстоящих от центра О на расстоянии г ^ 10Я напряженное состояние практически не отличается от одноосного, при котором ахх = д, ауу = аху = 0. Поэтому численное моделирование осуществлялось для прямоугольной области с размерами а = Ь = 10Я. При этом на краю области с х = 10Я ставились условия дх = д, ду = 0, а на краю с у = 10Я — условия дх = ду = 0.

Разбиение моделируемой четверти области на элементы осуществлялось следующим образом. Входящая в расчетную схему четверть контура отверстия разбивалась на п = 100 одинаковых элементарных дуг, порождая тем самым часть сетки прямоугольных элементов модели. Далее отрезки

О 15 30 45 60 75 в град

Рис. 3. Зависимость окружного напряжения ад на контуре отверстия от угла в в случаях бесконечной области и протяженной полосы

Я ^ х ^ 10Я на оси Ох и Я ^ у ^ 10Я на оси Оу разбивали на участки [Я, 1.2Я], [1.2Я, 1.4Я], [1.4Я, 2Я], [2Я, 4Я], [4Я, 6Я], [6Я, 10Я], которые в свою очередь разбивали соответственно на 30, 10, 20, 15, 10, 10 одинаковых отрезков. Через концы образованных элементарных отрезков проводили прямые, параллельные осям Ох и Оу, завершая формирование сетки расчетной модели. Отметим, что применяемый здесь способ генерации сетки модели приводит к существенной диспропорции размеров прямоугольных элементов, граничащих с отверстием и расположенных в окрестностях углов в = 0° и в = 90°. Этим, по-видимому, и объясняются отмеченные выше небольшие отклонения от точного решения в указанных зонах.

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

Подобным же образом генерировалась и сетка расчетной модели для квадратной области с отверстием. Как видно из рис. 4, результаты численного моделирования для случая Я/Ь = 0.7 (кривая 1) и случая Я/Ь = 0.8 (кривая 3) несколько отличаются от соответствующих результатов аналитического решения (кривые 2 и 4). Однако, учитывая приближенный характер данного аналитического решения, согласование результатов численного моделирования и аналитического решения здесь также можно считать хорошим.

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

m/q

20.0

0.0

-20.0

-40.0

0 20 40 60 в гРаб

Рис. 4. Зависимость окружного напряжения ag на контуре отверстия от угла в в случае квадратной области

Список литературы

1. Тимошенко С.П., Гудьер Дж. Теория упругости. М.: Наука, 1979. 560 с.

2. Савин Г.Н. Распределение напряжений около отверстий. Киев: Наукова думка, 1968. 888 с.

3. Зенкевич О. Метод конечных элементов в технике. М.: Мир, 1975. 541 с.

Мазин Василий Александрович (vasilm@yandex.ru), старший преподаватель, кафедра прикладной математики, Кубанский государственный университет.

Михайлова Виктория Львовна (tm@mami.ru), к.т.н., доцент, кафедра теоретической механики, Московский государственный технический университет «МАМИ».

Сухомлинов Лев Георгиевич (tm@mami.ru), д.т.н., профессор, кафедра теоретической механики, Московский государственный технический университет «МАМИ».

On the application of a variational finite-difference computational model to the analyses of stresses in rectangular domains with holes

V.A. Mazin, V.L. Mikhaylova, L.G. Sukhomlinov

Abstract. A variational finite-difference computational procedure for the analyses of stresses in rectangular domains with holes is presented. The results of

numerical simulations are compared with the existing analytical solutions and experimental data for stress distributions around a central circular hole in a strip and in a square plate.

Keywords : plane elastic problem, variational finite-difference procedure, stresses around holes.

Mazin Vasiliy (vasilm@yandex.ru), senior teacher, department of applied mathematics, Kuban State University.

Mikhaylova Virtoriya (tm@mami.ru), candidate of technical sciences, associate professor, department of theoretical mechanics, Moscow State Technical University «MAMI».

Sukhomlinov Lev (tm@mami.ru), doctor of technical sciences, professor, department of theoretical mechanics, Moscow State Technical University «MAMI».

Поступила 17.04.2010