НАУЧНО-ТЕХНИЧЕСКИИ ВЕСТНИК ИНФОРМАЦИОННЫХ ТЕХНОЛОГИИ, МЕХАНИКИ И ОПТИКИ июль-август 2020 Том 20 № 4 ISSN 2226-1494 http://ntv.itmo.ru/
SCIENTIFIC AND TECHNICAL JOURNAL OF INFORMATION TECHNOLOGIES, MECHANICS AND OPTICS July-August 2020 Vol. 20 No 4 ISSN 2226-1494 http://ntv.itmo.ru/en/
HHIIIDPMAPDHHhlX ТЕХНОЛОГИЙ, МЕХАНИКИ И ОПТИКИ
УДК 004.942 doi: 10.17586/2226-1494-2020-20-4-611-616
БЕССЕТОЧНОЕ МОДЕЛИРОВАНИЕ УПРУГИХ ДЕФОРМАЦИЙ ПОЛИМЕРНЫХ КОМПОЗИТНЫХ МАТЕРИАЛОВ ПРИ ИХ СТАТИЧЕСКОМ НАГРУЖЕНИИ А.В. Сизая, И.В. Цивильский
Казанский национальный исследовательский технический университет им. А.Н. Туполева (КАИ), Казань, 420111, Российская Федерация
Адрес для переписки: angel.sizaya@gmail.com Информация о статье
Поступила в редакцию 11.05.20, принята к печати 25.06.20 Язык статьи — русский
Ссылка для цитирования: Сизая А.В., Цивильский И.В. Бессеточное моделирование упругих деформаций полимерных композитных материалов при их статическом нагружении // Научно-технический вестник информационных технологий, механики и оптики. 2020. Т. 20. № 4. С. 611-616. doi: 10.17586/2226-1494-2020-20-4-611-616
Аннотация
Предмет исследования. Композиты — уникальные материалы, сочетающие в себе легкость и прочность, что делает их востребованными на рынке материалов для промышленного производства. Математическое моделирование данных материалов необходимо для прогнозирования их деформаций под различными нагрузками. Предложен бессеточный метод математического моделирования анизотропных композитных материалов, основанный на представлении их в виде метачастиц с эластичными связями. Особенностью предложенного метода является применение единых уравнений для мультимасштабного моделирования: на микромасштабе (одно волокно и элементарный объем наполнителя) и макромасштабе (сечение целого композитного элемента). Метод. Методика расчета заключается в применении интегрирования по времени методом Верле для определения смещений частиц с последующим решением системы уравнений для эластичных связей, используя вычисленные ранее значения масс и жесткостей. Метод позволяет оценить динамику деформации композитного материала во времени под действием прикладываемой нагрузки. Код для решения уравнений упругой динамики и визуализации результатов написан на JavaScript без сторонних библиотек. Основные результаты. С использованием предложенного метода смоделированы упругие деформации упрощенной двумерной модели углепластика (на основе смеси эпоксидной и полиэфирной смол в качестве наполнителя) под действием статичной механической нагрузки на прогиб. Расчеты для композита производились на микро- и макромасштабах с разным направлением укладки волокон: углы 0 и 90°. Выбран плоский тип укладки. Результаты расчетов верифицированы в программной среде ANSYS Composite PrepPost, реализующей метод конечных элементов. Совпадение результатов расчетов методом метачастиц и МКЭ составило 89 %. Практическая значимость. Разработанный решатель может использоваться для прогнозирования деформаций композитных материалов при моделировании на микро- и макромасштабах. Ключевые слова
бессеточные методы, композитные материалы, структурная динамика, прочность, математическое моделирование
doi: 10.17586/2226-1494-2020-20-4-611-616
MESHLESS MODELING OF ELASTIC DEFORMATIONS OF POLYMERIC COMPOSITE MATERIALS UNDER STATIC LOADING
A.V. Sizaya, I.V. Tsivilskiy
Kazan National Research Technical University named after A.N. Tupolev - KAI, Kazan, 420111, Russian Federation Corresponding author: angel.sizaya@gmail.com Article info
Received 11.05.20, accepted 25.06.20 Article in Russian
For citation: Sizaya A.V., Tsivilskiy I.V. Meshless modeling of elastic deformations of polymeric composite materials under static loading. Scientific and Technical Journal of Information Technologies, Mechanics and Optics, 2020, vol. 20, no. 4, pp. 611-616 (in Russian). doi: 10.17586/2226-1494-2020-20-4-611-616
Abstract
Subject of Research. Composites are unique materials combining lightness and strength, which makes them popular in the materials industry. Mathematical modeling of these materials is necessary to predict their behavior under
specific loads. The paper proposes a meshless method for mathematical modeling of anisotropic composite materials based on their representation as a set of meta-particles with spring-like bonds. Special feature of the presented work is the applicability of the same equations for the multi-scale modeling: on micro-scale (single fiber and an elementary volume of the compound) and macro-scale (cross-section of the whole composite part). Method. The method is based on Verlet-like time integration of particle displacements with consequent resolving of a system of elastic bonds with pre-calculated mass and stiffness. The method estimates temporal dynamics of composite deformations under applied mechanical load. The code for solution of motion equations and result visualization is written in pure JavaScript without any dependencies. Main Results. Elastic deformations of a simplified 2D model of carbon-fiber-reinforced plastic have been simulated (with a mixture of epoxy and polyester resins as a compound) under static mechanical load using the proposed method. Calculations on the micro-and macro-scales with different directions of fiber layering: at the angles of 0° and 90°. Flat layering type is selected. The results are verified via the ANSYS Composite PrepPost solver under equivalent conditions. The coincidence of the calculation results by the meta-particle method of the developed solver and the finite element method is 89 %. Practical Relevance. The results obtained can be used in the development of new types of composite materials at the modeling and forecasting stage, and can also allow simulations of composites taking into account micro-processes for exclusion of pore formation and other defects that cannot be resolved using macro-scale modeling. Keywords
meshless methods, composite materials, structural dynamics, strength, mathematical modeling
Введение
Композитные материалы (КМ) уникальным образом сочетают в себе свойства прочности и легкости, что делает их востребованными в промышленном производстве. Для того чтобы прогнозировать возможные нагрузки на элементы из композитов и соответственно учитывать эти данные при конструировании деталей, необходимо моделировать процессы, связанные с данными элементами, с помощью компьютера.
Композитный материал по своей структуре представляет однородный наполнитель, армированный прочным волокном. Изменяя код укладки и материалы компонентов, можно получить требуемые физико-механические свойства.
Традиционно структурная динамика КМ в большинстве существующих решателях моделируется классическим методом конечных элементов (МКЭ) [1-3]. В связи с анизотропностью композитов, уравнения для построения конечно-элементной сетки, основного элемента МКЭ [4-6], становятся сильно нелинейными, что затрудняет процессы моделирования и прогнозирования. Рассматривая композиты в виде классических симплексов МКЭ, даже уравнения на скалярные величины требуют дополнительной вычислительной нагрузки, в силу необходимости учета анизотропности теплопроводящих свойств и плотности материала [7].
На данный момент виртуальные испытания проводятся повсеместно, и уже существует несколько коммерческих решений для симуляций физических процессов, содержащих модули для моделирования слоистых композитов: COMSOL Multiphysics, ANSYS Composite PrepPost (ACP). Данные программные обеспечения основываются на классическом МКЭ. Данный метод включает в себя построение конечно-элементной сетки, для чего необходимо решить систему уравнений упругости, в частности закон Гука:
-V(C:8) = f,
где C — тензор жесткости; г — тензор деформаций; V — оператор набла (векторный дифференциальный оператор); f — внешняя сила.
Несмотря на преимущества КМ, они, в зависимости от геометрии армирующих компонентов и их взаимного расположения, могут быть не только изотропными, но и квазиизотропными, т. е. анизотропными в микрообъемах, но изотропными в объемах всего изделия, а также сильно анизотропными, что влечет за собой нелинейность уравнений, описывающих их механику. Матрица жесткости для композитов содержит 9 х 9 элементов, и при каждом изменении кода укладки, т. е. модификации ориентации волокон в наполнителе, ее необходимо определять повторно, а также проводить натурные эксперименты для верификации полученных значений характеристик упругих свойств. Этим объясняется существенное усложнение моделирования и прогнозирования физико-механических свойств КМ, а следовательно, и внедрение новых типов композитов.
Цель работы — разработать и верифицировать упрощенный и более естественный метод моделирования КМ для сокращения времени этапа симуляций, расширения их возможностей, следовательно, ускорения подготовительных процессов в производстве деталей из композитов.
Для данного исследования поставлены следующие задачи.
— Разработать бессеточный метод связанных дискретных метачастиц. Данный метод сосредоточен не на решении сеточных уравнений, а на расчете динамики частиц, связи между которыми могут удаляться при превышении предела прочности, что симулирует процесс разрушения.
— Сравнить полученные результаты с ANSYS Composites PrepPost1.
Выполнение необходимого решения разделено на четыре этапа:
Моделирование системы связанных метачастиц
Бессеточный метод связанных дискретных метачастиц, в отличие от классического МКЭ [8-11], используемого в большинстве существующих программных продуктах для моделирования композитных материалов, основывается на представлении каждого волокна композита и элемента наполнителя в виде метачастиц с эластичными связями.
Актуальные методы исследования процессов структурной динамики композитных материалов подробно рассмотрены в работах [12, 13], в которых показано, что макроскопические модели не всегда позволяют с достаточной точностью рассчитать количественные характеристики процессов структурной механики, определяющие поведение композитного материала под нагрузкой. С другой стороны, моделирование этих процессов на микромасштабе МКЭ требует серьезной вычислительной мощности оборудования и временных затрат.
Поскольку метод метачастиц не включает построение конечно-элементной сетки и основывается на реальной волокнистой композитной структуре, при больших деформациях и нагрузках он может быть более выгодным для моделирования на разных масштабах.
Благодаря автоматическому анализу жесткости композита на макроуровне на основе известных данных об упругих свойствах отдельно взятого волокна, при использовании предлагаемого решения ускоряется процесс подготовки модели. Эта особенность позволит автоматически определять упругие свойства целого ряда композитов, состоящих из заданных материалов волокна и наполнителя, но отличающихся способом укладки (в программах-аналогах для составления базы упругих свойств материалов требуются отдельные натурные испытания для каждого типа укладки).
Таким образом, авиапромышленное предприятие на этапе моделирования и прогнозирования сократит время расчета серии композитов с разнородной укладкой, но из одинаковых компонентов, пропорционально числу вариантов укладки (расчет свойств композитов с 10 различными укладками будет в 10 раз быстрее аналогов МКЭ). Изменение в структуре никак существенно не повлияет на ход расчета: точность будет не ниже, чем при сопоставительных расчетах на коммерческих аналогах, а стабильность при этом повысится. При усложнении конструкции разница в производительности метода увеличится в пользу предлагаемого бессеточного решения.
Основной алгоритм расчета условно делится на две части: расчет динамики положений частиц и динамики связей между ними.
К первой части относится вычисление таких базовых параметров, как масса частицы
где р — плотность материала; V — объем элементарной ячейки; N — число частиц в ячейке, а также жесткость связи
Д(1-у)
где Е — модуль Юнга для данного материала; V — коэффициент Пуассона.
Для расчета смещений используется интегрирование по времени методом Верле, имеющим второй порядок точности по времени. Для этого смещения г-ой частицы в моменты времени г - т и г + т раскладываются в ряд Тейлора:
йи, 1 сРщ
ф + т) = М,(0 + — X + -— т2 + о(т3) + оСт4) + ... (1) яг 2 аг
йи, 1 (12и,
иЦ-х) = «,(0 + 0(ТЗ) + •■■ (2)
где о — члены высших порядков относительно шага по времени.
Складывая уравнения (1) и (2), заменяя - на
— (из второго закона Ньютона) и перенеся и(г - т) тх
вправо, получаем выражение (3) для расчета смещения
на следующем шаге по времени через смещения на предыдущем шаге по времени и текущее, используя значение силы, действующей на частицу, содержащую гравитационную составляющую:
и (г + т) = 2и(г) - и (г - т) + —т2, (3)
где ¥ — сила, действующая на частицу; т — шаг по времени.
Ко второй части расчета относится решение уравнений для связей между частицами:
и(р2) = и(р2) - <1 уэт
I & и(р1) = и(р1)-й—ат
где d — единичный вектор вдоль связи; I — исходная длина недеформированной связи; ^ — жесткость связи; Р1 — первая частица; р2 — вторая частица.
Из данных уравнений вычисляется удлинение связи:
и12 = и(Р2) - и(Р1).
Данное удлинение характеризует жесткость связей волокна, наполнителя и связей между волокном и наполнителем (рис. 1).
В зависимости от материала, сопоставленного паре частиц, определяется жесткость связи между ними (либо материал наполнителя, либо материал волокна), т. е. введены индивидуальные жесткости.
Также существует несколько типов связей (рис. 2). Данные эластичные связи зависят от разных условий: температуры, давления, плотности и т. д. Эти связи
Исходная Сдвиг каждой частицы Связь
связь на половину разрешена
от общего перемещения вдоль связи
Рис. 1. Блок-схема для визуализации изменения симплекса сетки с указанием направляющих векторов
Рис. 2. Типы связей между частицами: крест-накрест (а); по диагонали (б); без внутренних связей (в)
определяют, как далеко частицы могут быть друг от друга в рамках одной системы.
Наиболее удачными связями для расчетов являются связи крест-накрест (рис. 2, а). Это обусловлено тем, что фигура считается жесткой, т. е. сохраняющей свою форму при деформациях, только если она состоит из треугольников. Фигуры по типу прямоугольника имеют тенденцию складываться при боковом сдвиге, поэтому квадраты хуже подходят для данного моделирования.
Результаты
Для моделирования динамики участка композитного материала проведены расчеты композита на прогиб в сечении слоя под статической механической нагрузкой в 10 Н, направленной вниз по оси У. Боковые стороны были жестко зафиксированы.
Композит рассматривался на разных масштабах: — микроуровень — отдельное волокно в элементарном объеме наполнителя;
— макроуровень — группа волокон, равномерно распределенная в наполнителе с разным направлением укладки волокон: с углами 0 и 90°. Толщина волокон — 0,7 мм, расстояние между ними — 15 мм.
Выбран плоский тип укладки. Выбранные материалы: углеродное волокно и смесь эпоксидной и полиэфирной смол.
Количество ячеек в моделях исследуемого композитного материала: микроуровень — 81, макроуровень — 100. Тип соединения точек в ячейках: крест-накрест.
Свойства данных материалов, в частности, модуль Юнга и коэффициент Пуассона, могут четко продемонстрировать основные характеристики элементов композита в отдельности и анизотропность композитного материала в целом (таблица).
Так как волокно имеет четкую ориентацию, его характеристики, в зависимости от рассматриваемой оси или плоскости, имеют разные значения. Наполнитель в свою очередь — однородная среда, характеристики которой не зависят от направления, следовательно постоянны во всем объеме.
Как на микроуровне (рис. 3), так и на макроуровне (рис. 4) можно отметить, что при ф = 0° прогиб больше, чем при ф = 90°. Это объясняется тем, что жесткость волокна во много раз превышает жесткость наполнителя (разница порядка 1000 МПа). Именно поэтому
►к** ►к ►!<►!<
►их* ***
Рис. 3. Моделирование прогиба в сечении композита бессеточным методом метачастиц с эластичными связями на микроуровне под разными углами: ф = 0° (а); ф = 90° (б), где ф — угол, под которым направлены волокна. Толщина волокон 0,7 мм. (Оранжевый цвет — наполнитель, синий — волокно, голубой — жестко зафиксированные точки)
Таблица. Свойства материалов
Волокно Наполнитель
модуль Юнга коэффициент Пуассона модуль Юнга коэффициент Пуассона
Ех = 11 500 МПа vXУ = 0,28 Е = 85 МПа V = 0,3
Еу = 6 430 МПа vXZ = 0,28
Е2 = 6 430 МПа vУZ = 0,34
Ех - модуль Юнга вдоль оси X; Еу — модуль Юнга вдоль оси У; Е2 — модуль Юнга вдоль оси 2; X — абсцисса; У — ордината; Z — аппликата.
Рис. 4. Моделирование прогиба в сечении композита бессеточным методом метачастиц с эластичными связями на макроуровне под разными углами: ф = 0° (а); ф = 90° (б), где ф — угол, под которым направлены волокна. Толщина волокон
Рис. 5. Моделирование прогиба в плоскости композита на макроуровне под углом ф = 0°. Размеры домена 7 х 7 см. Расчеты в ANSYS Composites PrepPost (а) и методом метачастиц на эквивалентной сетке (б)
прогиб в двух перпендикулярных направлениях так отличается.
Для верификации модели выполнены аналогичные расчеты МКЭ с помощью АСР [14] с эквивалентными нагрузками (рис. 5) — прогиб в плоскости композита на макроуровне под статической механической нагрузкой в 10 Н, направленной вниз по оси У. Боковые стороны были также жестко зафиксированы.
Выбран плоский тип укладки.
Размеры домена 7 х 7 см.
Вычислено совпадение результатов расчетов методом метачастиц и МКЭ:
/ и&»Р \\
\\ Wbottom & а ( Utop \\
\\ubottom &м
•100 % = 89 %,
где и0р — среднее смещение верхних граней; иЬоиот — смещение нижних граней; иге1 — относительное смещение; индекс «аср» указывает на то, что смещения получены в АСР; индекс «мч» указывает на то, что смещения получены методом частиц.
Заключение
Разработана математическая модель, описывающая динамику композитного материала на микро- и макроуровнях под статической нагрузкой с использованием бессеточного метода связных частиц.
Проведена симуляция процессов статической нагрузки для композита классическим методом конечных элементов в ANSYS Composites PrepPost. Получено совпадение результатов, равное 89 %.
Выявлены: преимущества метода метачастиц:
— расчет в масштабе микроуровня;
— расчет в сечении композита, а также его недостатки:
— при увеличении деформаций могут появляться вывернутые ячейки с отрицательным объемом.
В дальнейшем предстоит:
References
Авторы
Сизая Анжелика Владимировна — студент, Казанский национальный исследовательский технический университет им. А.Н. Туполева (КАИ), Казань, 420111, Российская Федерация, ORCID ID: 0000-0002-2027-3394, angel.sizaya@gmail.com
Цивильский Илья Владимирович — кандидат технических наук, доцент, Казанский национальный исследовательский технический университет им. А.Н. Туполева (КАИ), Казань, 420111, Российская Федерация, Scopus ID: 57008301500, 25722440800, ORCID ID: 0000-0002-5466-8992, ivtsivilskiy@kai.ru
Angelica V. Sizaya — Student, Kazan National Research Technical University named after A.N. Tupolev — KAI, Kazan, 420111, Russian Federation, ORCID ID: 0000-0002-2027-3394, angel.sizaya@gmail.com
Ilya V. Tsivilskiy — PhD, Associate Professor, Kazan National Research Technical University named after A.N. Tupolev — KAI, Kazan, 420111, Russian Federation, Scopus ID: 57008301500, 25722440800, ORCID ID: 0000-0002-5466-8992, ivtsivilskiy@kai.ru