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

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

1. Ландау, Л.Д. Теоретическая физика. В 10 т. Т. VII. Теория упругости [Текст]: Учеб. пос. / Л.Д. Ландау, Е.М. Лифшиц. — М.: Физматлит, 2007. — 264 с.

2. Жилин, П.А. Прикладная механика. Основы теории оболочек [Текст]: Учеб. пос. / П.А. Жилин. — СПб.: Изд-во Политехи. ун-та, 2006. — 167 с.

3. Карпов, В.В. Математические модели термоупругости оболочек переменной толщины при учете различных свойств материала [Текст] / В.В. Карпов, В.Н. Филатов // Вестник гражданских инженеров. — 2006. - № 3. - С. 42-45.

4. Карпов, В.В. Нелинейные математические модели деформирования оболочек переменной толщины и алгоритмы их исследования [Текст]: Учеб. пос. / В.В. Карпов, О.В. Игнатьев, А.Ю. Сальников. — М.: АСВ; СПб.: СПбГАСУ, 2002.- 420 с.

5. Жгутов, В.м. Математическая модель, алгоритм исследования и анализ устойчивости нелинейно-упругих ребристых оболочек при больших перемещениях [Текст] / В.М. Жгутов // Научно-технические ведомости СПбГПУ. Физико-математические науки. - 2009. - № 4. - С. 24-30.

6. Жгутов, В.м. Математические модели, алгоритм исследования и анализ устойчивости ребристых оболочек с учетом ползучести материала при конечных прогибах [Текст] / В.М. Жгутов // Научно-технические ведомости СПбГПУ. Физико-математические науки. - 2010. - № 2. - С. 53-59.

7. Жгутов, В.м. Математические модели, алгоритм исследования и анализ устойчивости упругих

ребристых оболочек при конечных прогибах [Текст] /

B.М. Жгутов // Научно-технические ведомости СПбГПУ. Физико-математические науки. - 2011. -№ 1. - С. 122 - 129.

8. Абдикаримов, Р.А. Математические модели задач нелинейной динамики вязкоупругих орто-тропных пластин и оболочек переменной толщины [Текст] / Р.А. Абдикаримов, В.М. Жгутов // Инженерно-строительный журнал. - 2010. - № 6. -

C. 38 - 47.

9. Абдикаримов, Р.А. Математические модели задач нелинейной динамики вязкоупругих изотропных пластин и оболочек гладко-переменной толщины (асимметричный случай) [Текст] / Р.А. Абдикаримов,

B.М. Жгутов // Инженерно-строительный журнал. -2010. - № 8. - С. 47 - 55.

10. Жгутов, В.м. Метод конструктивной анизотропии для ортотропных и изотропных ребристых оболочек [Текст] / В.М. Жгутов // Инженерно-строительный журнал. - 2009. - № 8. - С. 40 - 46.

11. Жгутов, В.м. Ответ профессору Карпову Владимиру Васильевичу (о научном приоритете в методе конструктивной анизотропии для ребристых оболочек и на функционал, описывающий ползучесть их материала) [Текст] / В.М. Жгутов // Инженерно-строительный журнал. - 2011. - № 3. -

C. 75 - 80.

12. Жилин, П.А. Теоретическая механика. Фундаментальные законы механики [Текст]: Учеб. пос. / П.А. Жилин. - СПб.: Изд-во СПбГТУ, 2003. - 339 с.

УДК 539.3, 537.226.4

Н.Г. Осипова, А.С. Семёнов

МОДЕЛИРОВАНИЕ НЕЛИНЕйНОГО ПОВЕДЕНИЯ ПЬЕЗОКЕРАМИКИ ТЕТРАГОНАЛЬНОй СТРУКТУРЫ методами конечно-элементной ГОМОГЕНИЗАЦИИ

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

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

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

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

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

В отличие от простейших подходов, не учитывающих взаимное влияние монокристаллов друг на друга и основанных на простом осреднении откликов в предположении равенства воздействий во всех монокристаллах [1—3] (приближение Рейсса) или упрощенно учитывающих это влияние при применении самосогласованных моделей [4, 5], используемый подход, основанный на непосредственном рассмотрении представительного объема поликристалла и применении МКЭ для его анализа, позволяет детально описать взаимное влияние отдельных кристаллитов друг на друга. Учет взаимного влияния и возникающих вследствие этого микронапряжений важен при анализе поведения при сложном непропорциональном и циклическом нагружении.

