Численное моделирование магнитных свойств композиционных материалов Текст научной статьи по специальности «Математика»

ПРИКЛАДНАЯ МАТЕМАТИКА И МЕТОДЫ I МАТЕМАТИЧЕСКОГО МОДЕЛИРОВАНИЯ |

УДК 539.3

Ю. И. Димитриенко, Е. С. Ничеговский

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ МАГНИТНЫХ СВОЙСТВ КОМПОЗИЦИОННЫХ МАТЕРИАЛОВ

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

E-mail: dimit@serv.bmstu.ru

Ключевые слова: композитные материалы, численное моделирование,

магнитная проницаемость, асимптотическое осреднение.

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

Задача "на ячейке периодичности". Рассмотрим композиционный материал, имеющий периодическую структуру (рис. 1), ячейка периодичности (ЯП) V которого

состоит из N фаз Va, а = 1,. . . , N. Рис. 1. Структура ячейки периодичности

Введем малый параметр к = 1/Ь ^ 1 как отношение характерного размера ЯП I к характерному размеру всего композита Ь, а также исходные декартовы координаты хк, глобальные координаты хк и локальные координаты £к. Для определенности будем полагать, что одна из фаз с номером а = N является связной областью и называется матрицей, а остальные представляют собой наполнитель (например, волокна).

Рассмотрим для такой структуры задачу магнитостатики, предполагая выполненными условия идеального контакта на поверхности раздела Еа(3 между компонентами:

дгЬ*а = 0 в К; Ь*а = р*аН*а;

И*а = дгу *а; (1)

V*а = , (ЬГ - Ь*н) Пг = 0 на ;

|Е1 = V:, пгЬги = ье,

где Ь*а — компоненты комплексной амплитуды вектора магнитной индукции в фазе Уа; Н*а — компоненты комплексной амплитуды вектора напряженности магнитного поля; V *а — комплексная амплитуда магнитного потенциала; р*а — комплексная амплитуда магнитной проницаемости а-й фазы композита; пг — компоненты вектора нормали; V*, Ь* — заданные значения магнитного потенциала и магнитного потока на внешней границе композита. Звездочка у всех величин означает, что рассматриваются амплитуды комплексных величин (например, р*а = р'а + гр"а, где г — мнимая единица), соответствующие затухающим электромагнитным колебаниям с известной частотой ш. Все компоненты векторов отнесены к прямоугольной декартовой системе координат, операция дг означает дифференцирование по г-й координате хг.

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

Следуя общей идее метода асимптотического осреднения, решение задачи (1) относительно потенциала V *а строим в виде асимптотических разложений по малому параметру к так:

Vю (хг) = V*(0) (хг) + ^ *а(1) (хг, £ £) +

+ к2УГ(2) (X) + ... в V« , (2)

где -у*а(2),... — неизвестные функции, подлежащие вычислению,

периодические по локальным координатам £к; функция -у*(0) зависит только от глобальных координат.

Подставляя разложение (2) в систему (1), используя правила дифференцирования периодических функций [1-3] и собирая члены при одинаковых степенях малого параметра к, для функций -у*а(1) получаем локальную задачу "на ячейке периодичности"

ь;;(0) = о в

7 *а(°) _ 11*атт*а(0).

Ьг - ^ Иг ,

TT*а(0) _ i* , „, * а(1).

H - ti* + V/t ,

(3)

v *а(1) - v * N (1)

b*а(0) _ b*N(0)

n* - 0

на E

«а

(v*a(1)) - 0; [[v*a(1)]] - 0.

Здесь Ь*а(0), V*а(1) — компоненты вектора магнитной индукции и магнитный потенциал фаз композита ; а — 1,..., N в рамках одной ЯП, а Й* — V,*(0) — компоненты эффективного вектора напряженности магнитного поля; также обозначены , I — д/дх1 и / — д/д£ производные по двум типам координат; £— поверхность контакта двух фаз композита в пределах отдельно взятой ЯП, т.е. £а — Е^ав П Р^.

На искомые функции V*а(1) задачи (3) дополнительно наложены условия периодичности [ [V*а(1)] ] — 0 на границах ЯП, а также условие нормировки (V* а(1)) — 0, обусловленное требованиями единственно-

сти решения задачи на ячейке [4, 6], где

,*а(1)\ _

операция осред-

нения на ЯП. Условие нормировки делает задачу (3) в общем случае интегродифференциальной.

Преобразование локальной задачи к задачам "классического"

типа. Решения задачи (3) ищем в виде сумм:

v

*а(1) _

Е v га»

p=1

где v*0) _ функции следующего вида:

- -й; е+к).

