Нестационарная задача прямого пьезоэффекта для асимметричной биморфной круглой пластины ступенчато переменной толщины

Механика деформируемого твердого тела

DOI: https://dx.doi.org/10.24866/2227-6858/2019-2-2 УДК 539.3

О. В. Ратманова

РАТМАНОВА ОЛЕСЯ ВИКТОРОВНА - ассистент кафедры строительной механики и сопротивления материалов, е-mail: olesya654@yandex.ru Самарский государственный технический университет Молодогвардейская ул., 244, Самара, 443001

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

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

Введение

Для оценки эффективности работы и исследования напряженно-деформированного состояния составных электроупругих систем применяется в основном классическая прикладная теория для тонких пластин, в которой кинематические гипотезы дополняются аналогичными допущениями о характере изменения напряженности электрического поля по аксиальной координате [2, 7, 10, 12]. В данной постановке скорость распространения сдвиговых деформаций является неограниченной. Для устранения этого недостатка в задачах электроупругости [9] используется уточненная двухмодовая теория Тимошенко [5].

При исследовании конструкций, несимметричных относительно срединной поверхности, необходимо рассматривать связанные мембранные и изгибные колебания. В такой постановке можно отметить ограниченное количество построенных решений. В частности, в работе [4] замкнутое решение получено для прямоугольных анизотропных пластин ступенчато переменной толщины. Анализ работы асимметричного биморфного преобразователя постоянной толщины в случае действия гармонической механической нагрузки проведен в [6]. Анализу продольно-поперечных гармонических колебаний пластины переменного сечения посвящено

© Ратманова О.В., 2019

О статье: поступила: 18.02.2019; финансирование: бюджет Самарского государственного технического университета.

исследование [11]. Наконец, математическая модель для расчета изотропных асимметричных пластин гладко переменной толщины разработана в [1].

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

Постановка задачи

Пусть круглая биморфная пластина состоит из металлической заземленной подложки радиусом Г = Ь и толщиной К2, а также пьезокерамического аксиально поляризованного элемента с гексагональной кристаллической решеткой класса 6шш (Г = а, а < Ь, К). Осесиммет-ричные изгибные колебания образуются в результате действия на верхней лицевой поверхности конструкции механической нагрузки , 0, являющейся произвольной функцией времени

и радиальной координаты Г (рис. 1). Подключение электродированных поверхностей пье-зокерамических пластин к измерительному прибору позволяет определить разность потенциалов V*(?,). Условия закрепления цилиндрической поверхности ( г = ь ) биморфной конструкции могут быть произвольными. Для определенности считаем ее жестко защемленной.

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

нижней и верхней поверхностей подложки, соответственно обозначены: г = К + е , г2 = е,

