Application of the grain boundary formulation and image processing-based algorithm in micro-mechanical analysis of piezoelectric ceramic

Abstract

The main subject of this paper is the micro-mechanical analysis of a piezoelectric ceramic. In micro-mechanical analyses, it is very important to have knowledge about the real and natural micro-structure of materials. Therefore, the barium titanate powder was prepared using the solid-state technique, and pellets and beams were manufactured by uniaxial and isostatic pressing. The boundary element method (BEM) is used in order to be combined with three different grain boundary formulations for investigation of micro-mechanics and crack nucleation and evaluation in piezoelectric ceramic. In order to develop a numerical programming algorithm, suitable models of polycrystalline aggregate have to be discretised for the BEM analysis. Hence, original comprehensive algorithms are designed on the basis of image processing methods. Several assumptions are made to model the grain boundary in micro-scale. In the first step, before having any cracks, the traction equilibrium and displacement compatibility are governing equations. When the onset micro-crack starts to initiate, one mixed-mode potential based cohesive law is applied to model grain boundaries and investigate the intergranular crack nucleation and evolution. Upon interface failure, a frictional law is utilised in order to study separation, sliding or sticking between micro-crack surfaces. Several numerical experiments on barium titanate polycrystalline aggregate are presented to show the effectiveness of image processing-based discretisation algorithms and grain boundary formulation in micro-mechanical analysis.

Keywords

Micro-mechanical analysis image processing-based algorithm piezoelectric ceramic grain boundary formulation boundary element method intergranular fracture

1. Introduction

Nowadays, very complex numerical simulations for nano-scale to macro-scale materials take a lot time and need powerful computational machines. For this reason, proposing new methods to improve computational procedures can lead to more efficient numerical tools.

The boundary element method (BEM) is one of the most diffused computational methods in many areas of engineering and science, including fracture mechanics, fluid mechanics, electromagnetism, acoustics and geology [1 –5]. In fact, this method allows in a very simple way to both save time and perform calculations: in order to use the BEM, one needs to discretise only the boundary of the system, so the computational efficiency can be increased significantly [6,7].

In order to utilise the BEM, it is essential to use the traction fundamental solution and the displacement fundamental solution, which can be obtained for anisotropic bodies by using two main methods. The first method is based on Stroh’s formulation [8], and it is utilised in the current study, and the second one refers to Lekhnitskii’s formalism [9].

The micro-mechanical analysis, via computational procedures, allows one to understand how much grain size and orientations can influence the crack initiation or fracture criteria in composites. Secondarily, it is useful to determine the effective mechanical properties of materials. To this aim, different methods have been proposed [9 –13]. The investigation of effective properties of materials is mostly used for composites that contain cracks, holes or rigid bodies. Some composites also consist of a matrix phase. So, it is necessary to consider a unit cell with periodically distributed inclusion and then analyse it by using the BEM [10] or finite element method (FEM) [14] for evaluating the overall material properties of the composite. This method has been used in the case of piezoelectric ceramics [11]. However, they did not consider the cohesive law and pre-existing cracks in their modelling, which play important roles in determining the effective behaviour of materials.

In order to predict the fracture behaviour in a composite, different cohesive law models or friction phenomena can be used [15 –25]: via cohesive law models, it is possible to predict both intergranular and transgranular crack initiations. In other words, a cohesive law is suitable in order to model interfaces between materials with different mechanical properties. Depending on the possibilities to predict fracture and crack initiation, there are two different kinds of cohesive law models, the potential-based law [16 –18] and the linear law [19 –22].

Crack nucleation and propagation in polycrystalline materials can be studied in cohesive models by different numerical tools such as the BEM or FEM [26,27]. Despite the fact that some of the cohesive laws consider fracture criteria for the first and second mode separately, in [28,29] mixed-mode relations to predict the location of the cohesive zone in polycrystalline are presented.

The BEM has been used also in micromechanics and multistate modelling [30 –32]. In some of these works, materials are modelled only in micro-scale and neither cohesive law nor averaging theory or nonlocal theories were used. However, some basic concepts of micro-mechanical analysis, such as representative volume element or grain boundaries, have been used in the earliest works in the area of polycrystalline analysis, where a periodic arrangement of hexagonal grains was used to model the polycrystalline material at the micro-mechanical level. Symmetry considerations allowed Fotiu et al. [31] to restrict the models to representative unit cells with special boundary conditions. Sfantos and Aliabadi [33], Benedetti and Aliabadi [34] and Gulizzi et al. [35] have done a comprehensive research in using the BEM in the micro-mechanical analysis. Their works contain modelling of grain boundaries as well as generating a numerical algorithm.

For studying crack initiation and evaluation in a microstructure, suitable models have to be generated and prepared. The Poisson–Voronoi tessellation method is extensively used in order to generate a suitable model of polycrystalline materials [36,37]. However, it is clear that this method provides artificial microstructures and, for studying the scanning electron microscopy (SEM) images of specimens that are prepared in the laboratory, it does not have suitable solutions.

In order to have a bridge between the micro and macro analysis, averaging theorems have been suggested. It is worth mentioning that averaging theorems can be considered as a bridge in order to transfer all of the quantities from micro-scale to macro-scale [38 –40]. In the current work the averaging theorems proposed in [41] are used. For a more rigorous approach, the reader can refer to, for example, [42 –45].

This paper deals with the micro-mechanical analysis of failure in piezoelectric ceramic. In order to have knowledge about the micro-structures of ceramics, barium titanate pellets and beams have been sintered in different temperatures. Afterward, the SEM images of pellets that were sintered at different temperatures were prepared to achieve a reliable model of micro-structures. However, in order to perform the numerical analysis, these SEM images have to be analysed and discretised precisely for having geometrical information about grains boundaries. It is worth stressing that an algorithm, based on image processing methods, is designed for the first time, and we establish if it is able to discretise both the natural and artificial microstructure of materials. The numerical algorithm for utilising the BEM for every grain is developed, taking into account mechanical, electrical and piezoelectric properties. Then, three different statuses are considered to have accurate information about crack nucleation and evaluation in the grain boundary in micro-aggregate. Hereafter, suitable algorithms are created to combine the BEM with the cohesive-frictional grain boundary micro-mechanical model. An adaptive matrix is based on a numerical programming code able to combine the mentioned assumptions and change its elements according to the status of grain boundaries via considering mechanical, electrical and piezoelectric properties.

For some interesting applications of damage theories that can be considered we refer to, for example, [46] for a new meso-scale model related to a multi-phase composite material, [47] for a formulation based on variational tools and [48,49] for a particle system formulation. These articles can be used to start different approaches to the introduced problem.

2. Experimental procedure

2.1 BaTiO₃ powder

The traditional solid-state method used in our research for BaTiO₃ synthesis has strict criteria with regard to particle size and purity of the raw material, while it is also easy to produce non-homogeneous powders [50]. The advantages of the solid-state reaction used in mass production are its simplicity, precise stoichiometric control and low cost, with the main disadvantage being that the high calcining temperature results in very large and non-uniform grain sizes [51]. As a consequence, this method does not allow the production of materials with a high dielectric constant [51]. In order to eliminate the disadvantages of the solid-state reaction, a variety of novel wet chemical synthesis methods have been developed [52].

Three steps of calcination were applied in order to obtain BaTiO₃ powder. The full report of these procedures was comprehensively explained in [53]. Figure 1 presents the X-ray diffraction (XRD) pattern of BaTiO₃ after the third activation and calcination. We remark in Figure 1 that the third mechanical activation and calcination provide single-phase powder considering that the diffraction pattern of BaTiO₃ powder after the third calcination and activation has revealed the tetragonal phase and fully overlaps with the standard.

Figure 1.

