以局部徑向基底函數佈點法分析三維功能梯度壓電半導體問題

(1)

國立臺灣大學工學院土木工程學系碩士論文

Department of Civil Engineering College of Engineering National Taiwan University

Master Thesis

以局部徑向基底函數佈點法分析三維功能梯度壓電半導體問題

Three-dimensional Analysis for Functionally Graded Piezoelectric Semiconductors by the Local Radial Basis

Function Collocation Method

陸學賢

Hubert Hsueh-Hsien Lu

指導教授：楊德良博士 Advisor: Der-Liang Young, Ph.D.

中華民國 105 年 7 月

July, 2016

(2)

審定書

(3)

致謝

首先感謝指導教授楊德良老師非常有耐心的指導我，從一個什麼都不懂的大學畢業生，到現在對自己研究生生涯感到驕傲，老師不僅在學術方面給予指導，

甚至在做人處事上都相當關心學生我，也改正了我不少壞習慣。起初，參與壓電材料的數值模擬計畫時相當的不開心，那時覺得我自己身為一個水利組的學生就應該做和流體力學或水利工程相關的研究，但在老師的循循善誘下，我漸漸地珍惜這個能跨領域學習的機會，在碩一升碩二的暑假甚至還讓研究生的我能到斯洛伐克參加學術交流，非常感謝老師讓我加入這個研究團隊，往後學生也將更努力地學習。

接著我想感謝 Sladek 教授的建議和指導，還有其團隊成員 Peter，非常有耐心積極的指正我對學理上的錯誤，一來一往百餘篇的電子郵件讓我非常感激。另外特別感謝口試委員們提供給我許多寶貴的意見，讓我能在研究和論文上更精進更完整。此外非常感謝學長吳清森博士不管在學術或是在遇到問題時該如何面對上都給我許多意見和幫助，也感謝研究室學長姐們在我剛進團隊時非常有耐心的教我，從一個程式都不會到 MATLAB、C++都相當熟悉。感謝研究室的學弟們、

助理，Mahmoud 讓我練習英文口說還有亂扯紓解壓力，也感謝我同屆的研究室好麻吉們，祝大家順利畢業。最後感謝我的家人們，總是支持我一切決定，就像隱形的翅膀一樣讓我走到今天，在我徬徨失落的時候給我依靠，希望以後我也能成為你們信任安心的依靠，並且實現自己的夢想。

(4)

摘要

本研究主旨為藉由局部徑向基底函數佈點法分析三維梯度功能壓電半導體問題。局部徑向基底函數無網格數值方法已經被廣泛運用在工程與科學領域上，

由於在處理空間尺度相差甚異與無須數值積分的優勢下，因此將此數值方法運用在壓電材料的工程問題上。

壓電材料可以分為介電體與半導體兩類，不同於壓電絕緣體，由電子密度及電流所組成的電守恆式需被額外用以描述壓電半導體的現象，這也加深了彈性位移及電場間相互作用分析的複雜度。此研究以局部徑向基底函數佈點法來求解存在非定常數的偏微分方程，在物理場中的空間變化以多元二次曲面徑向基底函數近似；時變性的微分方程系統問題以 Houbolt 有限差分法求解。

以有限元素法之相對應結果來驗證局部徑向基底函數佈點法之結果，且分析在不同載重情形下分析的樑所產生的力學反應、電場、電流場間互相的關係。此外，此研究也分析梯度參數及初始電子密度所產生之影響。最後，暫態分析也在此研究的範疇中。

關鍵詞: 局部徑向基底函數佈點法、功能梯度材料、壓電半導體、壓電效應、智能材料

(5)

Abstract

This thesis presents three-dimensional analysis of functionally graded piezoelectric semiconductor by the local radial basis function collocation method (LRBFCM). The LRBFCM is a commonly-used meshless numerical method in the field of engineering and sciences. On account of the advantages of addressing the problems with much different length scales in three dimensions and circumventing numerical quadrature, the LRBFCM is investigated and applied in the problems of piezoelectric materials.

Piezoelectric materials can be divided by dielectrics and semiconductors.

Unlike piezoelectric dielectric materials, the conservation of charge which is composed of electron density and electric current is additionally considered to depict the phenomenon for piezoelectric semiconductors. This will complicate our analyzing the mutual coupling of elastic displacements and electric fields. For the solution of the set of partial differential equations with non-constant coefficients the LRBFCM is proposed in this work. The spatial variations of all physical fields are approximated by the multiquadric radial basis function. For time dependent problems a resulting system of ordinary differential equations is solved by the Houbolt finite difference scheme as a time stepping method.