= е-Н2, е _ р(2Ч2-Р(1)К 2 (р(2 АЛ

г к

( р"1 ^, р"2 ^ — объемная плотность пьезокерамики и металла).

Ч*(г„, О

Рис. 1. Расчетная схема биморфной конструкции.

Дифференциальные уравнения радиальных и изгибных осесимметричных колебаний биморфной пластины относительно нормальных N(г,,^),Ыв(г, поперечных 0г(г*, О сил, а также изгибающих моментов мг(г„,г*,) , при использовании следующих кинематических гипотез

и *

имеют вид

ды. N - N д2и;

+ —-- - а —= 0

дг

1 а:

ВЕСТНИК ИНЖЕНЕРНОЙ ШКОЛЫ ДВФУ. 2019. № 2

дОг м а 52Ц* ч , V1 + — - = Я V* Ч,

дг г

дМг дг

1

м., - мс

+ а

д У

2 дг1

= 0

где и*,"*,цт- радиальные, аксиальные перемещения и угол поворота сечения в плоскости (г, ^); _ радиальные перемещения поверхности при z = 0;

«1 = р{2\Н (а - гт) + (р(1)к + р(2к) Н (г — а),

«2 =Р(2) 2 Н ( а-г,) +

/ зз з з Л

т ^ ~ ^2 + ^(2) ^2 ~ ^з

^ 3 3

V ^ ^ у

Н ( г,- а ):

Н (...) — единичная ступенчатая функция О. Хэвисайда.

Подключение электродного покрытия пьезокерамической пластины к измерительному прибору с большим входным сопротивлением позволяет, из условия отсутствия тока смещения на внешнем лицевом электроде =^ +е = 0 [3], определить аксиальную компоненту вектора напряженности электрического поля Ez (г, z, , т.е.

£ (г.,z,t.) = -—

'33

ди* и * --1--

V дг. г у

Интегральные характеристики напряженного состояния составной конструкции определяются в виде следующих соотношений:

% % ^

Мг(в) = I (Ш)^ , ^г = ^1(2) | ,

М

I

&гг ( Ш)^

где стгг,0'вв,0'Г1 — компоненты тензора механических напряжений, к2— коэффициенты поперечного сдвига.

В результате подстановки (4) в (2), с учетом гипотез (1), (3), получаем следующую систему дифференциальных уравнений в безразмерной форме:

Р1 -Р2 аз ^ = 0,

от ог от

Рз

д2Ц

д

Р + р з

о г

г дЦ л --^

Кдг

д_ д г

д >

- Р и 0-«44 = 0 , о t

где {Ц,ио,г,Я} = {ц*,ио*,г.,а}/Ь, t = Ь

4 С,

11

2) , Я = — Ь М ру ' т1

а1

аз = (2у

а0

а =

? 2> к2' 4 Р{%2 к2

Р1 т1 т5

Р2 -1 тзЬ1 т7Ь1

= т1

Рз т9 т10

Р4 т11Ь"2 т1зЬ 2

Н

Краевые условия в безразмерной форме имеют вид:

г = 1,0 и 0 (1, г) = о , W(1, ?) = 0, ^(1, г) = 0 ,

и0 (0, , w(0, ?)<<», ^(0, ,

? = 0 и0(г,0) = ц,0г), V (г,0) = Wо (г), ^(г,0) = ^0 (г),

ди0 . , V • / ч дц/ . , ч

-г

где ^2) - С11

+

2 Л

с (1) +

1(12 Г о

йзэ У

К

т7(8) _ тз(4) +

2 Л

с (1) + ез1 11(12) р

йзз У

^15 , ^ = , ^0 = К (С5(52)к2 + СК ) ,

- С(2) _ ^11(12) 11(12)^18 , Ш1з(14) _ Ш11(12) +

Г

с ез1

'11

т15

2 2 ¿1 - г2

2

т1б

з з

г -

з

т,

v

2 _ 2 %з

т,

16

'зз у

17

2

, т18

з _ 3

3

ез1, £зз — постоянные пьезокерамического материала; Ст1 — модули упругости электроупругого = 1) и упругого (.V = 2) материалов (т, к = 1,2,5); \У() , и{),!//,, , , ^о , '/А, - заданные

при ? = 0 перемещения и их скорости.

Для данного конструктивного решения электрическое напряжение холостого хода V*(?,) определяется следующим образом:

V(?,) = 5-1 \ф(п,¿1,?ф>5,

где ф(г„, г1,?») - потенциал электрического поля, генерируемый на нижней лицевой поверхности пьезокерамической пластины ( е =- ^^ ), 5 - площадь электроупругой пластины.

Окончательное выражение функции V (?„) имеет следующий вид:

V(0 = ^К[2и•(а,?.)-(К + 2е)¥(а,Г,)] .

^зз а

Построение общего решения

Для решения начально-краевой задачи (5)-(7) введем на сегменте [0,1] конечное интегральное преобразование (КИП) [8] с неизвестными компонентами Кх (Д, г),

К2 (Д., г), К3 (Д, г) вектор-функции ядра преобразования и весовыми функциями (Х,^,у:

О

да

{ц( г,?) V ( г,?) г,? )} = £ О (4? ){К1(^г ), ^.г ), Кз(^г )}|| К. ||

¿=1

= /0 (л, г+ (л, г)+^ (л, г,,

ВЕСТНИК ИНЖЕНЕРНОЙ ШКОЛЫ ДВФУ. 2019. № 2

где - безразмерные круговые частоты осесимметричных колебаний (/ = 1, связанные с размерными частотами (01 следующей зависимостью:

л

Щ = 4

И

C

р

.

Представляя кусочно-гладкие функции U0 ( r, t), W ( r, t) , r, t) в виде

{U0,W,¥} = U\W(a),¥(a))H(r-R)^,Wb,¥(b))H(R-r),

подвергаем систему уравнений (5) и начальные условия (7) преобразованиям КИП в соответствии со структурным алгоритмом [8]. В результате получаем начальную задачу для трансформанты G{Л¡, t), общий интеграл которой имеет вид

CÍA;, t)= Go cos(A¡t) + G0 sin(^it)/ Л - ^ £ F(Л,t)sin Л (t -r)dT. и, с учетом (6), однородную задачу для компонент ядра КИП

p ^ VK^-p2 A VK<V + = 0,

dr dr

Рз

r K j) Л

V -V jj

d

dr

p 4 VKT3 3 p3

dr

+а к2jj=0.

У

^ dK (j) Л dK2 _ K ( j j

dr3

d

- p 2 VKij) + «4 J' = 0.

dr

i2s-(j) -

r =

V У

1,0 к (Л ,1) = K (Л-Д) = K (Л,1) = 0, |Kb >(4,0), Kf >(Л,0), Kf >(Л,0)}<», r = *: K1(VR=Kf>(i,,R), ,) = ,

*3а)(; ) = K3b)(4,R ),

-

m0

dK'

rn1

U

dr

SU

dr

+ m-

+ m,

U0 m3 т4 щ

r b dr br

u0 m1 ^ m2 ^

r b dr b r

dr

\

r dU,

-K'

= m0

dK2b dr

--K

b

т

5

=R V

Sr

+ m.

U0 my тз щ r b dr b r

J\r=R

Ar=R

8U0 U0 m3 m4 m7—0 + ms — —13-

v

dr

b dr b r

Ar=R

Здесь K K2, K3} = {Kia>, K2a), K3a)} H (r-R) + {K(b>, K^, K3b>} H (R-r), F (4 ,t) = J0 а3 (,, rdr , ^ = ^ = 1 , / = ^ ,

а3

G0 0Л ) = Jo (^0K1 + ^W0K2 + mK3 3dr,

Go (Л) = Jo' («'цЛ, + pw0K2 + г'/^гЛз .

Равенства (16) определяют условия неразрывности деформаций и усилий в области нерегулярности структуру конструкции при г = Я, которые формируются при выполнении равенств:

и0а)(Я,) = и0Ь)(Я,), V(а)(Я,) = ^(Ь)(Я,?), а)(Я,) = <//Ь)(Я,), ы{а)(я, ) = ы{Ь)(я,), 0(а)(я, ) = е(Ь)(я,), м(а)( я,) = м (Ь)(я,).

Для решения краевой задачи (14)-(16) вводятся новые функции Я 1 , г) ,

Я21 (л„, г), яз11 (л„, г) по формулам:

¿К^1 ^(л г^

Ж,г) = ,г), Я^Я,,г) = , яЦ,г) = -Я11