(4)

(5)

Здесь §(р)(£к) — новые неизвестные функции от £к, уже не являющиеся периодическими. Производные от функций (4) по локальным координатам имеют вид

= Е v (РЖ = - Е H** + Е § (Рж = -H + Е § (Рж; (6)

;=1

;=1

;=1

;=1

тт*а(0) _ и* , „,*а(1)

Hi — Hi + v| i

3 3

H* - H* + £ ^*p)|i = £ 0

p=1 p=1

*a (p)|i .

Подставляя (6) и (7) в (3), находим, что функции являются решением следующей задачи:

а = 1,...,^;

i(p)/i °i(p)

= 0,

а £,

La _ ,. *a tj

°(p) - ^ Hi

i(p)'

Hi(p) = 0(p)|i'

0*а — ,a*N %) = %) ,

l*a _ h*N

°i (p) °i(p)

(

0(p)|q

— 0, — 0,

0p)

(

П — 0

|E, — 0,5Hp*;

— 0, q — p.

на £

(8)

Эта задача имеет "классический" тип, поскольку при фиксированном значении р она состоит из уравнений Лапласа (с комплексными неизвестными), записанных для областей (это пересечения областей Уа£ с 1/8 ЯП), с условиями идеального контакта на поверхностях контакта £и условиями первого или второго рода на координатных плоскостях Ер и торцевых плоскостях Ер = = 0,5}. Задачи (8) для 1/8 ЯП £ = У £^ будем далее обозначать как задачи Ьр.

а=1

Граничные условия для задачи Ь3 показаны на рис. 2.

Расчет тензора магнитной проницаемости композиционного материала. Для вычисления магнитного потенциала в нулевом приближении у*(0) из (1) и (2) получаем осредненную задачу теории

p

q

Рис. 2. Граничные условия для задачи L3

магнитостатики с эффективным тензором комплексных амплитуд магнитной проницаемости р*р, связывающим комплексную амплитуду вектора напряженности магнитного поля композита Й* и комплексную амплитуду вектора магнитной индукции композита (Ь*):

(Ь*) = Щ, (9)

где ( ) означает осреднение по У.

После решения серии задач Ьр для всех р проинтегрируем магнитную индукцию по областям, занятым волокнами и матрицей; тогда эффективный тензор комплексных амплитуд магнитной проницаемости р*р композита можно вычислить по следующим формулам:

Ь*

& = ; (10)

= О = Е / . (11)

а=1 v

Вариационная формулировка задач Ьр. Дадим вариационную формулировку задачи (8). Для этого рассмотрим класс комплексно-значных функций р*р) (возможные значения магнитного потенциала), которые определены во всей области у, являются гладкими в подобластях Уа£, удовлетворяют условиям непрерывности на границах фаз у, а также граничным условиям первого рода на торцевых поверхностях Ер, Ер:

К,)] = 0, на, ^ = 0, = 0,5Д*. (12)

Обозначим вариации возможных значений магнитного потенциала ^(*р), удовлетворяющие тем же условиям (12), но с нулеыми граничными условиями на торцевых поверхностях Ер. Тогда, следуя традиционной схеме вывода вариационного принципа [4], получаем, что истинный магнитный потенциал §(р), удовлетворяющий всем уравнениям задачи (8), отличается от всех возможных значений потенциала р*р) тем, что для него и только для него лагранжиан Ь имеет минимальное значение,

8Ь = 0, Ь = 11 р>(р)^(Р)^У. (13)

Соответствующее вариационное уравнение (вариационная формули-

ровка задачи ЬР) имеет вид

I ^ЧрЖН)!^У = 0. (14)

Метод конечного элемента для задач Ьр. Для решения вариационного уравнения (14) применим метод конечного элемента (МКЭ), согласно которому всю область У разбивают на конечные элементы (КЭ) Уе, = У Уе, и для каждого Уе записывают вариационное урав-