X-ray diffraction pattern of BaTiO₃ powder after the third activation and calcination.

2.2 BaTiO3 pellets

Based on the lattice constant value (see [54]) in BaTiO₃ pellets, the theoretical density is considered to be 6.02 g/cm³. Table 1 presents some properties of BaTiO₃ pellets in which the relative density is very high and can be compared with [54].

Table 1.

The apparent density, water absorbability and apparent porosity of sintered pellets.

Water absorbability, %	Apparent density, g/cm³	Apparent porosity, %	Relative density, %
0.05	5.82	0.28	96.68

2.3 Grain size

One of the important parameters that can influence the mechanical and electrical properties of piezoelectric materials is the grain size of a polycrystalline aggregate. Not only is the intergranular fracture is influenced by grain size, but also the electrical response can experience a huge effect via changing the porosity and grain size of the microstructure. One of the popular methods to determine this feature is analysing images that are prepared by the SEM machine. Figure 2 presents the pellet surface at two different magnifications, ×350 (Figure 2(a)) and ×1000 (Figure 2(b)), in which compact densification can be observed. In the microstructure of polycrystalline BaTiO₃ both large grains, which have sizes from 50 to 100 µm, and small ones, which have dimensions from 1 to 3 µm, can be seen. As shown in Figure 2, grains boundaries are very well fitted together and only slender pores are visible. These results are in agreement with the data for pellet densification presented in Table 1 and the investigations carried out by Yoon and Lee [55].

Figure 2.

Scanning electron microscopy microstructure of the pellet surface at two different magnifications: (a) ×500; (b) ×1000.

3. Preparing the boundary element method model in computer-aided design software and the numerical programming algorithm

In order to analyse the microstructure of piezoelectric ceramics, it is essential to have some information on such as points on boundaries. Some researchers used the artificial microstructure created by a quasi-random point generator; however, this method is not able to discretise all kinds of artificial grains and the SEM image form laboratory. In this study, SEM images of specimens as well as artificial grains have been discretised by generating a comprehensive numerical algorithm.

Firstly, one image processing tools has been applied to allocated special pixels in every grain in the micro-structure. In the case of SEM images, the computer-aided design (CAD) software has been used in order to clarify the boundary of the SEM images. Afterwards, as shown in Figure 3, the output from the CAD software has been used in the MATLAB software and image processing tools in order to allocate special pixels for every grain in the polycrystalline aggregate. However, for the artificial specimen, we ignored the use of CAD software and the image was directly inputted to the image processing algorithm, as illustrated in Figure 4.

Figure 3.

The procedure of allocating special pixels for every grain in the polycrystalline aggregate obtained by the scanning electron microscopy machine. CAD: computer-aided design.

Figure 4.

The procedure of allocating special pixels for every grain in the artificial polycrystalline aggregate.

In order to apply BEMs, the discretisation of the grain boundary is essential. So, some numerical algorithms were created to obtain the coordinates of the grain boundary via pixel information. Figure 5 presents the results of this process. It can be concluded that the generated algorithms were successful for both artificial grains and the SEM image and they are able to discretise the grains very well.

Figure 5.

Node points on grain boundaries that are separated into grains in specimens of the grain boundary in two different views: (a) artificial microstructure; (b) microstructure obtained by the scanning electron microscopy image.

In this study, due to the low ratio of grain boundaries compared with grain surfaces, the grain boundaries will not be considered as a separate part, that is, the study of behaviour of grain boundaries will be ignored. In Figure 5, it is clear that there are gaps between grain boundaries and each grain has its own boundaries. So, the distances between grain boundaries must be eliminated and the grains must have the same boundaries in the area that they have contact. To achieve this goal, more numerical algorithms were designed in order to unify the boundary of grains that are in contact. The output of these algorithms is shown in Figure 6. As it is obvious in Figure 6, these numerical algorithms act very well: not only were the distances between the grains eliminated and they have unified boundaries, but also in the corners, the grains have the same coordinates, which will increase the accuracy of numerical models.

Figure 6.

Discretised grain boundaries in two different views: (a) artificial microstructure; (b) microstructure obtained by the scanning electron microscopy image.

4. Application of the boundary element method in micro-scale

We consider the model of brittle materials for representation of the microstructure and the basic assumption that the whole microstructure consists of polycrystalline brittle materials. Due to a small area of grain boundaries compared with the grain, the modelling of them is ignored in the current research and only the intergranular fracture is the subject of the current study. To investigate crack nucleation, propagation and branching in brittle polycrystalline, the BEM (as suggested in [34]) is being used here. Let us write the constitutive relation for piezoelectric ceramics in two dimensions, namely the 1–3 plane. As is known, there is an extra field in piezoelectric materials, namely electric displacement $D_{i}$ and electrical potential $ϕ$ , which need to be considered

\begin{matrix} [\begin{matrix} σ_{11} \\ σ_{22} \\ σ_{12} \end{matrix}] = [\begin{matrix} c_{11} & c_{13} & 0 \\ c_{13} & c_{33} & 0 \\ 0 & 0 & c_{44} \end{matrix}] [\begin{matrix} u_{1, 1} \\ u_{2, 2} \\ u_{1, 2} + u_{2, 1} \end{matrix}] + [\begin{matrix} \begin{matrix} 0 \\ 0 \\ e_{15} \end{matrix} & \begin{matrix} e_{31} \\ e_{33} \\ 0 \end{matrix} \end{matrix}] [\begin{matrix} φ,_{1} \\ φ,_{2} \end{matrix}] \\ [\begin{matrix} D_{1} \\ D_{2} \end{matrix}] = [\begin{matrix} \begin{matrix} 0 & 0 & e_{15} \end{matrix} \\ \begin{matrix} e_{31} & e_{33} & 0 \end{matrix} \end{matrix}] [\begin{matrix} u_{1, 1} \\ u_{2, 2} \\ u_{1, 2} + u_{2, 1} \end{matrix}] + [\begin{matrix} κ_{11} & 0 \\ 0 & κ_{33} \end{matrix}] [\begin{matrix} φ,_{1} \\ φ,_{2} \end{matrix}] \end{matrix}

(1)

The displacement vector for piezoelectric materials can be written as