The presented LRBFCM method is verified by using the corresponding results

(6)

obtained by the finite element method. The effect of various loading scenarios is then considered in the numerical examples to analyze the mutual properties of the mechanical responses, electrical fields, and electrical current field. The influence of material parameter gradation and initial electron density is then investigated. The transient analysis is also analyzed.

Keywords: local radial basis function collocation method (LRBFCM), functionally

graded materials, piezoelectric semiconductors, piezoelectric effect, smart materials

(7)

List of Figures

Fig.2.2.1 The illustration of a computation point with its local points respect to the local influence area  within the computational _i domain .

10

Fig.2.4.1 The number of the nearest points with respect to every computation point by choosing NL 5 in the two-dimensional uniform point distribution case.

16

Fig.2.4.2 The selection of points within the fixed radius in the two- dimensional non-uniform point distribution case.

17

Fig.2.4.3 The cross-shaped selection 18

Fig.3.2.1 The illustration of the relations between the energies and the piezoelectric effects.

23

Fig.3.2.1 Configuration of the direct piezoelectric effect. (G:

galvanometer)

24

Fig.3.4.1 Configuration of two-layered composites and FGMs. 26

Fig.4.5.1 The geometry of the beam. 50

(11)

Fig.4.5.2 The distribution of the points. 51

Fig.4.5.3 The illustration of the boundary conditions. 52

Fig.4.5.1.1 Variation of vertical displacement for a line along x₁ and located at x ₂ 0_and x₃h/ 2 in the static analysis.

53

Fig.4.5.1.2 Variation of electric potential for a line along x₁ and located at

2 0

x  _and x₃ h/ 2 in the static analysis.

53

Fig.4.5.1.3 Variation of electron density for a line along x₁ and located at

2 0

x  _and x₃ h/ 2 in the static analysis.

54

Fig.4.5.2.1 Variation of vertical displacements along the line



x x1, 2 0,x3h/ 2



for different values of M₀ in the beam

under static load T ₃ 100Pa.

55

Fig.4.5.2.2 Variation of electric potential along the line



x x1, 2 0,x3h/ 2



55

Fig.4.5.2.3 Variation of electron density along the line



x x1, 2 0,x3h/ 2



56

(12)

Fig.4.5.3.1 Variation of vertical displacement along the line



x x1, 2 0,x3h/ 2



for different values of the gradation

parameter  in the FGM beam under the static load

3 100Pa T  .

57



x x1, 2 0,x3h/ 2



3 100Pa T  .

57



x x1, 2 0,x3h/ 2



3 100Pa T  .

58

Fig.4.5.4.1 Variation of vertical displacement along the line



x x1, 2 0,x3h/ 2



for the three combinations of the

gradation of material coefficients in the FGM beam under static load T ₃ 100Pa.

59



x x1, 2 0,x3h/ 2



for the three combinations of the 60

(13)



x x1, 2 0,x3h/ 2



60

Fig.4.5.5.1 Time evolution of vertical displacement along the line



x x1, 2 0,x3h/ 2



by the CT and BT methods for the first-

order temporal partial derivative under static and impact load

3 100Pa T  .

62

Fig.4.5.5.2 Time evolution of electric potential along the line



x x1, 2 0,x3h/ 2



by the CT and BT methods for the first-

order temporal partial derivative under static and impact load

3 100Pa T  .

62

Fig.4.5.5.3 Time evolution of vertical displacement along the line



x x1, 2 0,x3h/ 2



gradation of material coefficients in the FGM beam under static and impact load T ₃ 100Pa.

63

(14)

Fig.4.5.5.4 Time evolution of electric potential along the line



x x1, 2 0,x3h/ 2



gradation of material coefficients in the FGM beam under static and impact load T ₃ 100Pa.

64

(15)

List of Tables

Table.2.3.1 List of commonly-used RBFs 13

Table.2.3.2 The differential formulation of the MQ-RBF 14

(16)

Chapter 1 Introduction

1.1 Motivations and Objectives

On account of the development of computer science and technology, numerical analyses have been frequently utilized and even substituted for some experiments.