Определяющие уравнения микроэлектромеханической модели

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

был предложен в работе [4] и позволил учесть мультидоменную структуру и диссипативный характер движения доменных стенок.

Предполагая, что тензор деформации £ и вектор электрической индукции D допускают декомпозицию в виде суммы линейных (обратимых) £-, D'- и остаточных (аналог пластических) Pr- составляющих, и, используя для первых уравнения линейного пьезоэлектрического отклика, получаем определяющие уравнения в виде:

е = е' + ег = 48^ ..о + Е•3й + ег; Б = Б' + Рг = 3 й •ст + к0- Е + Рг,

(1)

где о — тензор напряжений, E — вектор напряженности электрического поля, £г — тензор остаточной деформации, Pr — вектор спонтанной поляризации, — тензор модулей упругой податливости кристалла (4-го ранга), М — тензор пьезоэлектрических модулей кристалла (3-го ранга) и К — тензор модулей диэлектрической проницаемости (2-го ранга).

Обращенная форма определяющих уравнений (1), используемая при проведении КЭ-расчетов, допускает следующее представление:

о = 4 С° •• (е-ег) - (Б - Рг) •3 Ь; Е = -3 Ь •• (е-ег) + рЕ-(Б - Рг),

(2)

где тензоры 4CD, %, Ре, характеризующие линейные электромеханические свойства монокристалла, определяются на основе инверсии М, к°. В дополнение к представлениям (1) и (2) необходимо ввести уравнения для £г и Pr.

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

= X С1 е I;

I=1

(3)

Рг = X С1Р,

I=1

где ег - концентрация (объемная доля) /-го домена в монокристалле, удовлетворяющая ограничениям

0 < с < 1, £ ^ = 1. (4)

I=1

Модули монокристалла 3d, к° определяются на основе модулей отдельных доменов &7 , к1 соотношениями, аналогичными (3):

4 = £ ^ 4 Б?,

I=1

3Й = £ С1 Ч> (5)

I=1

к° = £ С1К.

I=1

Изменение концентрации /-го домена монокристалла выражается через скорость переключения из состояния с ориентацией / в со-стяние с ориентацией /:

30

= Е л7 а /а,

(6)

а=1

где А1 а = 1, если /-й домен является реципиентом а-ой системы переключений (а: / ^ I); А1 а = -1, если /-й домен является донором а-й

системы переключений (а: I ^ /); А1 а = 0 в остальных случаях.

С учетом выражений (3) и (6) скорости остаточных деформаций и поляризации могут быть представлены в виде

6 30

ё" = Е ¿1ёI = Е /°КУа;

1=1 а=1

6 30

Р' = Е ¿IР1 = Е /а«аРа,

(7)

1=1 а=1

1

где ^а = Па + Па8а> - тензор ШмВДа; -единичный вектор в направлении изменения поляризации; уа, Ра - константы материала.

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

8 =

1 • 30

ст-ё + Е-1)--(о-в' +Е-Б') =£ Оа/а, (8)

2 а=1

на основе которого с учетом уравнений (1), (3),

(5)-(7) вычисляется движущая сила Оа , со-

/г/

:

О а =

+ Е-

1 6

уаца+ -£ А1 а(ст--4Sf +Е- 3й7) 21=1

1 6

