Full Text

Введение Развитие современных технологий позволяет создавать математические модели, реально отражающие картину распределения напряженного состояния пространственных конструкций. Особое внимание уделяется изучению влияния структурных особенностей материалов и конфигурации на напряженное состояние конструкций. Разработке математического моделирования и решению задач физически нелинейного деформирования изотропных тел с концентраторами напряжений посвящены исследования многих авторов, обзор которых представлен в [1]. В [2] излагаются основные положения (постулаты) механики сплошной среды. Наряду с классическими рассматриваются сравнительно новые модели: композит и модели, учитывающие связанность механических полей. В [3] рассмотрены алгоритм и математическое моделирование решения задачи с учетом физической нелинейности тел на основе теории малых упругопластических деформаций. Отмечается, что процесс решения задачи значительно ускоряется при использовании деформационной теории по сравнению с теорией течения. В [4] исследуются вопросы автоматизации построения моделей и визуализации результатов численного моделирования. Описание анизотропии механических свойств трансверсально-изотропных тел осуществляется на основе структурно-феноменологической модели, которая позволяет представить исходный материал в виде комплекса двух совместно работающих изотропных материалов: основного материала, рассматриваемого с позиций механики сплошной среды, и материала волокон, ориентированных вдоль направления анизотропии исходного материала. При этом предполагается, что волокна воспринимают лишь осевые усилия растяжения-сжатия и деформируются совместно с основным материалом. Рассматривается упругопластическая среда, которая представляет собой неоднородный сплошной материал, состоящий из двух компонентов: армирующих элементов и матрицы (или связующего), которая обеспечивает совместную работу армирующих элементов. В волокнистых материалах деформация упругопластической матрицы обеспечивает нагружение волокон высокой прочности. Известно, что волокнистый материал обладает свойствами трансверсально-изотропной среды. В связи с этим для решения задачи физически нелинейного деформирования волокнистых композитов применяется теория малых упругопластических деформаций для трансверсально-изотропной среды, предложенная проф. Б.Е. Победря [5]. В работе отмечается, что при рассмотрении армированного композита, жесткость армирующих элементов которого существенно превышает жесткость связующего, появляется возможность использования упрощенной деформационной теории пластичности. Упрощенная теория позволяет применить теорию малых упругопластических деформаций для решения конкретных прикладных задач. Суть упрощения заключается в предположении того, что при простом растяжении композита в направлении оси трансверсальной изотропии и направлении, перпендикулярном к ней, пластических деформаций не возникает. Вследствие этого интенсивность напряжений и деформаций определяется отдельно как по главной оси трансверсальной изотропии, так и по перпендикулярно расположенной плоскости. Применение упрощенной теории основывается на том, что у рассматриваемого армированного композита жесткость армирующих элементов существенно превышает жесткость связующего. Поэтому интенсивность напряжений и деформаций определяется отдельно как по главной оси трансверсальной изотропии, так и по перпендикулярно расположенной плоскости. В работе предложен метод решения задачи упругопластического деформирования трансверсально-изотропных многосвязных тел. Алгоритм решения задачи включает: 1) определение эффективных параметров трансверсально-изотропной среды; 2) построение конечно-элементной сетки конфигурации тела; 3) определение локально-минимального значения ширины ленты ненулевых коэффициентов системы уравнений; 4) построение коэффициентов матрицы жесткости и компоненты вектора узловых нагрузок уравнения состояния отдельного конечного элемента согласно теории малых упругопластических деформаций для трансверсально-изотропной среды; 5) формирование разрешающей системы симметрично-ленточных уравнений; 6) решение системы симметрично-ленточной системы уравнений модифицированным методом квадратных корней; 7) вычисление упругопластического напряженно-деформированного состояния трансверсально-изотропного тела методом упругих решений А.А. Ильюшина. Для каждого этапа решения задачи разработаны вычислительные алгоритмы, позволяющие уменьшить количество вычислительных операций за счет модификации существующих методов решения и учета структуры коэффициентов матрицы. В качестве примера приведено решение задачи о деформировании трансверсально-изотропного тела в форме прямоугольника с круговым вырезом в центре. 1. Постановка задачи и метод решения Исследуется упругопластическая среда из неоднородного сплошного материала, состоящего из двух компонент: армирующих элементов и матрицы (или связующего), которая обеспечивает совместную работу армирующих элементов. Общая постановка краевой задачи теории упругости для анизотропных тел включает: - уравнение равновесия (1) - обобщенный закон Гука (2) - соотношения Коши и краевые условия (3) (4) где - компоненты вектора перемещений; Xi, - объемные и поверхностные силы; S1, S2 - слагаемые части поверхности S объема V; nj - внешняя нормаль к поверхности S2 объема V; - тензор упругих констант. Для упрощенной теории малых упругопластических деформаций трансверсально-изотропной среды обобщенный закон Гука (2) принимает следующий вид: , (5) где (6) и - функции пластичности типа Ильюшина, значения которых в упругой зоне равны нулю ( и - коэффициенты упрочнения и пределы упругой деформации, соответственно по плоскости изотропии и оси изотропии ). В упругой области параметры σij определяются законом Гука, а в области пластических деформаций - на основе деформационной теории А.А. Ильюшина; li - упругие постоянные трансверсально-изотропных материалов; - составляющие разложения девиаторных частей трансверсально-изотропного тензора деформации. Интенсивность тензора напряжений и тензора деформаций соответственно по плоскости и по оси изотропии: причем Поскольку среда является однородной с эффективными механическими параметрами по оси изотропии и плоскости изотропии, для решения упругопластических задач используется итерационной процесс метода упругих решений А.А. Ильюшина [6]. Эффективные механические параметры трансверсально-изотропного материала связаны с модулями λi следующими соотношениями: Здесь μef и - эффективные коэффициенты Пуассона, Еef и - эффективные модули Юнга и Gef и - эффективные продольные модули сдвига соответственно по плоскости трансверсальной изотропии и оси трансверсальной изотропии. Для вычисления эффективных механических параметров волокнистых материалов в работе используются соотношения, позволяющие учесть внутреннюю структуру материала для расчета периодически неоднородных материалов на основе асимптотического метода осреднения и пригодные для любых значений свойств и объемных долей компонентов [7]: - эффективный продольный модуль Юнга: где - эффективный продольный модуль сдвига: где - эффективный продольный коэффициент Пуассона: В общем случае, посредством представления отношений между тензором напряжений и тензором деформаций в виде функции , соотношения Коши и вектора смещений каждой частицы в системе координат Оx1x2x3 как можно представить нелинейную связь между тензором напряжений и вектором смещений [5]: В этом случае уравнение равновесия (1) определяет систему из трех уравнений в частных производных относительно трех компонент вектора смещений. Для этой системы уравнений можно поставить три типа граничных условий: в перемещениях (3), напряжениях (4) или смешанный тип. Таким образом, процесс деформирования твердого тела, находящегося в равновесии под действием внешних сил, можно свести к определению вектора смещения . На основе решения краевых задач можно определить компоненты вектора смещений. По известным значениям компонент вектора смещений можно определить компоненты тензора деформаций и тензора напряжений. 2. Дискретизации многосвязной области Конечно-элементная сетка многосвязной области формируется посредством «сшивания» (объединения) канонических подобластей. Под канонической подразумевается область, для которой существует алгоритм построения конечно-элементной сетки. В качестве конечных элементов выбираются четырехугольник (в случае двумерного тела) и четырехугольная призма (для трехмерного тела), так как заполнение реальной области этими элементами является весьма эффективным. Вычислительный алгоритм построения конечно-элементной сетки многосвязной области состоит из последовательности следующих этапов: 1) формирование конечно-элементной сетки канонических подобластей; 2) «сшивание» (объединение) подобластей; 3) определение начального фронта; 4) упорядочение номеров узлов на основе фронтального метода; 5) минимизация ширины ленты системы уравнений. На начальном этапе решения задачи формируется библиотека конечно-элементных сеток канонических областей. Объединение областей выполняется на основе критерия совпадения граничных узлов посредством установления простой иерархии объемов, поверхностей, линий и точек. Формирование элементов множества начальных фронтов осуществляется посредством определения номеров вершин и граничных узлов конечно-элементного представления многосвязной области. Перенумерация узлов и конечных элементов проводится на основе фронтального метода [8], с учетом того, что в качестве начального фронта используются вершины и грани конструкции. Известно, что ширина ленты ненулевых коэффициентов системы разрешающих уравнений МКЭ зависит от упорядочения узлов и конечных элементов конечно-элементной сетки. Нумерация, при которой ширина ленты имеет наименьшее значение, определяется путем проведения вычислительного эксперимента. В качестве параметра варьирования используются элементы множества начальных фронтов. Разработанный алгоритм построения конечно-элементной сетки позволяет посредством упорядочения нумерации узлов определить такую последовательность, при которой ширина ленты ненулевых коэффициентов разрешающей системы уравнений МКЭ будет локально-минимальной [9]. 3. Построение разрешающей системы симметрично-ленточных уравнений Наиболее существенным и в значительной степени определяющим качество схемы расчета МКЭ является построение коэффициентов матрицы жесткости и формирование системы разрешающих уравнений. Каждый шестигранный элемент из семейства конечных элементов, на которые разбивается тело, можно отобразить в трехмерный элемент правильной формы, например в куб. Для линейного конечного элемента с 8 узлами функции формы Ni в локальной системе координат можно представить в виде заданных выражений, которые удовлетворяют всем необходимым критериям [10]: где - координаты i-го узла, С этой целью глобальные координаты, в данном случае декартовые, связываются с локальными координатами соотношением Пользуясь наличием такой связи, можно выразить производные от функции формы, заданные в глобальных координатах, через локальные координаты где - матрица Якоби. Здесь где точки с координатами xi, yi, zi по определению свойств функции формы совпадают с соответствующими точками границы элемента ( n - количество узлов в элементе). Отсюда можем записать: (9) Далее, учитывая взаимно-однозначное соответствие между локальной и глобальной системами координат, можно записать компоненты якобиевой матрицы в следующем виде: (10) Обозначим компоненты обратной матрицы следующим образом: (11) Используя полученные отношения (9)-(11), вектор деформации конечного элемента е можно записать в локальной системе координат: где - компоненты обратной матрицы 0 - нулевой вектор размерности 3; - нулевая матрица размерности 3´8. Вектор напряжений, связанный с вектором деформации законом Гука, также можно представить в виде Теперь, сделав соответствующую подстановку, можно записать выражение для вычисления матрицы жесткости: ´ ´ с учетом того, что Для вычисления интеграла от поверхностных нагрузок необходимо каждую поверхность конечного элемента рассматривать в отдельности. В зависимости от того, на какой поверхности задана нагрузка, можно записать шесть различных значений этого интеграла. Предположим, что на поверхности приложена равномерно распределенная нагрузка и на ней расположены четыре первых узла конечного элемента. Тогда вектор узловых сил можно представить в форме где Коэффициенты матрицы жесткости изопараметрических конечных элементов и компоненты вектора узловых нагрузок вычисляются посредством гауссовых квадратур, поскольку обеспечивается наибольшая точность при заданном числе точек интегрирования. Трижды применяя формулу Гаусса для интегрирования функции одной переменной, можно записать: где t - число точек интегрирования; - весовые коэффициенты. Такой подход особенно эффективен в трехмерном случае, так как позволяет примерно вдвое сократить число точек интегрирования по сравнению со стандартной - точечной формулой. Для интегрирования выражений и достаточно квадратуры второго порядка. Координаты узлов и весовые коэффициенты для квадратуры Гаусса приведены в работе [11]. Построенные коэффициенты жесткости для всех конечных элементов используются при формировании разрешающей системы уравнений. Поскольку в одном узле сходится обычно несколько элементов, сборка фактически состоит в том, чтобы для каждого узла в каждом из направлений просуммировать соответствующие коэффициенты матриц жесткости смежных элементов и расположить эту сумму в нужное место глобальной матрицы жесткости. Количественная оценка напряженно-деформированного состояния пространственных тел требует разбивки области, занятой телом, на большое число конечных элементов, что приводит к построению системы алгебраических уравнений высокого порядка и сопряжено с определенными трудностями при реализации на вычислительных системах. С этой целью в работе использован способ построчной подготовки данных для каждого узла в отдельности, что обеспечивает построение ленточной системы алгебраических уравнений высокого порядка с учетом симметричности ее коэффициентов. 4. Решения симметрично-ленточной системы уравнений высокого порядка Использование МКЭ приводит к разрешающей системе линейных алгебраических уравнений высокого порядка. Для решения системы уравнений применяется метод квадратных корней, модифицированный для симметрично-ленточной структуры матрицы коэффициентов [12]. Поскольку каждый конечный элемент связан с ограниченным числом других конечных элементов, матрица системы уравнений всегда получается редко заполненной. Расположение коэффициентов в матрице тесно связано с порядком нумерации узлов в теле. При случайном порядке нумерации ненулевые компоненты распределяются в матрице также случайным образом, а это, в свою очередь, ведет к увеличению времени расчета задачи. Поэтому необходимо подобрать такой порядок нумерации, при котором ненулевые элементы были бы сгруппированы в окрестностях главной диагонали матрицы, т.е. образовывали полосу (ленту). Процедура упорядочения номеров узлов в дискретной модели позволяет минимизировать ширину ленты ненулевых коэффициентов разрешающей системы уравнений. Поскольку коэффициенты матрицы такой системы уравнений и коэффициенты матрицы жесткости конечных элементов симметричны, на этапе решения системы уравнений целесообразно использовать только диагональные элементы и элементы, располагающиеся ниже ее, т.е. нижнюю треугольную матрицу. При этом коэффициенты последней матрицы также имеют ленточную структуру. В преобразованиях метода квадратных корней в основном применяется операция умножения матрицы на вектор, в связи с этим разработан алгоритм умножения матрицы на вектор для случая, когда заданы только коэффициенты нижней ленты треугольной матрицы. Если расположить эти коэффициенты построчно, формируется прямоугольная матрица Sij размерами , где n - порядок системы уравнений, l - половинная ширина ленты ненулевых коэффициентов, включая диагональные элементы. Причем диагональные элементы исходной матрицы располагаются на последнем, l-м столбце матрицы Sij. Для иллюстрации преобразований предположим, что n = 9, l = 4. В этом случае нижняя треугольная и прямоугольная Sij матрицы имеют нижеприведенные представления: , . Для формирования процесса умножения матрицы Sij на вектор xj в случае, когда за основу берутся преобразованные коэффициенты матрицы, расположенные по диагонали и ниже, разработано и используется следующее соотношение: , где , Приведенное соотношение позволяет в алгоритме метода квадратных корней, посредством использования соответствующих диагональных коэффициентов и коэффициентов нижней треугольной матрицы, обойтись без нулевых коэффициентов, расположенных вне ленты ненулевых элементов и коэффициентов верхней треугольной матрицы. На основе полученных решений системы разрешающих уравнений, которые соответствуют узловым перемещениям, по соотношениям (5), (6) вычисляются параметры напряженного состояния. 5. Тестовый пример Приводятся результаты численного моделирования решения задач упругопластического деформирования элементов конструкций из волокнистых материалов. Рассматривается трехмерная упругопластическая задача об одноосном растяжении (Рzz = 850 МПа) по направлению волокон прямоугольной пластины высотой 10 мм, с шириной 5 мм и толщиной 1 мм. Изолированное отверстие радиусом мм. Приводятся результаты численного моделирования решения задач физически нелинейного деформирования элементов конструкций из волокнистых материалов. В качестве материала матрицы используется алюминиевый сплав Д16 (бороалюминий) с параметрами: модуль упругости Е = 7,1·104 МПа, коэффициент Пуассона μ = 0,32, коэффициент упрочнения и предел упругости σs = 2,13·102 МПа. Для борного волокна Е′ = = 39,7·104 МПа, μ′ = 0,21, предел прочности при растяжении = = 2,5·103 МПа. Пластина разбивается на призматические конечные элементы. Поскольку в окрестности отверстия наблюдаются повышенные напряжения, она разбивается на более мелкие конечные элементы. На рисунке, а приведено поле распределения значений интенсивности деформаций pu по плоскости изотропии в окрестности деформированного отверстия (изображено деформированное состояние конфигурации конструкции). Зона пластических деформаций формируется в окрестности верхней и нижней части отверстия. Это связано с тем, что под действием растягивающей нагрузки центральная часть конструкции сжимается по оси вследствие этого конфигурация окрестности отверстия деформируется: растягивается по оси и сжимается по . При сжатии по оси сжимается матрица композита, а поскольку ее механические параметры намного меньше, чем у волокна, именно в этих местах и формируется зоны пластических деформаций матрицы волокнистого композита [13]. а б в Рис. Распределение интенсивности деформаций На рисунке, б приведено распределение значений интенсивности деформаций по главной оси трансверсальной изотропии. Повышенные значения упругих деформаций формируются по бокам отверстия, однако в окрестности точек пересечения горизонтального диаметра с контуром отверстия значения минимальны (тон закраски соответствует удвоенным значениям деформаций по плоскости изотропии, приведенным на рисунке, в). Заключение На основе метода конечных элементов и упрощенной теории малых упругопластических деформаций трансверсально-изотропных сред разработаны численная модель и алгоритм решения трехмерных задач упругопластического деформирования трансверсально-изотропных тел.

About the authors

A. M Polatov

National University of Uzbekistan named after Mirzo Ulugbek

References

Statistics

Views

Abstract - 37

PDF (Russian) - 12

Refbacks

• There are currently no refbacks.