There are advantages by numerical analyses such as cheap and enable some problems which do not exist analytical solutions to be solved and analyzed. The numerical methods generally can be categorized as two types, mesh-dependent methods and meshless methods. Mesh-dependent methods have been developed and commonly used in scientific research and engineering applications. The following four mesh-dependent methods are most commonly used methods, the finite difference method (FDM), the finite volume method (FVM), the finite element method (FEM), and the boundary element method (BEM). However, the mesh-dependent numerical schemes still exist challenges and problems on account of inevitable burdensome tasks such as mesh generation and numerical quadrature especially for multi-dimensional problems and irregular domains. In order to avoid those problems, various meshless methods have been developed in recent years. They have become more and more popular due to the ease of implementation and the flexibility of generation of computational nodes which can circumvent the problems of inaccuracy near where gradients of variables are high.

(17)

However, there are the stability and accuracy of meshless methods which should be concerned and analyzed. In this thesis, we focus on one of the meshless numerical methods, the local radial basis functions collocation method (LRBFCM), which will be introduced in subsection 1.1.2.

1.1.1 Mesh-dependent numerical methods

With the development of the computer technology in recent years, it has been more feasible and efficient to utilize numerical methods to simulate and analyze engineering problems. The finite difference method (FDM), the finite volume method (FVM), and the finite element method (FEM) have been developed to a robust and effective methods and widely used in engineering problems. The finite difference method is based on Taylor series expansion to approximate the derivative which accuracy is dependent on how many terms we utilize. However, it is complicated when we solve problems with irregular domain due to difficulty of construction of orthogonal mesh. Additionally, we should refine the mesh and add more terms of derivative to obtain more accurate results which is very time-consuming. In order to reduce the dependency of the meshes, researchers have developed a powerful alternative numerical scheme, the boundary

element method (BEM), to substitute for other mesh-dependent numerical methods. By the Green’s function, the BEM can reduce one dimensionality of the problems. The

BEM discretizes the computational 3D domain of the surface instead of the whole

(18)

domain. Therefore, it provides more flexibility and relatively decreases the dependency of meshes. However, the limitation of the fundamental solutions or free-space Green’s functions restricts engineers to apply the BEM to some problems. Due to the above difficulties of applying the mesh-dependent numerical methods, researchers have been paying attention on the development and improvement of meshless numerical methods.

1.1.2 Meshless numerical methods

In order to circumvent numerical quadrature and mesh generation, various meshless or mesh-free numerical methods have been developed such as the smoothed particle hydrodynamics (SPH) [1], the multiquadrics collocation method (MQ) [2]-[6], the method of fundamental solutions (MFS) [7]-[10], the method of particular solutions (MPS) [11][12], the method of approximate particular solutions (MAPS) [13][14], the differential quadrature method (DQ) [15][16], the boundary particle method [17], and the finite point method (FPM) [18]. The above numerical methods are classified as global-type methods. Global-type methods generally utilize the discretization of all collocation points within the global domain. Subsequently, the problems such as the ill- conditioned and dense resultant interpolation matrix are inevitable. In order to improve the efficiency of computation and deal with large-scale problems, researchers have been developing various localization methods for the corresponding global-type numerical methods.

(19)

The localization technique allows the global-type numerical methods to approximate the solution of partial differential equations (PDEs) with few local points and then the individual relations constitute the sparse resultant interpolation matrix instead of the dense matrix. That enables computation procedures to utilize the solvers of inversion of sparse matrix and save substantial computational time and memory loading. The MQ, the MAPS, and the DQ are respectively improved and localized as the localized multiquadric method (LMQ) [19], the localized method of approximate particular solutions (LMAPS) [20]-[22], and the localized differential quadrature (LDQ) [23][24]. In addition, there are other localized meshless numerical methods such as the compactly supported radial basis functions [25], the local radial basis function collocation method (LRBFCM) [26]-[29], and the meshless local Petrov–Galerkin (MLPG) method [30][31]. This research will focus on the LRBFCM, and utilize it to conduct three-dimensional analysis for functionally graded piezoelectric semiconductors.

1.2 Organization of the thesis

In order to analyze the functionally graded piezoelectric semiconductor problems by the local radial basis function method, we divide the thesis into five chapters and the brief introduction of them is shown as follows:

(20)

Chapter 1 Introduction

The motivation of this thesis is presented, and the difference between the mesh- dependent methods and meshless methods is also discussed in this chapter. In addition, the development of local meshless methods with respect to corresponding global-type meshless methods is mentioned.