и тогда система (14) приводится к следующему дифференциальному уравнению относительно

Яз-'ЧА,, г) :

где у2=т^-т, / , /=(ир )-1 аз^2« ,

¿г г ¿г и рз

/2 ^"Ч+^РзГ (аД2и2 +n5P() , И2=Р1(а4^^Рз)-Р22+азЛ2Р4,

и1 = Р1Р4, п3 = а3 Д2 (а4 Д2 - Р3), и4 = Р1 Рз, и5 = аз Л2 Рз. Использование порождающего уравнения

V2 Я1) =-Д2 Я1],

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

Отношение физических характеристик материалов, используемых при изготовлении данной конструкции, позволяет получить следующие корни уравнения (19): Дг, Дг - действительные положительные числа, _ мнимое число.

В этом случае общий интеграл дифференциального уравнения (18) записывается в виде 2

Я1ЧЛ, г н

х 4%/ 1

+ Щ]/1

где /„(•••),^и( - ),,К« (...)— обыкновенные и модифицированные функции Бесселя I и II рода порядка А = 0,1).

Выражения для компонент ядра преобразования К , г), к2, 1 (Д, г), К 1 (Д., г) определяются, с учетом (17) и связей, полученных в процессе приведения (14) к (18):

я2 1 = -[«^ + (« +И4)У2 + «]Я1) ,