u_{I} = {\begin{matrix} u_{i}, & I = 1, 2 \\ ϕ & I = 3 \end{matrix}

(2)

In a similar way, the traction vector is

σ_{iJ} = {\begin{matrix} σ_{ij}, & J = 1, 2 \\ D_{J} & J = 3 \end{matrix} \Rightarrow t_{I} = {\begin{matrix} t_{i} = \sum_{j = 1}^{2} σ_{ij} n_{j}, I = 1, 2 \\ D = \sum_{j = 1}^{2} D_{j} n_{j}, I = 3 \end{matrix}

(3)

In Figure 6(b), the microstructure subject of the current study is presented. In this microstructure, the basic assumptions are that (i) for two grains which are adjacent, namely grain numbers 9 and 10 (see Figure 6(b)) the value of traction is equal and (ii) these two grains are completely attached to each other as follows

\begin{matrix} t^{9} = - t^{10}, δ u^{I} = u^{9} - u^{10} = 0 \\ ϕ^{9} = ϕ^{10}, D^{9} = D^{10} \end{matrix}

(4)

where the values of $t^{I}$ and $δ u^{I}$ are the interface traction and the relative displacement jump in global coordinates. The values of traction and displacement can be represented in two components, the $x$ direction ( $t_{x}^{I}, u_{x}^{I}$ ) and the $y$ direction ( $t_{y}^{I}, u_{y}^{I}$ ).

The fundamental boundary element equation of displacement for each grain can be written as

\begin{matrix} C_{ij}^{k} (z') u_{j}^{k} (z') + \int_{S_{c}^{k}} T_{ij}^{k} (z', z) u_{j}^{k} (z) {dS}_{c}^{k} + \int_{S_{nc}^{k}} T_{ij}^{k} (z', z) u_{j}^{k} (z) {dS}_{nc}^{k} = \\ \int_{S_{c}^{k}} U_{ij}^{k} (z', z) t_{j}^{k} (z) {dS}_{c}^{k} + \int_{S_{nc}^{k}} U_{ij}^{k} (z', z) t_{j}^{k} (z) {dS}_{nc}^{k} k = 1, . . ., number of grains \end{matrix}

(5)

In Equation (5), $T_{ij}^{k} (z', z)$ and $U_{ij}^{k} (z', z)$ are fundamental solutions for traction and displacement, which can be obtained via the Stroh formulation for each grain. $S_{c}^{H}$ and $S_{nc}^{H}$ belong to the boundary of the internal grain (and contact boundary) and the free boundary. $C_{ij}^{k}$ refers to free terms. Equation (5) is discretised in unique global coordinates; however, in order to simplify the application of the boundary condition along with interface equations, it is better to convert Equation (5) into local coordinates.

It is worth mentioning that the transformation from global coordinate to local coordinate for the boundary condition or free boundary can be done by

{\tilde{t}}^{I} = [\begin{matrix} {\tilde{t}}_{s} \\ {\tilde{t}}_{n} \\ D \end{matrix}] = [\begin{matrix} - n_{y} & n_{x} & 0 \\ n_{x} & n_{y} & 0 \\ 0 & 0 & 1 \end{matrix}] [\begin{matrix} t_{x} \\ t_{y} \\ D \end{matrix}] = R t^{I}, {\tilde{u}}^{I} = [\begin{matrix} {\tilde{u}}_{s} \\ {\tilde{u}}_{n} \\ φ \end{matrix}] = [\begin{matrix} - n_{y} & n_{x} & 0 \\ n_{x} & n_{y} & 0 \\ 0 & 0 & 1 \end{matrix}] [\begin{matrix} u_{x} \\ u_{y} \\ φ \end{matrix}] = R u^{I}

(6)

where $R$ denotes the local rotation matrix from global coordinates to local coordinates. The same procedure can be used for the displacement and traction of internal grains. Now, Equation (5) can be rewritten as

\begin{matrix} {\tilde{C}}_{ij}^{k} (z') {\tilde{u}}_{j}^{k} (z') + \int_{S_{c}^{k}} {\tilde{T}}_{ij}^{k} (z', z) {\tilde{u}}_{j}^{k} (z) {dS}_{c}^{k} + \int_{S_{nc}^{k}} {\tilde{T}}_{ij}^{k} (z', z) {\tilde{u}}_{j}^{k} (z) {dS}_{nc}^{k} = \\ \int_{S_{c}^{k}} {\tilde{U}}_{ij}^{k} (z', z) {\tilde{t}}_{j}^{k} (z) {dS}_{c}^{k} + \int_{S_{nc}^{k}} {\tilde{U}}_{ij}^{k} (z', z) {\tilde{t}}_{j}^{k} (z) {dS}_{nc}^{k} k = 1, . . ., number of grains \end{matrix}

(7)

where matrices in local coordinates can be defined in the continuum form

{\tilde{C}}_{ij}^{k} = C_{in}^{k} R_{nj}^{k} {\tilde{U}}_{ij}^{k} = U_{in}^{k} R_{nj}^{k} {\tilde{T}}_{ij}^{k} = T_{in}^{k} R_{nj}^{k}

(8)

or in the matrix form

\begin{matrix} [\begin{matrix} {\tilde{C}}_{11} & {\tilde{C}}_{12} & {\tilde{C}}_{13} \\ {\tilde{C}}_{21} & {\tilde{C}}_{22} & {\tilde{C}}_{23} \\ {\tilde{C}}_{31} & {\tilde{C}}_{32} & {\tilde{C}}_{33} \end{matrix}] = [\begin{matrix} - n_{y} & n_{x} & 0 \\ n_{x} & n_{y} & 0 \\ 0 & 0 & 1 \end{matrix}] [\begin{matrix} C_{11} & C_{12} & C_{13} \\ C_{21} & C_{22} & C_{23} \\ C_{31} & C_{32} & C_{33} \end{matrix}], \\ [\begin{matrix} {\tilde{T}}_{11} & {\tilde{T}}_{12} & {\tilde{T}}_{13} \\ {\tilde{T}}_{21} & {\tilde{T}}_{22} & {\tilde{T}}_{23} \\ {\tilde{T}}_{31} & {\tilde{T}}_{32} & {\tilde{T}}_{33} \end{matrix}] = [\begin{matrix} - n_{y} & n_{x} & 0 \\ n_{x} & n_{y} & 0 \\ 0 & 0 & 1 \end{matrix}] [\begin{matrix} T_{11} & T_{12} & T_{13} \\ T_{21} & T_{22} & T_{23} \\ T_{31} & T_{32} & T_{33} \end{matrix}], \\ [\begin{matrix} {\tilde{U}}_{11} & {\tilde{U}}_{12} & {\tilde{U}}_{13} \\ {\tilde{U}}_{21} & {\tilde{U}}_{22} & {\tilde{U}}_{23} \\ {\tilde{U}}_{31} & {\tilde{U}}_{32} & {\tilde{U}}_{33} \end{matrix}] = [\begin{matrix} - n_{y} & n_{x} & 0 \\ n_{x} & n_{y} & 0 \\ 0 & 0 & 1 \end{matrix}] [\begin{matrix} U_{11} & U_{12} & U_{13} \\ U_{21} & U_{22} & U_{23} \\ U_{31} & {\tilde{U}}_{32} & U_{33} \end{matrix}] \end{matrix}

(9)

5. Cohesive law and crack initiation rules

5.1 Interface condition

In this section, the different conditions at the interfaces will be analysed. Two different statuses will be considered for node pairs as follows.

Undamaged or intact region $UR$ : in this region, the node pairs on the interface are completely bonded together and displacement compatibility and traction equilibrium equations are governed in this situation. Another essential assumption is that the traction in node pairs does not exceed threshold limit. The separation between node pairs is allowed in this region, that is, in Equation (4) we will have $δ {\tilde{u}}^{I} = 0$ .

Damaged or cohesive region $CR$ : if the traction exceeds the threshold limit in certain node pairs, then it can be claimed that node pairs experience a crack initiation mode and they are damaged. The displacement compatibility equation ( $δ {\tilde{u}}^{I} \neq 0$ ) is not valid anymore; however, the traction equilibrium holds, but it can be obtained via using some cohesive law equations and a displacement jump. After node pairs enter this region an irreversible cohesive law will determine whether they will enter to the next region, which is cracked one or not.

Cracked region $FR$ : in this region, free micro-cracks happen for node pairs and the micro-cracking analysis rules are governed here. In the current work, a frictional contact mechanism is modelled.

5.2 Interface equations

In the UR zone the displacement compatibility equations and traction equilibrium equation hold. It is worth mentioning that the traction equilibrium equation will be held in all of the status; however, displacement compatibility equations in the UR zone are different than the ones in CR, that is, for obtaining the displacement jump a cohesive traction-separation law will be introduced.

5.2.1 Undamaged interface continuity equations

When the node pairs are in the undamaged zone, the compatibility and equilibrium equations in local coordinates are the following

\begin{matrix} {\tilde{u}}_{i}^{a} + {\tilde{u}}_{i}^{b} = - δ {\tilde{u}}_{i}^{ab} = 0 \\ i = 1, 2 or i = t and n \\ {\tilde{t}}_{i}^{a} - {\tilde{t}}_{i}^{b} = 0 \\ {\tilde{ϕ}}^{a} = {\tilde{ϕ}}^{b}, {\tilde{D}}^{a} = {\tilde{D}}^{b} \end{matrix}

(10)

In this equation, $i = 1$ represents tangential components and $i = 2$ the normal components in local coordinates. Through components of normal and tangential traction, the effective traction can be defined as

t_{eff} = {[{〈 {\tilde{t}}_{n} 〉}^{2} + {(\frac{β}{α} {\tilde{t}}_{s})}^{2}]}^{\frac{1}{2}}

(11)

In this equation, 〈.〉 denotes Macauley brackets, which are defined by $〈 x 〉 = max (0, x)$ . Regarding Equation (11), it is clear that only the tensile value of normal traction contributes to the effective value of traction. $β$ and $α$ are suitable coefficients that determine the contribution of opening and sliding modes in the damage process. As long as $t_{eff} < T_{max}$ the interface has undamaged status, $T_{max}$ being the maximum value for traction.

5.2.2 Damaged interface continuity equations

Once effective traction exceeds the maximum value, $t_{eff} \geq T_{max}$ , the node pairs enter the damaged zone and crack initiation starts at those node pairs. In this status, the displacement jump happens for node pairs and one variable is needed to calculate the displacement jump in the cohesive zone and traction as well. In different works, some relations are applied in order to find variables. According to Ortiz and Pandolfi [21], the non-dimensional effective opening displacement can be introduced as

d = {[{〈 \frac{δ {\tilde{u}}_{n}}{δ {\tilde{u}}_{n}^{cr}} 〉}^{2} + β^{2} {(\frac{δ {\tilde{u}}_{s}}{δ {\tilde{u}}_{s}^{cr}})}^{2}]}^{\frac{1}{2}}

(12)

In this relation, $δ {\tilde{u}}_{n}^{cr}$ and $δ {\tilde{u}}_{s}^{cr}$ are the critical values of displacement in normal and sliding modes in which interface failure or micro-cracks can occur. $β$ is a weight that determines which mode has more influence on effective displacement in the mixed-mode equation. It is worth mentioning that the value of effective opening displacement varied between 0 and 1, $0 \leq d \leq 1$ , that is, at $d = 0$ , there is no displacement jump and the interface is completely intact and, at $d = 1$ , we have complete failure or cohesion in which micro-cracks can be defined.

Let ϕ be the scalar potential

ϕ (δ {\tilde{u}}_{n}, δ {\tilde{u}}_{s}) = ϕ

(13)

Through derivation of the last relation, a cohesive traction-separation law can be approached as follows

{\tilde{t}}_{n} = \frac{\partial ϕ}{\partial δ {\tilde{u}}_{n}}, {\tilde{t}}_{s} = \frac{\partial ϕ}{\partial δ {\tilde{u}}_{s}}

(14)

In the literature [28,29] this integral equation was suggested as scalar potential

ϕ (d) = δ {\tilde{u}}_{n}^{cr} \int_{0}^{d} t_{e} (ς) d ς

(15)

in which $d$ is the current effective opening displacement defined by Equation (12). The following equation can be considered for $t_{e} (ς)$

t_{e} (ς) = T_{max} \frac{1 - d^{*}}{d^{*}} ς with d^{*} = max {d} \in [0, 1]

(16)

The variable $d^{*}$ needs a little description. The value of $d$ in the current increment cannot be lower that its value in previous increments. In other words, another variable must be defined in order to hold the highest values of $d$ in all of the increments. This is why $d^{*}$ is defined by $0 \leq d \leq d^{*}$ . However, the question is how the amount of $d$ can be decreased during the process. Here, damage evaluation divides the process into two different ones, (i) the loading process in which $d = d^{*}$ and (ii) the unloading and reloading process in which $0 \leq d \leq d^{*}$ . More accurately, it can be defined as

d^{*} = d and \overset{\cdot}{d} > 0 loading; 0 \leq d \leq d^{*} and {\begin{matrix} unloading \overset{\cdot}{d} \leq 0 \\ reloading \overset{\cdot}{d} > 0 \end{matrix}}

(17)

Let consider Equations (15) and (16): the cohesive traction can be calculated via differentiation of the potential equation and by using the effective opening displacement. The normal and sliding amount of traction in local coordinates can be shown as

{\begin{matrix} {\tilde{t}}_{s} = α T_{max} \frac{1 - d^{*}}{d^{*}} \frac{δ {\tilde{u}}_{s}}{δ {\tilde{u}}_{s}^{cr}} \\ {\tilde{t}}_{n} = T_{max} \frac{1 - d^{*}}{d^{*}} \frac{δ {\tilde{u}}_{n}}{δ {\tilde{u}}_{n}^{cr}} \end{matrix}

(18)

The equation that relates α to $β$ is

α = β^{2} \frac{δ {\tilde{u}}_{n}^{cr}}{δ {\tilde{u}}_{s}^{cr}}

(19)

The variation of normal and tangential components of traction versus normal and tangential displacement openings is presented in Figure 7. This relation guarantees the accuracy of potential energy in Equation (15), that is, the potential energy or fracture energy are path independent and in normal and tangential modes the fracture energy is equal. However, because of microstructural details, such as interface roughness, the fracture energy cannot be considered independent [56], that is, the fracture energy in normal mode is different from in tangential mode. In order to prove this assumption, the displacement jump can be categorised into normal and tangential components via considering the angle $γ$ between the vector of the displacement jump in the cohesive region and the tangential to the cohesive surface, as follows

\begin{matrix} δ {\tilde{u}}_{n} = \tilde{u} \cos (γ) \\ δ {\tilde{u}}_{s} = \tilde{u} \sin (γ) \end{matrix}

(20)

Figure 7.

Variation of normal and tangential components of cohesive traction versus opening displacement: (a) three-dimensional plot where the traction change is due to both components of the opening displacement; (b) assigning $δ {\tilde{u}}_{s} = 0.1$ and $δ {\tilde{u}}_{n} = 0.1$ for plotting normal and tangential components separately.

Now, let us write an integration relation to show how the total strain energy release rate can be computed

G (γ) = \int_{0}^{u_{cr}} t_{n} (\tilde{u}) d δ {\tilde{u}}_{n} + t_{s} (\tilde{u}) d δ {\tilde{u}}_{s} = \frac{T_{max}}{2} \frac{\frac{\cos^{2} (γ)}{δ {\tilde{u}}_{n}^{cr}} + α \frac{\sin^{2} (γ)}{δ {\tilde{u}}_{s}^{cr}}}{{(\frac{\cos (γ)}{δ {\tilde{u}}_{n}^{cr}})}^{2} + β^{2} {(\frac{\sin (γ)}{δ {\tilde{u}}_{s}^{cr}})}^{2}}

(21)

The strain energy release rate for modes one and two can be shown by

{\begin{matrix} G_{n} = G (0) = \frac{1}{2} T_{max} δ {\tilde{u}}_{s}^{cr} \\ G_{s} = G (\frac{π}{2}) = \frac{α}{2 β^{2}} T_{max} δ {\tilde{u}}_{s}^{cr} \end{matrix}

(22)

It is worth mentioning that the values of $α$ and $β$ are independent from each other and can be chosen arbitrarily. However, according to the amount of $\frac{G_{s}}{G_{n}}$ , Equation (22) is not always valid except when the ratio of strain energies in both modes is equal to one, or according to microstructural properties, the fracture energy in sliding and normal modes is unequal. Hence, in this work, Equation (21) is used instead of potential energy (Equation (15)). The main advantage of Equation (21) over Equation (15) is that it is easy to distinguish between fracture energy in normal and tangential modes. This idea was developed by van den Bosch et al. [57], and Snozzi and Molinari [56] modified it by providing an example about a masonry wall.

5.2.3 Cracked region

In order to model free micro-cracks, a frictional mechanism is considered. In Tables 2 –4 the frictional model, which contains either contact or separation, is presented. Here $μ$ is a friction coefficient.

Table 2.

The constraints for different frictional mechanisms of the micro-crack surface.

Separation	Contact
	Stick	Slip
$\begin{matrix} t^{n} = 0 \\ t^{s} = 0 \end{matrix}$	$\begin{matrix} δ u^{n} = 0 \\ δ u^{s} = 0 \end{matrix}$	$\begin{matrix} δ u^{n} = 0 \\ t^{s} \pm μ t^{n} = 0 \end{matrix}$

Table 3.

Violation check for the contact mode.

Assumption/decision	Separation	Contact
Separation	$δ u^{n} > 0$	$δ u^{n} \leq 0$
Contact	$t^{n} \leq 0$	$t^{n} > 0$

Table 4.

Violation check for the contact status.

Assumption/decision	Stick	Slip
Stick	$t^{s} > μ t^{n}$	$t^{s} \leq μ t^{n}$
Slip	$t^{s} δ u^{t} < 0$	$t^{s} δ u^{t} \geq 0$

6. The evaluation algorithm

After the image processing procedure in aggregate, the BEM must be applied to all the grains separately. Then, the grain boundary formulation according to different types of regions can be applied in the assembly of the matrix. As was mentioned in the previous section, each grain in the polycrystalline aggregate is supposed to have unique mechanical and electrical properties depending on the different orientation that every grain has. Actually, the orientation of grains in the polycrystalline aggregate was determined via allocating random orientations. The status equations at interface between two adjacent grains namely, $a$ and $b$ , are

{\begin{matrix} Ω_{ab}^{i} [{\tilde{u}}_{j}^{a}, {\tilde{u}}_{j}^{b}, φ^{a}, φ^{b}, {\tilde{t}}_{j}^{a}, {\tilde{t}}_{j}^{b}, D^{a}, D^{b}] = 0 \\ Θ_{ab}^{i} [{\tilde{t}}_{j}^{a}, {\tilde{t}}_{j}^{b}, D^{a}, D^{b}] = 0 \end{matrix}

(23)

In this system of equations $Ω_{ab}^{i}$ depends on the interface, while $Θ_{ab}^{i}$ represents for the equilibrium equation in all regions. In the undamaged region, these two equations can be written as

{\begin{matrix} Ω_{ab}^{i} = {\tilde{u}}_{i}^{a} + {\tilde{u}}_{i}^{b} = - δ {\tilde{u}}_{i}^{ab} \\ Ω_{ab} = {\tilde{ϕ}}^{a} - {\tilde{ϕ}}^{b} = δ {\tilde{ϕ}}^{ab} \\ i = 1, 2 \\ Θ_{ab}^{i} = {\tilde{t}}_{i}^{a} - {\tilde{t}}_{i}^{b} \\ Θ_{ab} = {\tilde{D}}^{a} - {\tilde{D}}^{b} \end{matrix} or s and n

(24)

When node pairs enter the cohesive region, these equations can be rewritten by taking into account the cohesive law

\begin{matrix} {\begin{matrix} Ω_{ab}^{s} = C_{s} (d^{*}) δ {\tilde{u}}_{s}^{ab} - {\tilde{t}}_{s} \\ Ω_{ab}^{n} = C_{n} (d^{*}) δ {\tilde{u}}_{n}^{ab} - {\tilde{t}}_{n} \\ Ω_{ab} = {\tilde{ϕ}}^{a} - {\tilde{ϕ}}^{b} \\ Θ_{ab}^{i} = {\tilde{t}}_{i}^{a} - {\tilde{t}}_{i}^{b} i = 1, 2 or s and n \\ Θ_{ab} = {\tilde{D}}^{a} - {\tilde{D}}^{b} \end{matrix} \\ \begin{matrix} C_{s} (d^{*}) = \frac{1 - d^{*}}{d^{*}} \frac{α T_{max}}{δ {\tilde{u}}_{s}^{cr}} \\ C_{n} (d^{*}) = \frac{1 - d^{*}}{d^{*}} \frac{T_{max}}{δ {\tilde{u}}_{n}^{cr}} \end{matrix} \end{matrix}

(25)

We remark that all the variables and equations are defined in local coordinates. In this paper, the uniform displacement boundary condition applied to the microstructure is represented as

u_{i}^{m} = ε_{ij}^{M} x_{j}^{m}

(26)

A linear relation for one single grain, namely grain number 10 in Figure 6, can be written as follows

[\begin{matrix} H_{nc}^{10} & H_{c}^{10} \end{matrix}] [\begin{matrix} u_{nc}^{10} \\ u_{c}^{10} \end{matrix}] = [\begin{matrix} G_{nc}^{10} & G_{c}^{10} \end{matrix}] [\begin{matrix} t_{nc}^{10} \\ t_{c}^{10} \end{matrix}]

(27)

Via applying the interface and boundary condition and by using Equation (27) for all the grains in the polycrystalline aggregate, the adaptive system of equations that contains BEM and grain interface equations can be written in a matrix form

[\begin{matrix} \begin{matrix} A^{1} B^{1} & 0 & \dots & 0 \\ 0 & A^{2} B^{2} & \dots & 0 \\ ⋮ & ⋮ & ⋱ & ⋮ \\ 0 & \dots & 0 & A^{n_{g}} B^{n_{g}} \\ Interface equation & Ω \\ Interface equation & Θ \end{matrix} \end{matrix}] [\begin{matrix} x_{nc}^{1} \\ x_{c}^{1} \\ ⋮ \\ x_{nc}^{n_{g}} \\ x_{c}^{n_{g}} \end{matrix}] = [[\begin{matrix} \begin{matrix} D^{1} & 0 & \dots & 0 \\ 0 & D^{2} & \dots & 0 \\ ⋮ & ⋮ & ⋱ & ⋮ \\ 0 & \dots & 0 & D^{n_{g}} \\ 0 \\ 0 \end{matrix} \end{matrix}]] [\begin{matrix} y_{nc}^{1} \\ y_{nc}^{2} \\ ⋮ \\ 0 \\ y_{nc}^{n_{g}} \end{matrix}]

(28)

To get this linear equation, the uniform displacement boundary condition has been applied. In Equation (28), $n_{g}$ represents the number of grains. $D$ contains the matrix $H_{nc}$ in Equation (27) and $A$ the matrix $G_{nc}$ in Equation (27); however, contrary to Equation (27), all these matrices here are defined in local coordinates. $B$ is a matrix that has the form $[\begin{matrix} G_{c} & - H_{c} \end{matrix}]$ . The vector $x_{nc}$ contains the unknown value of $\tilde{t}$ in the non-contact or boundary condition area and $x_{c}$ contains ${[\begin{matrix} \tilde{t} & \tilde{u} \end{matrix}]}^{T}$ , which are traction and displacement in the contact area, while the vector $y_{nc}$ contains the prescribed boundary conditions in the non-contact area. $Ω$ and $Θ$ represent the interface status in Equation (23) and contain many zeros, some 1 and −1 in the undamaged areas; however, they can be described through Equation (23) in the cohesive region too. It is clear that the total matrix of Equation (28) is linear and it contains the nonlinear equations of cohesive law that can be solved in every iteration.

Hence, the linear system of Equation (28) can be extended in another form

[\begin{matrix} W (BCs) \\ Ω \\ Θ \end{matrix}] X = [\begin{matrix} y (BCs, load increments) \\ 0 \\ 0 \end{matrix}] or AX = Y

(29)

In this matrix, $W$ on the left-hand side contains the information about boundary conditions, shape functions and Jacobians. $y$ contains the information about boundary conditions, shape functions and Jacobians, and load increments. $Ω$ and $Θ$ are interface conditions.

In order to obtain macro-stress and macro-strain we use averaging theory

{\bar{σ}}_{ij}^{m} = \frac{1}{V^{t}} \sum_{k = 1}^{n_{g}} \int_{S_{nc}^{k}} x_{i}^{m} t_{j}^{m} d S^{t} {\bar{ε}}_{ij}^{m} = \frac{1}{2 V^{t}} \sum_{k = 1}^{n_{g}} \int_{S_{nc}^{k}} (u_{i}^{k} n_{j}^{k} + u_{j}^{k} n_{i}^{k}) {dS}_{nc}^{k} + \int_{S_{c}^{k}} (u_{i}^{k} n_{j}^{k} + u_{j}^{k} n_{i}^{k}) {dS}_{c}^{k}

(30)

where $S_{nc}^{k}$ is a non-contact boundary condition, $S_{c}^{k}$ a contact boundary condition and $V^{t}$ is the volume of the microstructure.

7. Numerical algorithm and damage evaluations

Let us to explain the numerical algorithms that are used in order to apply the boundary element formulations to the interface equations, damage initiations and propagation procedures. These algorithms have been developed in order to analyse the piezoelectric ceramics from both a micro-mechanical and multi-scale point of view.

7.1 Micro-structural and damage evaluation algorithm

Thanks to the system of linear equation in (28) or (29), it is possible to apply a load to the system. Actually, the external load can be increased incrementally. In the cases where all of the interfaces are in the undamaged status, the micro-structure shows the same behaviour for all the external loads and there are no notable problems in applying the external load; however, this situation is not always available in on the polycrystalline aggregates and, after some particular conditions, the damage will start to be evaluated. More precisely, we define the load increment $Δ λ$ as $λ_{k + 1} = λ_{k} + Δ λ$ . As long as $t_{eff}$ does not exceed the value of $T_{max}$ in all the grains, the choice of $Δ λ$ does not need much attention, but, after the first node pairs enter into the cohesive region, $Δ λ$ must be chosen carefully in order to have precise results.

We developed a search algorithm in order to find the $i$ th node pair that enters the cohesive zone, via considering indicators like $t_{eff} > T_{max}$ . This event happens at the $k$ th increment, in which the load is $λ_{k}$ , and an arbitrary amount is allocated for $d_{start}^{*}$ in the $i$ th node pairs, while matrix $A$ in Equation (29) is updated to solve the linear equation at the $k + 1$ th increment. In other words, $A_{k + 1}$ is prepared by using the cohesive law’s equations (Equations (18), (23) and (24)).

In the $k + 1$ th increment, $λ_{k + 1}$ can be obtained via choosing $Δ λ$ and then new load $λ_{k + 1}$ is applied to the right-hand side of Equation (29) by considering $Y_{k + 1}$ . In the next step, solution $X_{k + 1}$ is approached for the $k + 1$ th increment. For the node pairs, which are in the cohesive region, the amount of $d_{k + 1}^{1}$ is evaluated by using Equation (12) for the $k + 1$ th increment and the first iteration. The number of iterations can be $1, 2, 3, . . ., j$ . Then the amount of $d_{k + 1}^{1}$ in the first iteration will be compared with $d_{k}^{*}$ . Status two can happen if $d_{k + 1}^{1}$ is lower than $d_{k}^{*}$ , that is, the unloading–reloading situation on the node pair where the amount of $d_{k}^{*}$ will be the same for the next step $d_{k + 1}^{*} = d_{k}^{*}$ . On the other hand, if $d_{k + 1}^{1}$ is higher than $d_{k}^{*}$ , we have the loading status. In this situation, the amount of $d_{k + 1}$ will be evaluated for the next increments $(d_{k + 1}^{1}, d_{k + 1}^{2}, d_{k + 1}^{3}, . . ., d_{k + 1}^{j})$ and the evaluation will be stopped when the following convergence condition will be joined

| d_{k + 1}^{j + 1} - d_{k + 1}^{j} | \leq ε_{d} | d_{k + 1}^{j} |

(31)

$ε_{d}$ being a tolerance parameter. Now, the amount of $d_{k + 1}^{*}$ can be updated by $d_{k + 1}^{*} = d_{k + 1}^{j + 1}$ and matrix $A$ will be updated via new amount of $d_{k + 1}^{*}$ for the next incremental load.

7.2 Numerical simulation of micro-damage and the cohesive law without considering contact and frictional effects

The efficiency of cohesive law and the developed grain boundary formulation will be proved via some numerical examples. We consider a practical example based on the electrical and mechanical properties of barium titanate (see Table 5).

Table 5.

The property of Barium Titanate piezoelectric ceramic.

$c_{11} (\frac{N}{m^{2}})$	$c_{13} (\frac{N}{m^{2}})$	$c_{33} (\frac{N}{m^{2}})$	$c_{44} (\frac{N}{m^{2}})$	$e_{31} (\frac{C}{m^{2}})$	$e_{33} (\frac{C}{m^{2}})$	$e_{15} (\frac{C}{m^{2}})$	$κ_{11} (\frac{C}{Vm})$	$κ_{33} (\frac{C}{Vm})$
$9.54 \times 10^{10}$	$3.26 \times 10^{10}$	$7.57 \times 10^{10}$	$2.56 \times 10^{10}$	$- 6.98$	$13.8$	$13.4$	$60 \times 10^{- 10}$	$54.7 \times 10^{- 10}$

The specified value for parameters was introduced during previous sections and subsections, and can be found in Table 6. The fracture toughness (see [58]) is considered here $G_{s} = 2.34 \times 10^{- 3} \frac{N}{mm}$ for both normal and tangential modes. In other words, we assume $\frac{G_{s}}{G_{n}} = 1$ except when another value is specified.

Table 6.

The interface property of piezoelectric ceramic.

$T_{max} (MPa)$	$α$	$β$	$d_{start}^{*}$	$ε_{d}$
80	1	$\sqrt{2}$	$0.5 \times 10^{- 5}$	$0.5 \times 10^{- 3}$

The uniform displacement boundary condition is applied on the microstructure in the horizontal direction. Figure 8 shows the relation between stress and strain in the microstructure of homogenous piezoelectric ceramic of Figure 9, which is subjected to tensile mechanical load in the horizontal direction (i) by considering the piezoelectric effect and (ii) without the piezoelectric effect. It can be concluded that the piezoelectric coupling has a significant influence on the fracture criteria of aggregates, that is, through considering the piezoelectricity and for a specified amount of strain, the polycrystalline aggregate without the piezoelectric effect can experience more stress. It could be deduced that for a specific strain energy that is entered to the system, the piezoelectric aggregate should be divided between mechanical and electric effects, while in polycrystalline without piezoelectricity this energy can only have influence on mechanical properties. The intergranular fracture criteria of this aggregate, which are subjected to tensile mechanical load in the horizontal direction by considering the piezoelectric effect and without it, are shown in Figure 9. It can be concluded in this figure that by increasing the strain, the crack grows faster and faster.

Figure 8.

Normal strain and normal stress of polycrystalline aggregate shown in Figure 9 and subjected to tensile mechanical load in the horizontal direction by considering piezoelectric effects and without piezoelectric effects.

Figure 9.

The fracture pattern for polycrystalline aggregate under horizontal tensile mechanical load: (a) without piezoelectric effects; (b) with piezoelectric effects.

In this work, the node pairs enter the cohesive zone when the effective traction exceeds the ultimate value of traction. Hence, it could be interesting to know how the traction changes in the node pairs. Figure 10 represents the variation of normal and tangential traction in three different zones of barium titanate polycrystalline subjected to horizontal tensile load. It could be concluded that all of the zones in this figure enter the cohesive zone. As can be seen, these figures are similar to Figure 7, where the normal displacement jump and normal traction have positive values as expected. Moreover, for negative values of tangential traction, the value of the tangential displacement jump is negative and vice versa.

Figure 10.

Variation of normal and tangential traction versus normal and tangential displacement jump in three different zones of barium titanate polycrystalline subjected to horizontal tensile load: (a) area in which the components are shown; (b) normal component; (c) tangential component.

7.3 The influence of frictional contact on the mechanical and electrical response of piezoelectric ceramic

The presented results were based on separation and no contact was assumed. In this section, slip and stick conditions are assumed in micro-mechanical grain boundary formulations as well as separation. Let us to compare the results according to the assumption of frictional contact and contact without friction. Figure 11 represents the tension test of piezoelectric polycrystalline aggregate in Figure 9 with and without considering frictional effects. It can be concluded that in the curve with frictional effect the stress has a greater value for the same amount of strain between $0.2 \times 10^{- 3}$ and $1 \times 10^{- 3}$ , and that can be justified by considering the frictional effect and slip or stick mode: the interface traction has a non-zero value as a result of the higher amount of average stress.

Figure 11.

Normal strain and normal stress of piezoelectric polycrystalline aggregate of Figure 9, which is subjected to horizontal tensile load by considering frictional effects and without frictional effects.

Figure 12 demonstrates the modelling in the stress–strain micro-structure of piezoelectric polycrystalline aggregate of Figure 13, which is subjected to a horizontal compression mechanical load, and Figure 13 demonstrates the micro-fracture path. It is obvious in Figure 13 that for a higher value of friction coefficient the micro-crack propagates slower, that is, a higher value of frictional loads opposes the sliding and carrying more loads. In Figure 12, as can be observed, while the polycrystalline aggregate experiences failure and crack propagations, the amount of stress still increases because of inter-lock at the interface and the internal fraction results in intergranular traction higher than zero.

Figure 12.

The stress–strain curve of the micro-structure in Figure 13 subjected to compression horizontal mechanical load for different values of friction coefficient.

Figure 13.

The micro-fracture path of piezoelectric polycrystalline subjected to compression horizontal load for different values of friction coefficient at $ε_{11}^{m} = 0.3$ .

8. Conclusions

The crack initiation and evaluation in piezoelectric ceramic according to image processing-based discretisation algorithms have been presented. The micro-structure was obtained through doing some experimental procedures. Then, some comprehensive numerical algorithms have been developed to utilise the geometrical information of pixels, which are allocated to every grain in the micro-aggregate in order to prepare the final model for numerical analysis. The nonlinear grain boundary formulations in different zones, namely the undamaged, cohesive and failure zones, have been combined with the BEM in an adaptive system of linear equations. The image processing-based discretisation algorithms and developed grain boundary formulation were found suitable for investigation of the intergranular micro-fracture in piezoelectric ceramic. Furthermore, the influences of different parameters, such as piezoelectricity effects, friction effects and mechanical loading, have been analysed by using the developed methods. The results of numerical analyses indicate that the image processing-based numerical algorithm is a suitable method for discretising the natural and artificial microstructure of materials. Moreover, considering friction effects in analysis can increase the strength of the materials.

The problem here presented can be conveniently used to study damage and fracture at the micro-level in some mechanical systems, for example pantographic lattices [59]. Numerical tools able to model fracture in this case can be found in [60,61].

Footnotes

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This work was supported by the People Programme (Marie Curie Actions) of the European Union’s Seventh Framework Programme FP7/2007-2013/, Grant Number PITN-GA-2013-606878.

References

Aliabadi

. The boundary element method, vol II: application to solids and structures. Chichester: Wiley, 2001.

Atluri

Pipkins

. Large deformation analysis of plates and shells. In: Beskos

(ed.) Boundary element analysis of plate and shells. 1991, 141–166. Berlin Heidelberg: Springer-Verlag.

Banerjee

Butterfield

. Boundary element methods in engineering science. London: McGraw-Hill, 1981.

Ben Mariem

Hamdi

. A new boundary element method for fluid-structure interaction problems. Int J Numer Meth Eng 1987; 24: 1251–1267.

Bernhard

Gardner

Mollo

, et al. Prediction of sound fields in cavities using boundary-element methods. AIAA J 1987; 25: 1176–1183.

Fenner

. Boundary element methods for engineers: part II plane elastic problems. London: Bookboon, 2014.

Beskos

Maier

(eds). Boundary element advances in solid mechanics. Wien: Springer, 2014.

Pasternak

. Boundary integral equations and the boundary element method for fracture mechanics analysis in 2D anisotropic thermoelasticity. Eng Anal Bound Elem 2012; 36: 1931–1941.

Yin

. Deconstructing plane anisotropic elasticity: Part I: the latent structure of Lekhnitskii’s formalism. Int J Solids Struct 2000; 37: 5257–5276.

10.

Yang

Qin

. Micro-mechanical analysis of composite materials by BEM. Eng Anal Bound Elem 2004; 28: 919–926.

11.

Qin

. Micromechanics-BE solution for properties of piezoelectric materials with defects. Eng Anal Bound Elem 2004; 28: 809–814.

12.

Placidi

. A variational approach for a nonlinear 1-dimensional second gradient continuum damage model. Continuum Mech Thermodynam 2015; 27: 623-638.

13.

Placidi

. A variational approach for a nonlinear one-dimensional damage-elasto-plastic second-gradient continuum model. Continuum Mech Thermodynam 2016; 28: 119–137.

14.

Trindade

Benjeddou

. Finite element characterisation of multilayer d31 piezoelectric macro-fibre composites. Compos Struct 2016; 151: 47–57.

15.

Rahali

Assidi

Goda

, et al. Computation of the effective mechanical properties including nonclassical moduli of 2.5 D and 3D interlocks by micromechanical approaches. Composites Part B 2016; 98: 194–212.

16.

Jafari

Saeed

Nima

. Modeling of the plastic deformation of polycrystalline materials in micro and nano level using finite element method. J Nano Res 2013; 22: 41-60.

17.

Needleman

. Numerical simulations of dynamic interfacial crack growth allowing for crack growth away from the bond line. Int J Fract 1996; 74: 253–275.

18.

Tvergaard

. Effect of fibre debonding in a whisker-reinforced metal. Mater Sci Eng A 1990; 125: 203–213.

19.

Espinosa

Zavattieri

Gordon

. Adaptive FEM computation of geometric and material nonlinearities with application to brittle failure. Mech Mater 1998; 29: 275–305.

20.

Camacho

Ortiz

. Computational modeling of impact damage in brittle materials. Int J Solids Struct 1996; 33: 2899–2938.

21.

Ortiz

Pandolfi

. Finite-deformation irreversible cohesive elements for three-dimensional crack-propagation analysis. Int J Numer Meth Eng 1999; 44: 1267–1282.

22.

Espinosa

Zavattieri

Dwivedi

. A finite deformation continuum\discrete model for the description of fragmentation and damage in brittle materials. J Mech Phys Solids 1998; 46: 1909–1942.

23.

Giorgio

Scerrato

. Multi-scale concrete model with rate dependent internal friction. Eur J Environ Civ Eng 2017; 21: 821-839.

24.

Scerrato

Giorgio

Della Corte

, et al. Towards the design of an enriched concrete with enhanced dissipation performances. Cement Concrete Res 2016; 84: 48–61.

25.

Scerrato

Giorgio

Della Corte

, et al. A micro-structural model for dissipation phenomena in the concrete. Int J Numer Anal Meth Geomech 2015; 39: 2037–2052.

26.

Zavattieri

Raghuram

Espinosa

. A computational model of ceramic microstructures subjected to multi-axial dynamic loading. J Mech Phys Solids 2001; 49: 27–68.

27.

Verhoosel

Gutiérrez

. Modeling inter-and transgranular fracture in piezoelectric polycrystals. Eng Fract Mech 2009; 76: 742–760.

28.

Galvis

Sollero

. Boundary element analysis of crack problems in polycrystalline materials. Procedia Mater Sci 2014; 3: 1928–1933.

29.

Benedetti

Aliabadi

. Multiscale modeling of polycrystalline materials: a boundary element approach to material degradation and fracture. Comput Meth Appl Mech Eng 2015; 289: 429–453.

30.

Goldberg

Hopkins

. Application of the boundary element method to the micromechanical analysis of composite materials. Comput Struct 1995; 56: 721–731.

31.

Fotiu

Heuer

Ziegler

. BEM analysis of grain boundary sliding in polycrystals. Eng Anal Bound Elem 1995; 15: 349–358.

32.

Benedetti

Aliabadi

. A three-dimensional grain boundary formulation for microstructural modeling of polycrystalline materials. Comput Mater Sci 2013; 67: 249–260.

33.

Sfantos

Aliabadi

. A boundary cohesive grain element formulation for modeling intergranular microfracture in polycrystalline brittle materials. Int J Numer Meth Eng 2007; 69: 1590–1626.

34.

Benedetti

Aliabadi

. A three-dimensional cohesive-frictional grain-boundary micromechanical model for intergranular degradation and failure in polycrystalline materials. Comput Meth Appl Mech Eng 2013; 265: 36–62.

35.

Gulizzi

Milazzo

Benedetti

. An enhanced grain-boundary framework for computational homogenization and micro-cracking simulations of polycrystalline materials. Comput Mech 2015; 56: 631-651.

36.

Crocker

Flewitt

PEJ

Smith

. Computational modelling of fracture in polycrystalline materials. Int Mater Rev 2005; 50: 99–125.

37.

Okabe

Boots

Sugihara

, et al. Spatial tessellations: concepts and applications of Voronoi diagrams. Chichester: John Wiley & Sons, 2009, 501.

38.

Nemat-Nasser

Hori

. Micromechanics overall properties of heterogeneous materials. Amsterdam: Elsevier Science, 1999.

39.

Kouznetsova

Brekelmans

WAM

Baaijens

FPT

. An approach to micro-macro modeling of heterogeneous materials. Comput Mech 2001; 27: 37–48.

40.

Ren

Chen

, et al. Micro-cracks informed damage models for brittle solids. Int J Solids Struct 2011; 48: 1560–1571.

41.

Kouznetsova

Geers

MGD

Brekelmans

WAM

. Multi-scale constitutive modelling of heterogeneous materials with a gradient enhanced computational homogenization scheme. Int J Numer Meth Eng 2002; 54: 1235–1260.

42.

Carcaterra

dell’Isola

Esposito

, et al. Macroscopic description of microscopically strongly inhomogenous systems: a mathematical basis for the synthesis of higher gradients metamaterials. Arch Ration Mech Anal 2015; 218: 1239–1262.

43.

Pideri

Seppecher

. A second gradient material resulting from the homogenization of an heterogeneous linear elastic medium, Continuum Mech Thermodynam 1997; 9: 241–257.

44.

Altenbach

Eremeyev

. On the linear theory of micropolar plates. ZAMM 2009; 89: 242–256.

45.

Misra

Poorsolhjouy

. Granular micromechanics based micromorphic model predicts frequency band gaps. Continuum Mech Thermodynam 2016; 28: 215–234.

46.

Contrafatto

Cuomo

Gazzo

. A concrete homogenisation technique at meso-scale level accounting for damaging behaviour of cement paste and aggregates. Comput Struct 2016; 173: 1–18.

47.

Cuomo

Contrafatto

Greco

. A variational model based on isogeometric interpolation for the analysis of cracked bodies. Int J Eng Sci 2014; 80: 173–188.

48.

Della Corte

Battista

dell’Isola

, et al. Modeling deformable bodies using discrete systems with centroid-based propagating interaction: fracture and crack evolution. In: dell’Isola

Sofonea

Steigmann

(eds) Mathematical modelling in solid mechanics. Advanced structured materials, vol 69. Singapore: Springer, 2017.

49.

Battista

Rosa

dell’Erba

, et al. Numerical investigation of a particle system compared with first and second gradient continua: deformation and fracture phenomena. Math Mech Solids 2017; 22(11): 2120–2134.

50.

Wang

Cui

. Preparation and characterization of BaTiO₃ powders and ceramics by the sol–gel process using organic monoacid as surfactant. Scr Mater 2007; 57: 623–626.

51.

Kao

Yang

. Preparation of barium strontium titanate powder from citrate precursor. Appl Organomet Chem 1999; 13: 383–397.

52.

Chen

Shen

Liu

, et al. Preparation and properties of barium titanate nanopowder by conventional and high-gravity reactive precipitation methods. Scr Mater 2003; 49: 509–514.

53.

Trzepieciński

Ryzińska

Biglar

, et al. Modelling of multilayer actuator layers by homogenisation technique using Digimat software. Ceram Int 2017; 43: 3259–3266.

54.

Duran

Gutierrez

Tartaj

, et al. Densification behaviour, microstructure development and dielectric properties of pure BaTiO3 prepared by thermal decomposition of (Ba, Ti)-citrate polyester resins. Ceram Int 2002; 28: 283–292.

55.

Yoon

Lee

. Processing of barium titanate tapes with different binders for MLCC applications - part I: optimization using design of experiments. J Eur Ceram Soc 2004; 24: 739–752.

56.

Snozzi

Molinari

. A cohesive element model for mixed mode loading with frictional contact capability. Int J Numer Meth Eng 2013; 93: 510–526.

57.

van den Bosch

Schreurs

PJG

Geers

MGD

. A critical evaluation of the exponential Xu and Needleman cohesive zone law for mixed mode decohesion. Eng Fract Mech 2005; 72: 2247–2267.

58.

Park

Sun

. Fracture criteria for piezoelectric ceramics. J Am Ceram Soc 1995; 78: 1475–1480.

59.

Spagnuolo

Barcz

Pfaff

, et al. Qualitative pivot damage analysis in aluminum printed pantographic sheets: numerics and experiments. Mech Res Comm 2017; 83: 47–52.

60.

Turco

dell’Isola

Rizzi

, et al. Fiber rupture in sheared planar pantographic sheets: Numerical and experimental evidence. Mech Res Commun 2016; 76: 86-90.

61.

Turco

Rizzi

. Pantographic structures presenting statistically distributed defects: Numerical investigations of the effects on deformation fields. Mech Res Commun Volume 2016; 77: 65-69.

Application of the grain boundary formulation and image processing-based algorithm in micro-mechanical analysis of piezoelectric ceramic

Abstract

Keywords

1. Introduction

2. Experimental procedure

2.1 BaTiO3 powder

2.2 BaTiO3 pellets

2.3 Grain size

3. Preparing the boundary element method model in computer-aided design software and the numerical programming algorithm

4. Application of the boundary element method in micro-scale

5. Cohesive law and crack initiation rules

5.1 Interface condition

5.2 Interface equations

5.2.1 Undamaged interface continuity equations

5.2.2 Damaged interface continuity equations

5.2.3 Cracked region

6. The evaluation algorithm

7. Numerical algorithm and damage evaluations

7.1 Micro-structural and damage evaluation algorithm

7.2 Numerical simulation of micro-damage and the cohesive law without considering contact and frictional effects

7.3 The influence of frictional contact on the mechanical and electrical response of piezoelectric ceramic

8. Conclusions

Footnotes

Funding

References

2.1 BaTiO₃ powder