Chapter 2 The Local Radial Basis Function Collocation Method

The meshless methods, the radial basis function collocation method (RBFCM) and the local radial basis function collocation method (LRBFCM), are introduced in this chapter.

Chapter 3 Piezoelectricity

The piezoelectricity is introduced and can be divided by five sections, historical overview, principles of piezoelectric effect, applications in civil engineering, functionally graded materials, and the constitutive equations of piezoelectric materials.

Chapter 4 The Local Radial Basis Function Collocation Method for

Functionally Graded Piezoelectric Semiconductor

The LRBFCM is utilized to analyze the functionally graded piezoelectric

(21)

semiconductor problems in this chapter. The results of the LRBFCM and those of the finite element method by the commercial software, COMSOL, are compared and in good agreement. The effect of functionally graded properties is also illustrated and compared. The transient analyses are also analyzed in this chapter.

Chapter 5 Conclusions and Future works

Conclusions and future works are summarized and proposed in this chapter, respectively.

(22)

Chapter 2 The Local Radial Basis Function Collocation Method

This chapter will report the numerical tool, LRBFCM, which is utilized to deal with the piezoelectric sensor problem in this thesis. In the first section, we will illustrate the concept and the approximation procedure of the LRBFCM thoroughly. Then, we will introduce the RBFs which play an important role in the LRBFCM in the second section. Although the LRBFCM can overcome the drawbacks of the RBFCM, it is still a developing meshless numerical method whose stability and the accuracy should be further investigated. The accuracy and the stability in the LRBFCM is strongly depended on the selection of the supported local nodes and the value of the shape parameters. For the purpose of promoting the accuracy and the stability of the approximation results, we will introduce some kinds of methods for choosing the supported local nodes in the section 3 and the shape parameters in the section 4, respectively. In the last section, we provide the normalization technique to improve the multi-scale domain problems.

(23)

2.1 The radial basis function collocation method

The multiquadric (MQ) scheme was first proposed by Hardy [2] in 1971. Then, the radial basis function collocation method (RBFCM), which is a modified MQ scheme, was developed by Kansa. The RBFCM [3][4] requires a linear combination of radial basis functions with regard to all the computation points within the computational domain  to approximate any given variable denoted by . Let all the computation points be defined as x_i ,i [1,N], where N is the total number of global points.

Then the given variable at any computation point within the computational domain can be approximated by the following term

   

1

, ,

N

k

k k

t f

 





 

x x x x (2.1-1)

where f is the radial basis function,  is the weighting coefficient to be determined, and x x is the Euclidean distance between the global points. In order to evaluate  _k the weighting coefficients, Eq. (2.1-1) at each computational point should be enforced.

Then, the system of the algebraic equations becomes

 f ,

  (2.1-2)

where 

 

x1,t , ,



x_N,t



^T,



1, ,_N



^T,_and

(24)

     

1 1 1 2 1

2 1 2 2 2

1 2

.

N N

N N N N

f f f

    

 

  

 

  

 

    

 

f

x x x x x x

To implement any given or operator ^

 

, the approximate summation equation of the variable in Eq. (2.1-1) becomes



^

 

^,^t



^

 

^,

 x  f  (2.1-3)

By solving Eq. (2.1-3), we obtain the weighting coefficient with respect to every computation point within the global domain. Then substituting the weighting coefficients into Eq. (2.1-1), we get the approximate values of the variable within the computational domain.

2.2 The local radial basis function collocation method

Due to time-consuming computation, large memory loading, dense and ill- conditioned resultant interpolation matrix, and sensitive shape parameter, the localization technique, the local radial basis function collocation method, was developed [26]. The local radial basis function collocation method (LRBFCM) enables engineers to analyze large-scale realistic problems and more efficiently compute partial differential equations by sparse solvers. Contrary to the RBFCM, LRBFCM utilizes the approximation of ^

 

^x^,t with respect to the point x_i,i[1,N] which is supported

(25)

by NLfrom the local influence area  as _i

 

^,



^,



^,

1

, ,

NL

ik i i k i

k

t f k

  











x x x x (2.2-1)

where _{i k}, (i[1,N];k[1,NL]) are the weighting coefficients, f is the RBF, and

,

 i k

x x is the Euclidean distance between two points. In order to evaluate the weighting coefficients, Eq. (2.2-1) at each computational point from  should be _i enforced. Then, the system of the algebraic equations becomes

i  fi i,

  (2.2-2)