азЛ«4 - Рз«5

Подстановка К 1 ,К21 ,К1 в граничные условия (15), (16) позволяет определить посто-

Я1 } = Р21 [_РзЯ] 1 + Р4^лС } + (а4Л|2 - Рз) Я!1 };

[одстановка К^1 янные ^• ••Д. ,^•••Дз.^ и собственные значения Л1.

ВЕСТНИК ИНЖЕНЕРНОЙ ШКОЛЫ ДВФУ. 2019. № 2

Применяя к выражению трансформанты (13) формулы обращения (10), получаем общие выражения для безразмерных динамических компонент перемещения ио (г, *), Ш (г,*) и проекции угла поворота сечения г, *) .

Пример расчета

В качестве примера рассматривается биморфная конструкция, состоящая из металлической подложки из латуни ( ^ = 0.5 х103 м, Ь = 3 х 10~2 м) и пьезокерамической пластины состава ЦТС-19. Используются следующие физические характеристики материалов: латунь -

2 3

Е = 9.8х 1010 Н/м , V = 0.35 , р = 8600 кг/м , е, V,р- модуль упругости, коэффициент Пуассона

и плотность изотропного материала ( с(2) = Е С1 ~ у) , с(2) =

11 "(1 + V )(1 - 2 V) 12

Ev

., с(2)=_Е_); пье-

зокерамика

- £33 = 7.26х109 Ф/м , е31 =-4.9 Кл/м

31

С ,С5(5)} = {10.9, 9.1, 6.1, 5.4, 2.4}х10

10 , 2 Н/м .

-11 >^55 у

Рассмотрим случай действия гармонической нагрузки:

4 ( ^ * ) = 40 С1" г ) 81П №, где 40, %— амплитуда нагрузки и частота вынужденных колебаний в безразмерной форме.

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

в анными поверхностями пьезокерамической пластины V ((пунктирная линия), а на рисунках 3, 4 показаны зависимости изменения амплитудных значений разности потенциалов V (*) при различных значениях геометрических размеров пьезокерамической пластины.

{ш ( 0, * ), V (* )х10 -4} / 4 0

1x10

500

500

1x10

/ \ ' * / / \ ' * / У %/ \ г V 7 V*^

4 * \ / * Л А> * \ / • ♦ \ / Л > -* ^^

100

200

300

400

Рис. 2. Изменение Ш(0,*), V (*) по времени * (^ = 0 5 х103 м, а = 1.8 х 102 м, %- 0.Ц).

Численный анализ результатов позволяет сделать следующие выводы: 1) при высокочастотном внешнем воздействии / = 0.8^ эффект «биения», наблюдаемый в конструкции при вычислении вертикальных перемещений Ш (0, , также распространяется на характер изменения V (*) по времени;

2

0

0

2) анализ результатов расчета, полученных для пластины постоянной и ступенчато переменной толщины, показывает, что уменьшение жесткости конструкции при a < b приводит к росту перемещений W (0, t) и разности потенциалов V (t);

