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

МЕХАНИКА

УДК 539.3

ЧИСЛЕННО-АНАЛИТИЧЕСКОЕ ПОСТРОЕНИЕ МАТРИЦ ГРИНА ТРЕХМЕРНЫХ ТЕОРИЙ УПРУГОСТИ И ЭЛЕКТРОУПРУГОСТИ

© 2010 г. Л.А. Игумнов, С.Ю. Литвинчук, В.П. Пазин, А.Н. Петров

НИИ механики Нижегородского госуниверситета им. Н.И. Лобачевского

[email protected]

Поступила в редакцию 29.03.2010

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

Ключевые слова: матрица Грина, трехмерная постановка, анизотропная упругость, электроупругость, численное моделирование.

Введение

Матрица Грина важна как фундаментальное решение соответствующей системы дифференциальных уравнений, а также для построения граничных интегральных уравнений соответствующих краевых задач. Получение упругостатической матрицы Грина в трехмерной постановке имеет свою историю. Из соответствующего корпуса публикаций выделим работы следующих авторов: Freedholm (1900), Lifshitz и Rozenzweig (1947), Kroner (1953), Synge (1957), Willis (1965), Mura и Kinoshita (1971), Pan и Chou (1976), Wang (1997). В итоге, в трехмерной постановке строятся выражения для упругостатических матриц Грина анизотропных материалов при помощи преобразования Радона и контурного интегрирования. Ссылки на первые работы по детальному обсуждению упругостатических матриц Грина для анизотропного случая можно найти в обзоре Bacon и др. (1978), а также в работах Mura (1987) и Ting (1996). Исследования упругостатических матриц Грина электроупругости для бесконечных и полубес-конечных трансверсальных изотропных тел в замкнутой форме были проведены Wang и Zheng (1995), Dunn и Wienecke (1996, 1998), Ding и др. (1996-2000) и Pan (2002). Определенные итоги подведены в [1].

1. Постановка задач

Приведем базовые уравнения трехмерной линейной теории электроупругости, которые можно найти, к примеру, в [2, 3]:

А,] + = 0, А,; - б = 0, г'>.] = 1,3,

= Су1тг1т - еку Ек, А = еукУ]к + Е],

I, у, к, I, т = 1,3,

еу = 1(и,-,у + иу,), Е, = -ф„-, г, у = 1,3, где Сту, Д - соответственно компоненты тензора напряжения и вектора электрического смещения; ^, б - компоненты плотностей объемной силы и плотности электрического заряда; - компоненты тензора деформации,

Еу - компоненты вектора электрического поля;

су1ш , еук и - упругие модули электр°уп-

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

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

единых позиций, группируем базовые компоненты следующим образом [4]:

Гм., г = 1, 2, 3, к;, г = 1,2,3,

М=\ €=\ 9

(ф, г = 4,

I- Е;, г = 4,

[ст„, 7 = 1,2,3,

€ =\ 4 4 [А, 7 = 4,

С/И, 7, к = 1,2,3,

е/г/,

7 = 1,2,3, к = 4,

, 7 = 4 к =1,2,3

-еа, ], к = 4.

3,г 3

где

°Ук = А I(г,- ))-18(г,-г Жг,) =

8п *2

= -^ I(М% (г, ))-18(г,г0 )Ж(г,),

8п2г ^2

г° = “ , 8(г,-г) = 8(гг,г°) =18(г,-г°), г г

*«■ .

2

Интеграл по сфере £ : = 1 преобразу

ется к интегралу по окружности:

Единичная окружность является пересечением единичной сферы с плоскостью 2г- Г° = 0.

Тогда физические уравнения и уравнения равновесия запишутся следующим образом:

Г^., 7 = 1, 2,3,

= 1- е, 7=4.

2. Получение формул представлений матриц Грина

Математически матрица Грина может быть определена следующим образом:

ёцнбьт» (х) = -5 >8( х),

— (1) 1,7, к, 1, т = 1,3.

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

Матрица О^ с использованием обратного

преобразования Радона может быть записана следующим образом [5]:

Рис. 1. Пересечение единичной сферы с плоскостью

^ г0 = 0

Продемонстрируем второй метод. Обратную матрицу [м22 (г)) представим через матрицу алгебраических дополнений Л^ (г) и определитель г), тогда матрица Грина примет вид

[3, 6]:

в = _2_ Г (г,~) 5(2.г. (7,) =

7 8п2 5{ Б(7,) 11> 1>

1 +7 А#(Р + ^ ^ ^

К- - - ^•

(3)

4п г Ар + Сq)

Применение теории вычетов позволяет записать интеграл (3) в следующем виде:

Сд (X) = -: 1

4

аш X-

2 пг А/к (Р + С т^

(4)

(С т -С! ) П (^ т -С к )(С т Чк )

где ^ т - корни многочлена восьмой степени от

^ -0(^ + = 0 ; а9 - коэффициент при ;

1т^т > 0, да = 1,4 ; - сопряженное к ^т .

Итак, в отличие от интеграла (2), второй способ дает явную формулу (4), что гарантирует эффективность использования С(х) и точность получаемых результатов. Отметим, что при построении (4) все полюсы предполагались простыми. В случае кратных полюсов небольшие изменения в константах позволяют свести задачу к простым полюсам с незначительными ошибками в вычислениях расширенной матрицы Грина [7, 8]. При вычислении важно, что

матрица Мсимметрична, как и ее алгебраическое дополнение А, а значит, расширенная

матрица Грина симметрична [9] и необходимо вычислить лишь 10 элементов из 16.

т=1

а

к =1 ,к Фт

0

3. Численные эксперименты

Для удобства тензор упругих модулей Сук/ запишем через следующие константы

сар (а = 1,6; Р = 1,6) [10]. Электроупругие константы вщ могут быть записаны с двумя индексами екр (к = 1, 3;р = 1, 6) следующим образом:

ек1 = ек11; ек2 = ек22;

ек 3 = ек 33; ек 4 = ек 23 = ек 32;

ек 5 = ек13 = ек 31; ек 6 = ек12 = ек 21В качестве первого примера рассмотрим изотропный материал со следующими параметрами: ц = 493.6 Н/м2; v = 0.25 [6].

Для изотропного упругого материала матрица Грина может быть найдена аналитически и численно. Вычисления компонент упругой матрицы Грина были проведены в точке х = (1,1,1), и результаты представлены в табл. 1 в сравнении с аналитическим решением.

В качестве другого примера рассмотрим трансверсальный изотропный электроупругий материал, для которого известны точные решения в замкнутой форме [11]. Предположим, что ось симметрии материала параллельна оси Х3, тогда отличные от нуля элементы констант материала следующие:

c11, c22 (= c11), c33, c13, c23 (= c13 ), c44,C55(= C44),C66,CX2(= cn - 2c66),

Таблица 1

Значения индексов компонент ( j, к ) Аналитическое решение (м) Численное решение (м)

1,1 7.2395217310 х 10-5 7.2402710926 х 10-5

1,2; 2,1 1.0342173901 х 10-5 1.0344001002 х 10 -5

1,3; 3,1 1.0342173901 х 10-5 1.0341648562 х 10-5

2,2 7.2395217310 х 10-5 7.2392019103 х 10-5

2,3; 3,2 1.0342173901 х 10-5 1.0340308316 х 10-5

3,3 7.2395217310 х 10-5 7.2390925138 х 10-5

Таблица 2

cii С13 C12 c33 C44 c66 e3i e33 ei5 eii Б33

139.0 74.3 77.8 115.0 25.6 30.6 -5.2 15.1 12.7 6.4605 5.61975

Таблица 3

Значения индексов компонент ( j, к ) Аналитическое решение (м) Численное решение (м)

1,1 1.1507372825 х 10-12 1.1507373432 х 10-12

1,2; 2,1 1.9428530900 х 10-13 1.9428530549 х 10-13

1,3; 3,1 1.7241444282 х 10-13 1.7241445986 х 10-13

1,4; 4,1 2.0825816731 х 10-4 2.0825819231 х 10-4

2,2 1.1507372825 х 10-12 1.1507373215 х 10-12

2,3; 3,2 1.7241444282 х 10-13 1.7241445725 х 10-13

2,4; 4,2 2.0825816731 х 10-4 2.0825818427 х 10-4

3,3 7.9729216465 х 10-13 7.9729219583 х 10-13

3,4; 4,3 1.5239325664 х 10-3 1.5239326228 х 10-3

4,4 - 4.3915928465 х 106 - 4.3915929246 х 106

Таблица 4

cii С22 С33 C12 С13 С23 C44 C55 C45 c66 ci6 C26 C36

95.5 25.9 16.3 28.9 4.03 4.65 4.40 6.45 -1.78 32.7 44.7 15.6 0.54

е31, е33, е15, е32(- е3і), е24(- е15 ) є11,є22(- є11),є33 •

Для трансверсально изотропной пьезокерамики (Р2Т-4) [12] соответствующие константы материала приведены в табл. 2, упругие константы сар даны с точностью до множителя 109 Н/м2,

пьезоэлектрические коэффициенты екр - Кл/м2 и диэлектрические постоянн^іе ~у - 10-9 Кл/В^м.

о о

Рис. 2. Компонента Сц

Рис. 4. Компонента С13

Для точки X = (1,1,1) компоненты электро-упругой матрицы Грина представлены в табл. 3 в сравнении с точным решением [11].

4. Визуализация матриц Грина

Опишем построение матриц Грина в виде двумерных поверхностей на основе интерполяционной схемы вычисления. Функции Оук (х)

о о

Рис. 3. Компонента

о о

Рис. 5. Компонента С 22

о о

Рис. 6. Компонента С 23

Рис. 7. Компонента С33

Рис. 8. Компонента Сц

-13

х 10

5

о о

Рис. 10. Компонента С13

-12

и 10

3

0 о

Рис. 12. Компонента С22

зависят от трех пространственных переменных *1, Х2 и Х3. Перейдем к сферическим координатам (г, 01,02).

Для единичной сферы Сук (01,02) - функции только двух переменных - полярного угла 0 <01 < п и азимутного угла 0 <02 < 2п. Зна-

Рис. 9. Компонента С12

-4

о о

Рис. 11. Компонента С14

з< 10

5

О о

Рис. 13. Компонента С 23

чение Сд от конкретных 01 и 02 строится с

помощью интерполяции Лагранжа [5].

На рис. 2-7 представлены компоненты анизотропной матрицы Грина для графито-эпок-сидного моноклинического материала, параметры материала приведены в табл. 4. На рис. 8-17 представлены компоненты электроупругой

о о

Рис. 14. Компонента G24

-3

х 10

о о

Рис. 1б. Компонента G34

матрицы Грина для трансверсальной изотропной пьезокерамики (PZT-4).

Заключение

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

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

1. Qing-Hua Qin Green’s function and boundary elements of multifield materials. Elsevier, 2007. 254 p.

2. Партон В. З., Кудрявцев Б. А. Электромагнито-упругость пьезоэлектрических и электропроводных тел. М.: Наука, 1988. 470 с.

17-,.-""

0 0

Рис. 15. Компонента G33

к 10

о о

Рис. 17. Компонента G44

3. Pan E., Tonon F. Three-dimensional Green’s functions in anisotropic piezoelectric solids // International Journal of Solids Structures. 2000. Num. 37. P. 943-958.

4. Li X., Wang M. Three-dimensional Green’s functions for infinite anisotropic piezoelectric media // International Journal of Solids Structures. 2007. Num. 44. P. 1680-1684.

5. Баженов В.Г., Игумнов Л.А. Методы граничных интегральных уравнений и граничных элементов в решении задач трехмерной динамической теории упругости с сопряженными полями. М.: Физмат-лит, 2008. 352 с.

6. Gaul L., Kogl M., Wagner M. Boundary Element Methods for Engineers and Scientists. Berlin Springer, 2003. 488 p.

7. Pan E., Amadei B. Fracture mechanics analysis of cracked 2-D anisotropic media with a new formulation of the boundary element method // Int. J. Fracture.

1996. Num. 77. P. 161-174.

8. Pan E.A general boundary element analysis of 2-D linear elastic fracture mechanics // Int. J. Fracture.

1997. Num. 88. P. 41-59.

9. Chen T. Green's functions and the non-uniform transformation problem in a piezoelectric medium // Mech. Res. Comm. 1993. Num. 20. P. 271-278.

10. Ting T.C.T. Anisotropic Elasticity - Theory and Applications. New York: Oxford University Press, 1996.

11. Dunn M.L., Wienecke H.A. Green's functions for transversely isotropic piezoelectric solids // International Journal of Solids Structures. 1996. Num. 33. P. 4571-4581.

12. Dunn M.L., Taya M. An analysis of piezoelectric composite materials containing ellipsoidal in homogeneities // Proc. R. Soc. Lond. 1993. Num. A443. P. 265-287.

THE NUMERICAL-ANALYTICAL CONSTRUCTION OF GREENS MATRICES OF 3-D ELASTICITY AND ELECTRO-ELASTICITY THEORIES

L.A. Igumnov, S.Yu. Litvinchuk, V.P. Pazin, A.N. Petrov

Two methods to construct Green’s matrices on the basis of a unified formalization for Green’s matrices of 3-D theories of anisotropic elasticity and electro-elasticity have been presented. The specific feature of the first method is getting the Green’s matrices in the integral form while the second one employs an analytical formula which uses Green’s matrix root sets. Components of the obtained representations of Green’s matrices have been numerically constructed for particular material properties. The results obtained are compared with those of other authors. Green’s matrix visualization has been presented in the form of surfaces.

Keywords: Green’s matrix, 3-D formulation, anisotropic elasticity, electro-elasticity, numerical modeling.