where _i ^



x_i,1, ,t



,



x_{i NL}, ,t



^^T,_i  _i,1, ,_i,_N_L^T , and

     

,1 ,1 ,1 ,2 ,1 ,

,2 ,1 ,2 ,2 ,2 ,

, ,1 , ,2 , ,

.

i i i i i i NL

i

i NL i i NL i i NL i NL

f f f

    

 

    

 

 

 

    

 

f

x x x x x x

Fig.2.2.1 The illustration of a computation point with its local points respect to the local influence area  within the computational domain_i ^.

x

,

xi kⁱ

x





i

(26)

Since f_i is invertible, the weighting coefficient can be determined as

1 , [1, ].

i fi^ i i N

  (2.2-3)

Additionally, Eq. (2.2-1) can be rewritten as

 

^,t i

 

i

 x  F x  (2.2-4)

by introducing the row-vector ^Fⁱ

 

^x ^^_^f



^x^^xⁱ^,1



^, ^,^f



^x^^xⁱ^,^NL



^_^. ^To

implement any given or operator ^

 

, the approximate summation equation of the variable in Eq. (2.2-1) gives

        ^ ^{ } ^

1

,

, ,

NL

i k i k i i

k

t f

    









 ^F ^

x x x x (2.2-5)

where

       

i

i i i

 F  F _

x x x x .

For convenience, we define the row vector (1NL)

   

¹ ,1, ,2, , ,NL

i  Fi i fi^  si si si 

s x , and then Eq. (2.2-5) becomes

 



i^,t



i i i ^[1,N^].

  x s   (2.2-6)

Subsequently, we transform the local system to the global system as shown in Eq. (2.2- 6).

 



i^,t



i i i i ^[1,N^],

  x s S  (2.2-7)

where

(27)

     

,1 ,2 ,N , 1 2

,

, , , , 0, , and , , , , , , .

,

j i T

i i i i i j N

i j j i

S S S S t t t

s

   



 

 

      

S x

x x x

x 

Finally,

 

 

 

 

 

 

 

1 1

2 2

N

, ,

N, t t

t

 

  

 

   

   

      S S S



 

 x

x

(2.2-8)

where  

1,1 1,2 1,N

2,1 2,2 2,N

1 2 N

N,1 N,2 N,N

, , , ^T ,

S S S

 

 

  

 

 

S S S S

and we can obtain the solution of  by solving the linear system in Eq. (2.2-1).

Additionally, we define S^ij as the S whose given operator is the second-order partial spatial derivatives with respect to x_i and xj in this thesis; that is,

 

²^/ x xi j

 

     where ,i j represent the indices of Cartesian spatial dimension

number.

2.3 Radial basis function

The radial basis functions (RBFs) have been frequently utilized and applied to approximate scattered data in recent years. By the Euclidean distance r, the LRBFCM can form the discretization equation. However, there are several types shown in the Table.2.3.1 where the shape parameter is denoted by c. The decision of appropriate shape parameter is very cumbersome especially in practical engineering problems. We will elucidate that in the following section.

(28)

Table.2.3.1. List of commonly-used RBFs

Name of RBFs Formulation

Multiquadric (MQ) r²+c²,c  0

Inverse multiquadric (IMQ) 1/ r²+c²,c  0 General Multiquadric (GMQ)



^r²⁺^c²



^1/ⁿ^,^c^^{0 (}ⁿ^^1,3,5...)

Gaussian (GA) ^exp



^^cr²



^,^c^⁰

Polyharmonic Splines (PS) of order m in 2D

 

2^mln

r r

Polyharmonic Splines (PS) of order m in 3D

2m 1

r ^

Franke [32] compared several methods in 1982 from the characteristics, such as

accuracy, sensitivity to parameters, timing, storage requirements, and so on. Among the

methods, the multiquadric RBF (MQ-RBF), that is Hardy’s multiquadric method, is one of the most frequently-used methods due to the characteristics of accuracy and stability in the LRBFCM. Therefore, we adopt the MQ-RBF and use it to form the approximate equations in this thesis. To utilize the radial basis function f as the MQ-RBF with

respect to computation points x and x^k, the radial basis function can be expressed as

(29)

 

k k

f  f r (2.3-1)

where r xx^k  (x₁x₁^k)² (x_d x_d^k)² is the Euclidean distance in the Cartesian coordinate system and dis the number of spatial dimension. The differential formulation in the Cartesian coordinate system is shown in Table.2.3.2.