3) для заданной металлической подложки (b = 3 х 102 м, h2 = 0.5 х103 м ) можно подобрать геометрические размеры пьезокерамической пластины, позволяющие наиболее эффективно преобраз овать механическую нагрузку в электрический импульс, в частности, в данном примере необходимо использовать следующие параметры: R = 0.62 ( a = 1.86 х 10"2 м), h1= 0.9х10~3 м.

V(t)/^

80

75

70 65

0.5 0.6 0.7 0.8 0.9

Рис. 3. Зависимость амплитудного значения V (t) от радиуса пьезокерамической пластины R (% = 0.2Д).

V(t)/q,

76 75 74 73 72

71

- 4 - 4 - 3 - 3

8x1 0 4 9x1 0 4 1 x1 0 3 1.1 x1 0 3

Рис. 4. Зависимость амплитудного значения V (t) от толщины пьезокерамической пластины h1 (^ = 0.2^).

Заключение

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

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

1. Абдикаримов Р.А., Жгутов В.М. Математические модели задач нелинейной динамики вязкоуп-ругих изотропных пластин и оболочек гладко-переменной толщины (асимметричные случаи) // Инженерно-строительный журнал. 2010. № 8. С. 47-55.

R

h, м

2. Ватульян А.О., Рынкова А.А. Об одной модели изгибных колебаний пьезоэлектрических би-морфов с разрезными электродами и ее приложениях // Изв. РАН. МТТ. 2007. № 4. С. 114-122.

3. Гринченко В.Т., Улитко А.Ф., Шульга Н.А. Механика связанных полей в элементах конструкций. Киев: Наукова думка, 1989. 279 с.

4. Еленицкий Э.Я., Дьяченко Ю.П. Применение метода начальных параметров к решению нестационарной задачи динамики для прямоугольной пластины ступенчатого сечения // Изв. вузов. Строительство. 1997. № 11. С. 13-18.

5. Паймушин В.Н. Соотношение теории тонких оболочек типа теории Тимошенко при произвольных перемещениях и деформациях // Прикладная механика и техническая физика. 2014. Т. 55, № 5. С. 135-149.

6. Рудницкий С.И., Шарапов В.М., Шульга Н.А. Колебания дискового биморфного преобразователя типа металл-пьезокерамика // Прикл. механика. 1990. Т. 26, № 10. С. 64-71.

7. Савин В.Г., Бабаев А.Э. Действие акустического импульса на плоскую электроупругую систему из биморфов // 1нформацшш системи мехашка та керуания: сб. Киев: Нац. техн. ун-т Украины. 2009. Вып. 3. С. 30-39.

8. Сеницкий Ю.Э. Метод конечных интегральных преобразований - обобщение классической процедуры разложения по собственным вектор-функциям // Изв. Саратовского ун-та. Сер. Матем., механика, информатика. 2011. № 3(1). С. 61-89.

9. Шляхин Д.А. Вынужденные осесимметричные колебания пьезокерамической тонкой биморф-ной пластины // Изв. РАН. МТТ. 2013. № 2. С. 77-85.

10. Янчевский И.В. Управление колебаниями изгиба круглого асимметричного биморфного пьезо-преобразователя с разрезными электродами // Проблемы машиностроения. 2012. Т. 15, № 2. С.37-43.

11. Dzyuba V.A., Steblyanko P.O. Use of splines in the calculation of deflections for plates of variable thickness. Science and Education a New Dimension. Natural and Technical Sciences, II(4). 2014; 32:41-47.

12. Wang Y., Xu R.Q., Ding H.J. Analytical solutions of functionally graded piezoelectric circular plates subjected to axisymmetric loads. Acta Mechanica. 2010(215);1-4:287-305.

FEFU: SCHOOL of ENGINEERING BULLETIN. 2019. N 2/39

Mechanics of Deformable Solids www. dvfu. ru/en/vestnikis

DOI: https://dx.doi.org/10.24866/2227-6858/2019-2-2

Ratmanova O.

OLESYA RATMANOVA, Assistant, Department of Building Mechanics and Structural Resistance, e-mail: olesya654@yandex.ru Samara State Technical University 244 Molodogvardeyskaya St., Samara, Russia, 443001

Non-stationary direct piezoelectric effect task for an asymmetric bimorph round plate with stepped variable thickness

Abstract: Cosideration is being given to the dynamic axisymmetric problem of direct piezoelectric effect for a circular two-layer transducer consisting of a metal substrate and a piezoceramic plate of smaller diameter. The oscillations of the electroelastic system are carried out due to the action of a mechanical load on its front surface, which is an arbitrary function of time and radial coordinate. The presented mathematical model of calculation allows to take into consideration the coherence of longitudinal-transverse oscillations of a structure that is asymmetric in height. The new analytical solution is obtained within the framework of the refined Tymoshenko theory by the method of finite integral transformations. The calculated ratios make it possible to choose the geometric dimensions of the bi-morphic plate, which most effectively allow to convert the mechanical load into the generated electri-

cal pulse. The considered designs are widely used as touch switches that regulate the operation of valves for various purposes.

Keywords: asymmetric bimorphic plate, direct piezoelectric effect, electroelasticity, dynamic load, closed solution.

REFERENCES

1. Abdikarimov R.A., Zhgutov V.M. Mathematical models of problems of nonlinear dynamics of viscoe-lastic isotropic plates and shells of smooth variable thickness (asymmetric cases). Civil Engineering Journal. 2010;8:47-55.

2. Vatulyan A.O., Rynkov A.A. A model of Flexural vibrations of piezoelectric bimorphs with split electrodes and its applications. Izv. Russian Academy of Sciences. MTT. 2007;4:114-122.

3. Grinchenko V.T., Ulitko A.F., Shulga N.Ah. Mechanics of related fields in design Elements. Kyiv, Sciences Dumka, 1989, 279 p.

4. Elnicki E.J., Dyachenko P.Y. Application of the method of initial parameters to the solution of unsteady problems of dynamics for a rectangular plate stepped section. Izv. Higher Educational. Construction. 1997; 11:13-18.

5. Paymushin V.N. The ratio of the theory of thin shells such as the theory of Tymoshenko with arbitrary displacements and deformations. Applied mechanics and technical physics. 2014(55);5:P.135-149.

6. Rudnitsky S.I., Sharapov V.M., Shulga N.A. Vibrations of disk bimorph transducer of the type metal-piezoceramics. Prikl. Mechanics. 1990(26);10:64-71.

7. Savin V.G., Babaev A.E. The Effect of acoustic pulse with flat elastic and electric system of bimorph Information system mechani Caruana. Kyiv, National Technica Un-t Ukraine. 2009;3:30-39.

8. Senitskiy Yu.E. Method of finite integral transformations is a generalization of the classical procedure of expansion in eigen vector-functions. Izv. Saratov Un-ta. Series. Mathematics, Mechanics, Informatics. 2011;3:61-89.

9. Shlyakhin D.A. Forced axisymmetric vibrations of a piezoceramic thin bimorph plate. Izv. Russian Academy of Sciences. MTT. 2013;2:77-85.

10. Yanchevskii I.V. Control of oscillation of the bending of circular piezoelectric transducer asymmetrical bimorph with split electrodes . Problems. Mechanical Engineering, 2012(15); 2:37-43.

11. Dzyuba V.A., Steblyanko P.O. Use of splines in the calculation of deflections for plates of variable thickness. Science and Education a New Dimension. Natural and Technical Sciences, II(4). 2014; 32:41-47.

12. Wang Y., Xu R.Q., Ding H.J. Analytical solutions of functionally graded piezoelectric circular plates subjected to axisymmetric loads. Acta Mechanica. 2010(215);1-4:287-305.