е

нение (14). При этом в правой части появляется ненулевое слагаемое, обусловленное потоком магнитной индукции Ь(р) = через

поверхность отдельного КЭ:

J = / Ь^и^Е. (15)

Зададим аппроксимации для магнитного потенциала в виде р(р) =

= {Ф}К)}Т, где Н)} = Н)(<Й) ...^(Р))} - координатная

строка, составленная из значений магнитного потенциала ) в

узлах КЭ; — координаты ]-го узла, ] = 1... т, т — число узлов; {Ф} = {Ф1(^г) ...Фт (<Т)} - строка функций формы. Вычисляя вариацию магнитного потенциала 5к(р) = {Ф}^{^(<р)}т, затем представляя производные от потенциала в виде координатного столбца Н)|1... ^(р)|з} = где {Ь} = {д/д£1 ...д/д£3} — строка

операторов дифференцирования, и вводя матрицу [В] = {Ь}Т{Ф} ({^(р)|1 ...^(р)|3}Т = [В]{^(р)}), переписываем вариационное уравнение (15) в виде

I ¿Н)}Т[В]Тр* [В]{^У = 18{1У{Р)}{Ф}ТЩР)<1Е. (16)

Вынося вариацию за знак интеграла, из (16) получаем разрешающую систему линейных алгебраических уравнений (СЛАУ)

[К *]К„)}Т = {/}Т, (17)

>)J

где [K*] = f [Б]T^[B]dV; {f }T = f {Ф}тb^p)d£ - локальная матри-

ца жесткости и столбец правых частей.

При численной реализации МКЭ были выбраны конечные элементы в форме 4-узловых тетраэдров, функции формы которых имеют линейный вид. Из локальной СЛАУ (17) составили глобальную СЛАУ для всей рассматриваемой области у, решение которой находили путем разделения комплексных переменных на действительную

и мнимую части с привлечением QMR-методов решения СЛАУ для возникающих в таком подходе несимметричных матриц.

Расчет для 3Б ортогонально-армированного композиционного материала. С помощью разработанного метода решены задачи Ьр для 3Б ортогонально-армированного композита (число N = 4, а = 1, 2, 3 — волокна, а = 4 — матрица). Число КЭ при решении задач Ьр составляло 6343, число степеней свободы — 1502. При численной реализации учитывалось, что армирующие в трех направлениях волокна могут состоять из различных материалов.

Проведено 10 тестовых расчетов для 3Б ортогонально-армированного композиционного материала с одинаковыми волокнами (см. рис. 1 и 2) с модельными значениями коэффициентов магнитной проницаемости волокон и матрицы -/,-т. Были выбраны следующие соотношения между -/,-т: Яе-т/ Ке-/ = 0,1, 1т-т/ Ке-т = 0,1, 1т-// Яе-/ = 0,1. По формулам (10), (11) были вычислены эффективные значения коэффициента магнитной проницаемости композита. Для сравнения проводились также вычисления вилки Фойгта-Рейсса (линейная и обратно-линейная зависимости)