Table.2.3.2. The differential formulation of the MQ-RBF

2-dimensional formulation 3-dimensional formulation

f



^r²^^c²



¹²



^r²^^c²



¹²

1

f x



 

1 1

1

2 2 2

x xk

r c





 

1 1

1

2 2 2

x xk

r c





2

f x



 

2 2

1

2 2

k

2

x x

r +c



 

2 2

1

2 2 2

x xk

r c





3

f x



 0

 

3 3

1

2 2 2

x xk

r c





2

2 1

f x



 

2

1 1

1 3

2 2 2 2 2 2

1 x x^k

r c r c

 

 

 

2

1 1

1 3

2 2 2 2 2 2

1 x x^k

r c r c

 

 

2

2 2

f x



 

2

2 2

1 3

2 2 2 2 2 2

1 x x^k

r c r c

 

 

 

2

2 2

1 3

2 2 2 2 2 2

1 x x^k

r c r c

 

 

2

2 3

f x



 0

 

2

3 3

1 3

2 2 2 2 2 2

1 x x^k

r c r c

 

 

2

1 2

f x x



 

  

 

1 1 2 2

3

2 2 2

k k

x x x x

r c

 





  

 

1 1 2 2

3

2 2 2

k k

x x x x

r c

 





(30)

2

1 3

f x x



  0

  

 

1 3

2

k k

1 3

3 2 2

x x x x

r +c

 



2

2 3

f x x



  0

  

 

2 3

2

k k

2 3

3 2 2

x x x x

r +c

 



2.4 Local influence area

The numerical accuracy and stability of LRBFCM are highly dependent on the local influence area. Inappropriate local influence will possibly cause the ill- conditioned resultant interpolation matrix. Therefore, there have been many fashions of selection of local influence area developed to deal with the above problems. The most popular three fashions are the selection of fixed number of the nearest points, the selection of points within the fixed radius, and the cross-shaped selection. Due to the easy implementation, we adopt the selection of fixed number of the nearest points in this research.

First, the selection of fixed number of the nearest points is the easiest and most common selection method due to its characteristic of simple coding. By utilizing this method of selection, we should determine the number of the nearest points with respect to every computation point. To exemplify, the two-dimensional uniform point distribution case is shown in Fig.2.4.1. It is obvious that the selection in the boundary and corner is asymmetric and unbalanced. This phenomenon will become more

(31)

significant in the non-uniform point.

Fig.2.4.1 The number of the nearest points with respect to every computation point by choosing NL 5 in the two-dimensional uniform point distribution case.

Second, the selection of points within the fixed radius is also a commonly-used method which should be given a radius to determine the local points with respect to every computation point. From Fig.2.4.2, computational inefficiency possibly occurred in the relatively dense region and insufficient local points are illustrated. If the non- uniform point distribution cases or the problems with the much different length scales in three dimensions are applied by the selection of fixed number of the nearest points, this phenomenon will be more significant which increases the opportunity of

(32)

occurrence of numerical instability and ill-conditioned matrix.

Fig.2.4.2 The selection of points within the fixed radius in the two-dimensional non- uniform point distribution case.

Third, the cross-shaped selection [33] as shown in Fig.2.4.3, which utilizes the provided shape of local influence area and provided number of local points with respect to computation point to determine the local point within each local domain, provides more efficient computation and lowers the opportunity of occurence of numerical instability and ill-conditioned matrix.

(33)

Fig.2.4.3 The cross-shaped selection.

2.5 Shape parameter

The decision of optimal shape parameter c plays a crucial role in the stability and accuracy of the LRBFCM for different numerical applications. It strongly depends on geometry of global domain and types of point distributions.

In general, shape parameter is located on the specific range. If c is too small, it could result in the singularity; however, if c is too big, the influence of radius will be decreased by the shape parameter. It is very difficult to determine optimal shape parameter when we deal with the problems which the length scales in three dimensions

(34)

are much different. We will introduce and elucidate a normalization technique to assure optimal shape parameter to be in the same order in the next section.

2.6 Normalization technique

Many normalization techniques have been developed due to the difficulty of determination of optimal shape parameter. Generally, normalization techniques can be divided by two parts, normalized distance and normalized shape parameter.

2.6.1 Normalized distance

One of the most popular normalization techniques is to normalize the distance as

2

, , , ,

, ,

1 ,

d i k d i m

Nd i k m