Ра8а+ -X А1 а(3й/ "Ст + к?Е)

СТ +

(9)

Уравнение эволюции /а, удовлетворяющее условию неотрицательности диссипации 8 > 0 , принимаем в следующем виде:

ра

& а _ в а о_

_ оа

о0

п-1 /

,йопог(а) \т

С

(10)

где Оа > 0,Ва > 0, п > 0, т > 0,С0 > 0 - константы материала, определяющие форму гистере-зисных кривых; (а) - концентрация /-го домена (донора а-й системы переключений).

Для описания склерономного поведения следует выбрать п >> 1, как обычно и поступают при анализе пластичности кристаллов. Введение последнего сомножителя в уравнение (10) позволяет описать эффект насыщения и удовлетворить неравенствам (4).

Представительный объем поликристалла

В качестве представительного объема поликристалла рассматривался куб с различным числом разбиений (рис. 1). Каждая Гауссова точка конечного элемента рассматривалась как отдельный монокристалл; ориентация таких монокристаллов генерировалась случайным образом (нормальное распределение для каждого их трех углов Эйлера). При проведении вычислительных экспериментов с использованием КЭ комплекса РАЭТОСЯАТОЯ [6]

Рис. 1. Примеры представительных объемов поликристалла сегнетоэластика с различным числом разбиений: 2 х 2 х 2 (а), 5 х 5 х 5 (б), 10 х 10 х 10 (в). Показаны поля распределения угловых отклонений ориентаций (001) кристаллов от вертикали

Рис. 2. Распределение вертикальных компонент электрического поля Е' (а) и тензора напряжений а (б). Деформированное состояние для наглядности увеличено в 300 раз

исследовалось поведение различных вариантов поликристаллов с разбиениями от 1 х 1 х 1 до 10 х 10 х 10 элементов, содержащих от 8 до 8 000 монокристаллов и имеющих от 48 до 48 000 доменов соответственно.

На рис. 2 представлены характерные результаты КЭ-вычислительных экспериментов с представительным объемом поликристалла при задании граничных условий п -о|^ = 0, ф|^ = . Полученные в результате расчетов поля обладают ярко выраженной неоднородностью.

Целью проведения расчетов является определение эффективных значений полей ё, а, Б, Е:

е = — \гйУ; ст = — [стйУ;

У-1 у■>

у У у У

- 1 г - 1 ,

Б = — ( Б йУ; Е = — ( Е йУ УУ

УУ УУ

(11)

и установления связей между ними при различных способах задания граничных условий.

В данной работе основное внимание уделяется исследованию нелинейных характеристик, определяющих форму механического и диэлектрического видов гистерезиса. Определение линейных эффективных характеристик 3ё, к ? выполнялось стандартными методами гомогенизации для трансверсально-изотропного материала.

Формулировка краевой задачи для представительного объема поликристалла включает уравнения в объеме V

V-a = 0; е = (Vu)S ; V-D = 0; E = -Vy,

(12)

определяющие уравнения (2) и граничные условия на внешней поверхности S:

n-о S = n о; u S = e-г; n - D| S = n - D; ф| S =E- г.

(13)

При получении КЭ-решений связанных электромеханических краевых задач (12), (13) использовалась векторно-потенциальная формулировка [7, 8], основанная на вариационной формулировке с выбором в качестве варьируемых величин вектора перемещения и и электрического векторного потенциала ш [8]:

f(ст-- ôe + E-ÔD) dV =

V

= f (n - ст-ôu + n xE -ôy) dS.

(14)

Сравнение прогнозов моделей конечно-элементной гомогенизации с экспериментами

На рис. 3 представлено сравнение результатов КЭ-моделирования поведения поликристалла (представительный объем с разбиением 3 х 3 х 3) при циклическом нагружении с результатами экспериментов для трех различных пьезокерами-ческих материалов: PIC 151 [9], PZT-5H [10], PLZT 8/65/35 [11]. Во всех случаях наблюдаются лишь незначительные отклонения результатов расчета от эксперимента, что указывает на возможность использования рассматриваемой модели материала для адекватного описания по-

D,. Кл/м

-1,5 -1,0 -0,5

Dz, Кл/м 0,30

0.15

-0.15

-0,30

-1,0

1,0 1,5

Е... Ш/м

1,0

F,.,- MR/m

Рис. 3. Сравнение результатов КЭ-расчета (линии ) с экспериментом (символы ) для трех материалов: PIC 151 (a), PZT-5H (б), PLZT 8/65/35 (в)

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

Константы материалов, используемые в КЭ расчетах

Величина Обозначение Размерность Материал

Р1С 151 Р2Т-5Н РЬ2Т 8/65/35

Упругие модули С1111 С3333 С1133 МПа 1,26105 1,18105 0,53105 0,75-105 0,34105 0,21105 1,65-10п 1,32-Ю11 0,95-Ю11

Диэлектрическая проницаемость К33 Ф/м 2,2010-8 8,810-8 5,6310-8

Пьезоэлектрические модули ¿33 ¿31 ¿15 м/В 3,151010 —1,28-10-10 4,82-10-10 2,2-10-10 —1,11010 1,810-10 5,9410-10 —3,9510-10 3,70-10-10

Спонтанная поляризация и деформация р0 £0 Кл/м2 % 0,57 0,3 0,49 0,3 0,38 0,05

Критическое значение движущей силы В-Кл/м3 0,44106 0,38-106 0,18106

Коэффициент вязкости В90 — 0,00015 0,05 1

Показатели степени уравнения эволюции п т — 12 2,5 15 1,5 21 1,4

Параметр насыщения Со — 0,01 0,01 0,01

Сравнение прогнозов моделей КЭ-гомогенизации с прогнозами моделей прямого осреднения и монокристалла

Основным преимуществом рассматриваемой модели сегнетоэлектроупругого материала, основанной на КЭ-гомогенизации, является учет взаимного влияния отдельных кристал-

литов друг на друга. Это наглядно проявляется в сравнении с прогнозами других моделей. На рис. 4 показано сравнение кривых диэлектрического гистерезиса, полученных на основе модели гомогенизации (представительный объем с разбиением 2 х 2 х 2), модели прямого осреднения (60 монокристаллов) и модели одиночного монокристалла с результатами эксперимента

Кл/м2

/ ^ггг--

й- 3 А>/ 1

/ Г-2 -'-

-1-

Еъ, МВ/м

Рис. 4. Сравнение результатов прогнозов модели КЭ-гомогенизации (1) с моделью прямого осреднения (2) и с моделью монокристалла (3)

для материала PIC 151 [9]. Неучет взаимного влияния приводит к увеличению высоты гистерезиса и резкому переходу от его верхних и нижних сторон к боковым (соответствующих процессам переключений).

Недостатком рассматриваемой модели является значительное увеличение времени счета, связанное с необходимостью КЭ-решения задачи. Для получения решения задачи гомогенизации (64 монокристалла, 162 степени свободы в КЭ-модели, 1800 инкрементов нагрузки) потребовалось в 16 раз больше времени, чем для задачи прямого осреднения (60 монокристаллов, 500 инкре-

ментов), и в 69 раз больше, чем в расчете с использованием модели монокристалла (500 инкрементов).

Влияние числа монокристаллов в представительном объеме поликристалла

Исследование влияний числа разбиений представительного объема поликристалла на форму кривой диэлектрического гистерезиса показало, что начиная с числа разбиений 3 х 3 х 3 (216 монокристаллов), кривые становятся визуально неотличимы (рис. 5). В данной серии расчетов разбиение 1 х 1 х 1 дает границу сверху, а 2 х 2 х 2 является нижней границей.

Dz, Кл/м2

i -5

Е,, МВ/м

Рис. 5. Влияние числа разбиений представительного объема поликристалла на кривые гистерезиса: 1 — 1х1х1 разбиений (8 кристаллов), 2 — 2х2х2 разбиений (64 кристалла), 3—5 — 3х3х3—5х5х5 разбиений (216—1000 кристаллов)

Анализ влияния механических напряжений на диэлектрический гистерезис

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

меняющемся электрическом поле. На рис. 6 показаны расчетные кривые гистерезиса при различных значениях сжимающих напряжений. Наблюдается удовлетворительное совпадение кривых,полученных на основе применения мо -дели КЭ-гомогенизации, с экспериментальными результатами для материала PIC 151 [9]. Несколько заниженный уровень деформаций объясняется неучетом ромбоэдрической фазы в расчетах.

а)

Dz, Кл/м2

0,4

б)

02

-ОД.

-0.4,

и

[У^З 4

1 2 Е,, МВ/М

%

0,1

-0,1 -ОД -0,3 -0,4

-2

1

^•3

1 2 Е7, МВ/М

Рис. 6. Зависимость формы петли диэлектрического (а) и электромеханического (б) гистерезисов от уровня механических сжимающих напряжений ст, МПа: 0 (1), —50 (2), —200 (3), —400 (4)

Моделирование диэлектрического гистерезиса с выдержками

В условиях нагружения с различными скоростями и программах нагружения с промежуточными выдержками проявляются вязкие свойства сегнетоэлектроупругих материалов. В целях проверки возможности предсказания данных эффектов рассматриваемой моделью были проведены расчеты цикла нагружения с амплитудой 2 МВ/м со скоростью нагружения 0,08 МВ/(мс) при наличии промежуточных выдержек в течении 300 с через каждые 0,5 МВ/м

а)

Б, Кл/м2

0,40

(т. е. при Ez =±0,5; ± 1,0; ± 1,5; ±2.0 МВ/м). Результаты подобных испытаний представлены в работе для материала PIC 151 [9].

На рис. 7 показаны расчетные кривые диэлектрического и электромеханического гистерезисов в форме бабочки для рассматриваемой программы нагружения с выдержками. При фиксированных постоянных уровнях нагрузки с течением времени наблюдается ярко выраженная ползучесть. Полученные результаты демонстрируют качественное совпадение с результатами опытов [9].

0,090

0,060

0,030

1,0 2,0 Ez, МВ/м

1,0 2,0 Ez, МВ/м

Рис. 7. Расчетные кривые диэлектрического гистерезиса при наличии выдержек в процессе нагружения

е„. %

Итак, результаты проведенных многочисленных вычислительных экспериментов с использованием микромеханически мотивированной модели сегнетоэлектроупругого материала, основанной на определении эффективных свойств поликристалла методами конечно-элементной гомогенизации, проде-

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

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

1. Hwang, S.C. Ferroelectric/ferroelastic interactions and a polarization switching model [Text] / S.C. Hwang, C.S. Lynch, R.M. McMeeking // Acta Met. - 1995. -Vol. 43. - P. 2073-2084.

2. Huber, J.E. Multi-axial electrical switching of a ferroelectric: Theory versus experiment [Text] / J.E. Huber, N.A. Fleck // J. Mech. Phys. Solids. - 2001. - Vol. 49. -P. 785-811.

3. Belov, A.Y. Viscoplastic models for ferroelectric ceramics [Text] / A.Y. Belov, W. S. Kreher // J. Eur. Ceram. Soc. - 2005. - Vol. 25. - P. 2567-2571.

4. Huber, J.E. A constitutive model for ferroelectric polycrystals [Text] / J.E. Huber, N.A. Fleck, C.M. Landis, R.M. McMeeking // J. Mech. Phys. Solids. - 1999. -Vol. 47. - P. 1663-1697.

5. Landis, C.M. A self-consistent constitutive model for switching in polycrystalline barium titanate [Text] / C.M. Landis, R.M. McMeeking // Ferroelectrics. -2001. - Vol. 255. -P. 13-34.

6. Семёнов, А.С. PANTOCRATOR - конечно-элементный программный комплекс, ориентированный на решение нелинейных задач механики [Текст] / А.С. Семёнов // Тр. V-й Междунар. конф. «Научно-

технические проблемы прогнозирования надежности и долговечности конструкций и методы их решения». СПб.: Изд-во СПбГПУ, 2003. - С. 466-480.

7. Семёнов, А.С. Эффективные методы решения нелинейных краевых задач сегнетоэлектроупругости [Текст] / А.С. Семёнов, А.Ч. Лисковски, П. Ноймай-стер [и др.] // Морские интелектуальные технологии. - 2010. - № 1. - С. 55-61.

8. Semenov, A.S. Return mapping algorithms and consistent tangent operators in ferroelectroelasticity [Text] / A.S. Semenov, A.C. Liskowsky, H. Balke // International Journal for Numerical Methods in Engineering. - 2010. -Vol. 81. - Р. 1298-1340.

9. Zhou, D. Experimental investigation of non-linear constitutive behavior of PZT piezoceramics [Text] / D. Zhou // Dissertation. Karlsruhe, 2003. - 139 p.

10. Huber, J. Ferroelectrics: models and applications [Text] / J. Huber // Dissertation. Cambridge, 1998. -177 p.

11. Lynch, C.S. The effect of uniaxial stress on the electro-mechanical response of 8/65/35 PLZT [Text] / C.S. Lynch // Acta Materialia. - 1996. - Vol. 44. -Р. 4137-4148.

УДК 519.62: 537.8

М.М. Корсун, М.Э. Рояк

УЧЕТ ШИХТОВАННОСТИ МАТЕРИАЛОВ ПРИ МОДЕЛИРОВАНИИ ЭЛЕКТРОМАГНИТНЫХ ПРОЦЕССОВ

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

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