* * , * (л \ 1 V/ , (1 - V/) /юч

- = VI + -т(1 - VI); — = + --*—, (18)

- -/ -т

где V/ — коэффициент армирования.

На рис.3,а изображена действительная часть коэффициента магнитной проницаемости композита Яе-* в зависимости от V/, на рис. 3, б — мнимая часть 1т-*, а на рис. 4 — абсолютная величина |-* |. Результаты расчетов показывают, что значения Яе-*, 1т-* и |-* | для 3Б ортогонального композита укладываются в вилку Фойгта-Рейсса,

Рис. 3. Зависимости действительной (а) и мнимой (б) частей магнитной проницаемости 3Б композита от коэффициента армирования, рассчитанные по методу Фойгта-Рейсса (Ф, Р) и методом асимптотического осреднения (звездочки)

что свидетельствует о хорошей точности предложенного метода расчета эффективных характеристик. В то же время сами соотношения (18) Фойг-та или Рейсса не могут быть выбраны даже в качестве приближенных выражений для эффективного коэффициента магнитной проницаемости ввиду того, что вилка Фойгта-Рейсса является слишком широкой. Различие между значениями (|, вычисленными по предложенному методу и по соотношениям Фойгта или Рейсса, составляет около 60% для коэффициента армирования ^^ = 0,55, что является весьма плохим результатом с точки зрения прогнозирования эффективных свойств композитов. Таким образом, для расчета эффективных магнитных свойств композитов со сложными структурами армирования необходимо использовать более точные методы, чем метод Фойгта-Рейсса, например метод асимптотического осреднения, рассмотренный в данной работе, который позволяет вычислять более точные в математическом смысле значения эффективных характеристик композитов.

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

СПИСОК ЛИТЕРАТУРЫ

1. Санчес-Пал енсия Э. Теория колебаний и неоднородные среды. - М.: Мир, 1984.-472 с.

2. Бахвалов Н. С., Панасенко Г. П. Осреднение процессов в периодических средах. - М.: Наука, 1984. - 352 с.

3. Победря Б. Е. Механика композиционных материалов. - М.: Изд-во МГУ, 1984.- 336 с.

4. Димитриенко Ю. И. Механика композиционных материалов при высоких температурах. - М.: Машиностроение, 1997. - 367 с.

5. Димитриенко Ю. И., Кашкаров А. И. Конечно-элементный метод для вычисления эффективных характеристик пространственно-армированных композитов // Вестник МГТУ им. Н.Э. Баумана. Сер. "Естественные науки". -2002.-№ 2.-С. 95-108.

\p*\/R е/4

0 0,2 0,4 0,6

Рис. 4. Зависимость абсолютного значения магнитной проницаемости от коэффициента армирования для 3Б композита, рассчитанная по методу Фойгта-Рейсса (Ф, Р) и методом асимптотического осреднения (звездочки)

6. ДимитриенкоЮ. И., КашкаровА. И., МакашовА. А. Разработка конечно-элементного метода решения локальных задач теории упругости "на ячейке периодичности" для композитов с периодической пространственной структурой // Математика в современном мире / Под ред. Ю.А. Дробышева. -Калуга.: Изд-во КГПУ, 2004. - С. 177-191.

7. ДимитриенкоЮ. И., Соколов А. П., КашкаровА. И. Разработка конечно-элементного метода решения задач расчета эффективных характеристик композиционных материалов на многопроцессорных вычислительных системах // Аэрокосмические технологии. - М.: Изд-во МГТУ им. Н.Э. Баумана, 2004. - С. 113-114.

8. Димитриенко Ю. И., С о к о л о в А. П. Разработка системы автоматизированного вычисления эффективных упругих характеристик композитов. // Вестник МГТУ им. Н.Э. Баумана. Сер. "Естественные науки". - 2008. - № 2. -C. 57-67.

Статья поступила в редакцию 13.10.2009

Юрий Иванович Димитриенко родился в 1962 г., окончил в 1984 г. МГУ им. М.В.Ломоносова. Д-р физ.-мат. наук, профессор, заведующий кафедрой "Вычислительная математика и математическая физика" МГТУ им. Н.Э. Баумана, действительный член Академии инженерных наук. Автор более 160 научных работ в области вычислительной механики, нелинейного тензорного анализа, термомеханики композитов, математического моделирования в материаловедении.

Yu.I. Dimitrienko (b.1962) graduated from the Lomonosov Moscow State University in 1984. D. Sc. (Phys.-Math.), professor, head of "Computational Mathematics and Mathematical Physics" department of the Bauman Moscow State Technical University. Full member of the Russian Academy of Engineering Sciences. Author of more than 140 publications in the field of computational mechanics, nonlinear tensor analysis, thermomechanics of composite materials, mathematical simulation in science of materials.

Евгений Сергеевич Ничеговский родился в 1985 г., окончил в 2007 г. МГТУ им. Н.Э. Баумана. Аспирант кафедры "Вычислительная математика и математическая физика" МГТУ им. Н.Э. Баумана. Автор ряда научных работ в области численных методов в механике композитов.

E.S.Nichegovskiy (b. 1985) graduated from the Bauman Moscow State Technical University in 2007. Post-graduatу of "Computational Mathematics and Mathematical Physics" department of the Bauman Moscow State Technical University. Author of some publications in the field of the numerical methods in composite mechanics.