d d i

r L

x x



 

  



 (2.6.1-1)

where r_{i k m}_{, ,} is the distance between the collocation points within the local influence area  , _i x_{d i k}_{, ,} denotes the position of the kth local point within the local influence area  in the_i ^d^th dimension, Nd is the total number of spatial dimension, and the maximum distance between collocation points in all the dimensions within the local influence area  is defined as _i L . _{d i}_,

2.6.2 Normalized shape parameter

Another normalization technique is developed in [26] and it is to normalize the given shape parameter c associated with each local influence area _i as

(35)

i i,

c  c L (2.6.2-1)

where c is the normalized shape parameter and_i L is the maximum distance between _i the computation points within the local influence area _i. By normalizing the shape parameter, the range of optimal shape parameter could be narrowed.

Due to the efficiency and improvement, we adopt the second normalization technique, normalized shape parameter, for simulation of piezoelectric problems in this thesis. Additionally, there are the cross partial derivative terms in our problems (see Eqs. (4.2-19)-(4.2-23) in Chapter 4), so the local influence must cover more local points rather than local points of cross-shaped selection. As a result, normalization technique of shape parameter is necessary to solve piezoelectric problems.

(36)

Chapter 3 Piezoelectricity

3.1 Historical overview

The story of the piezoelectric materials starts in 1880, when the Pierre and Jacques Curie discovered [34] that several natural materials, including quartz and Rochelle salt, exhibited a special property. The Curie brothers demonstrated that if the specially prepared materials were imposed a mechanical stress, an electric output was produced.

They showed this coupling by measuring the charge induced across electrodes placed on the material when it was imposed an applied mechanical deformation. They defined this effect the piezoelectric effect. The name comes from a Greek word for squeeze – piezein. Few years later it was demonstrated by Gabriel Lippmann [35] that

piezoelectric materials also exhibited the reciprocal property; namely, a mechanical strain was induced when an electric field was applied to the materials.

However, the coupling weak, which means the amount of electrical signal produced by applied mechanical deformation was small, limited the application due to the lack of precision instrumentation. The first engineering application was developed to locate submarines, which is the basis of sonar, until World War I. The piezoelectric materials were widely used in sonar during World War II and developments in electronics also stimulate different uses of piezoelectric materials, such as electronic

(37)

oscillators and filters. On account of the increasing need for better piezoelectric materials, the synthetic materials were developed to exhibit better piezoelectric properties. To exemplify, the early synthetic piezoelectric material, Barium titanate

(BaTiO ) , is superior to quartz crystals in piezoelectric and thermal properties. In the 3

1950s and 1960s the most widely used piezoelectric material, lead-zirconate-titanate (PZT), was developed and motivated more applications. Nowadays, piezoelectricity is utilized everywhere. For example, motion and force sensors, the airbag, accelerometers, and atomic force microscopes (AFMs). The application in civil engineering will be presented in section 3.3.

3.2 Principles of piezoelectric effect

The piezoelectric effect can be divided by two types as shown in Fig.3.2.1. The first is the direct piezoelectric effect which depicts piezoelectric materials transform the applied mechanical strain into the electric output. The second is the converse effect which describe mechanical strain energy is produced by an applied electrical potential on piezoelectric materials.

(38)

Fig.3.2.1. The illustration of the relations between the energies and the piezoelectric effects.

In general, we will utilize the direct piezoelectric effect to be a sensor and utilize the converse piezoelectric effect to be an actuator. From Fig.3.2.1, if the piezoelectric material generates a positive electric field by the applied tensile stress (see (a)), the applied compressive load will generate a negative electric field (see (b)). Furthermore, the converse piezoelectric effect also exhibits this phenomenon. If the positive electric filed is imposed on the piezoelectric material and generates a contraction of the material, an expansion of the material will be generated by the applied negative electric field.

Electrical Energy

Mechanical Energy

Direct Piezoelectric

Effect

Converse Piezoelectric

Effect

(39)

Fig.3.2.1 Configuration of the direct piezoelectric effect. (G: galvanometer)

3.3 Applications in civil engineering

In recent decades, there have been many piezoelectric applications in civil engineering due to advances of science and technology. The so-called smart materials, which can be significantly changed by the applied stimuli such as stress or electric output, were presented [36]. This thesis also introduces the related applications of smart materials including piezoelectric materials to civil and mechanical infrastructure systems. Piezoelectric materials have been playing an important role in the vibration control of structures [37]. In addition, piezoelectric materials are also applied [38] in structural damping mechanism by passive electrical circuits, while piezoceramics [39]

are used in various forms for active control of structural vibration and are applied in civil structures such as beams and steel frames. The embedded piezoelectric wafer active sensors (PWAS) [40] perform an important function in structural health

(40)

monitoring (SHM) by exciting and detecting tuned Lamb waves. Since oil resources have been gradually depleted, more and more researchers devote themselves to alternative technology of power harvesting and piezoelectric materials are considered a feasible and renewable resources. Erturk [41] introduced and analyzed the energy harvesting of piezoelectric materials from moving load excitations and surface strain fluctuations in civil infrastructure system.

3.4 Functionally graded materials

In recent decades, a novel advanced materials, functionally graded materials (FGMs), have been gradually attached importance in various engineering applications.

The characteristic of FGMs is that continuous and gradual variation of material properties over the spatial coordinates. The difference of structure between two-layered composites and FGMs is shown in Fig.3.4.1. FGMs have shown advantage of better performance over multilayered composites. In general, conventional multilayered composites suffer from abrupt changes of material properties [42] at the interface between contiguous layers of composites which results in problems such as delamination and large residual thermal stresses. In contrast, FGMs are utilized to reduce the stress concentration and the fracture toughness [43][44].

From the above advantages, the properties of FGMs provide prospect for

(41)

applications of piezoelectric materials. The attributes of low thermal expansion coefficient, low dielectric constant, high toughness, high strength, and increase of gradient of the material properties [45] can be applied to extend the lifetime, and improve reliability of piezoelectric structures [46]. In general, the grading variation is defined by power-law, sigmoid, or exponential function. In this work the exponential function is adopted and we will elucidate the definition of grading properties in Eq.

(4.2-24).

Fig.3.4.1 Configuration of two-layered composites and FGMs.

3.5 The constitutive equations of piezoelectric materials

From thermodynamic considerations, the constitutive equations of piezoelectric materials can be categorized into the following four forms [47] with

respect to their independent variables:

  D form:

c

D

kD

   

E    k  

^

D

^(3.5-1)

(a) Two-layered composites (b) FGMs

(42)

  E form:

s

E

dE

   

D d    h E

^ ^(3.5-2)

  D form:

s

D

gD

   

E   g   

^

D

^(3.5-3)

  E form:

c

E

eE

   

D e    h E

^ ^(3.5-4)

where  , , D, and E denote the stress, strain, electric displacement, and electric field, respectively. The stiffness constants at a constant electric displacement and a constant electric field are c^D c^E, respectively. The compliance constants s^D s^E are at a constant electric displacement and a constant electric field, respectively. The dielectric constants (permittivity) at a constant strain and a constant stress are h^ h^ , respectively. The impermittivity constants ^  are at a constant strain and a ^ constant stress, respectively. The constants k , d , g , and e are piezoelectric constants. To exemplify, the  ^{ E} form of the constitutive equation in Eq. (3.5-4) can be rewritten in tensor as

以局部徑向基底函數佈點法分析三維功能梯度壓電半導體問題

國立臺灣大學工學院土木工程學系 碩士論文

Department of Civil Engineering College of Engineering National Taiwan University

Master Thesis

以局部徑向基底函數佈點法分析 三維功能梯度壓電半導體問題

Three-dimensional Analysis for Functionally Graded Piezoelectric Semiconductors by the Local Radial Basis

Function Collocation Method

陸學賢

Hubert Hsueh-Hsien Lu

指導教授：楊德良 博士 Advisor: Der-Liang Young, Ph.D.

中華民國 105 年 7 月

July, 2016

審定書

致謝

摘要

Abstract

Table of Contents

List of Figures





















































List of Tables

Chapter 1 Introduction

1.1 Motivations and Objectives

1.1.1 Mesh-dependent numerical methods

1.1.2 Meshless numerical methods

1.2 Organization of the thesis

Chapter 2 The Local Radial Basis Function Collocation Method

2.1 The radial basis function collocation method

   



 









     

     

     

 



 



 

2.2 The local radial basis function collocation method

 

 















     

     

     

國立臺灣大學工學院土木工程學系碩士論文

以局部徑向基底函數佈點法分析三維功能梯度壓電半導體問題

指導教授：楊德良博士 Advisor: Der-Liang Young, Ph.D.

        ^ ^{ } ^