The Meshless Analysis of Scale-Dependent Problems for Coupled Fields

The meshless local Petrov–Galerkin (MLPG) method was developed to analyze 2D problems for flexoelectricity and higher-grade thermoelectricity. Both problems were multiphysical and scale-dependent. The size effect was considered by the strain and electric field gradients in the flexoelectricity, and higher-grade heat flux in the thermoelectricity. The variational principle was applied to derive the governing equations within the higher-grade theory of considered continuous media. The order of derivatives in the governing equations was higher than in their counterparts in classical theory. In the numerical treatment, the coupled governing partial differential equations (PDE) were satisfied in a local weak-form on small fictitious subdomains with a simple test function. Physical fields were approximated by the moving least-squares (MLS) scheme. Applying the spatial approximations in local integral equations and to boundary conditions, a system of algebraic equations was obtained for the nodal unknowns.


Introduction
In nanocomposites the nano-sized particles are incorporated into a matrix. In these materials the ratio of surface to volume is significantly larger than in their bulk-sized equivalent. Then, their properties can be many times improved with respect to those known for the component parts.
The mechanical strength, toughness and electrical or thermal conductivity can be drastically improved in nanocomposites. Therefore, they are starting to be intensively utilized in many engineering fields. However, in nano-sized structures it is observed the size-effect phenomena if the characteristic length of material structure is compared with the size of the analyzed body [1][2][3][4][5]. Even some new phenomena are observed in nano-sized structures, like electric polarization in centro-symmetric crystals. It is explained by the direct flexoelectricity effect [6][7][8]. With respect to piezoelectricity, it can be viewed as a higher order effect [9]. If stresses are proportional to the gradients of electric intensity vector, we are talking about the converse flexoelectricity [10][11][12]. The strain-and electric intensity vector-gradients are very large in nano-sized dielectrics and they should be considered in constitutive equations.
Nanotechnology is also successfully utilized in improving of thermoelectric properties [13].
Thermoelectric materials have a potential to convert waste heat directly into electricity if the thermal conductivity is reduced without affecting the high electrical conductivity [14]. The thermal conductivity is reduced significantly in nanostructures only. Scattering of phonons is observed only on interfaces of nanostructures. Due to scattering of phonons the thermal conductivity is reduced. Since the electrons are smaller, they are not scattered and the electric conductivity is not reduced. It requires to develop a theory for heat conduction, where size effect is considered.
The microstructural material characteristics are not considered in the classical continuum mechanics and results are size-independent. Atomistic models are able to consider the size-effect in structural elements. However, extremely high requirements are put on computer memory in this approach. It is more convenient to develop advanced continuum models, where the intrinsic length parameter (characteristic of material microstructure) is considered [15][16][17]. Former gradient theory is very complicated due to many length scales as material parameters, which are not available. Therefore, Aifantis [18] simplified former theory by introducing only one length parameter. The nonlocality should be considered in the heat conduction problems if the temperature gradients are large [19]. For a special weight function in nonlocal integral expression of the heat flux vector it is possible to transform the integral form into a differential equation with a characteristic length parameter representing the nonlocality. In both flexoelectric and thermoelectric problems, we have constitutive equations with the intrinsic material parameter representing microstructure and higher derivatives of physical fields than in corresponding problems described by classical theory. In both problems there are interactions of several physical fields.
It is needed to have a reliable and accurate computational tool to solve these multiphysical problems described by gradient theories with intrinsic material parameters. Higher order derivatives in governing equations require the C 1 -continuous elements in numerical domain discretization methods to guarantee the continuity of variables and their derivatives on interfaces of elements. It is a difficult task to obey such a requirement. It is more convenient to develop the mixed formulation in the finite element method (FEM) [20,21]. The order of continuity of the Moving Least-square (MLS) approximation in the Meshless Local Petrov-Galerkin method (MLPG) can be tuned to a desired value very easily [22][23][24].
In the present paper, the authors have developed a meshless method based on the MLPG weakform to solve multiphysical problems in dielectric and thermoelectric solids. Both the direct and converse flexoelectricity is considered in dielectric solids. Nodal points are spread on the analyzed domain and each node is surrounded by a small circle for simplicity, but without loss of generality. The spatial variations of primary physical fields are approximated by the moving least-squares (MLS) scheme. After performing the spatial integrations, a system of algebraic or ordinary differential equations for unknown nodal values is obtained. The essential boundary conditions on the global boundary are satisfied by the collocation. Numerical examples are presented and discussed to compare the results obtained by the gradient theory with those obtained by classical theory.

The direct and converse flexoelectricity
The electric enthalpy density for piezoelectric solids with strain and electric intensity vector gradients and can be written as [25,26] where a and c denote the second-order permittivity and the fourth-order elastic constant tensors, respectively. Symbols e and f are used for the piezoelectric and flexoelectric coefficients, respectively. Symbol g denotes the higher order elastic coefficients. Finally, the symbols b and h are used for the converse flexoelectric coefficients and higher-order electric parameters, respectively.
The strain tensor ij  and the electric field vector j E are related to the displacements i u and the electric potential  by The strain-gradients are denoted by symbol η ( ) The constitutive equations are obtained from the electric enthalpy density expression (1) as where ij  , k D , jkl  and ij Q are the stress tensor, electric displacements, higher order stress and electric quadrupole, respectively.
In the simplified gradient theory [18,27,28] only one internal length-scale parameter l is present.
Then, the higher-order elastic parameters with li  being the Kronecker delta. A similar idea of simplification has been applied to the higher-order electric parameters where q is another length-scale parameter.
In the simplified direct flexoelectricity there are considered two independent parameters 1 f and 2 f for the direct flexoelectric coefficient 20]. Taking into account the above simplifications, the electric enthalpy density can be written as Finally, we reduce the number of independent converse flexoelectric coefficients ijkl b . The stresses induced by electric intensity vector in the orthotropic piezoelectric material with poling direction along the transversal isotropy 3 x -axis can be written as 11  A similar form is considered for induced stresses induced by the converse flexoelectricity ( ) with three independent converse flexoelectric coefficients Then, the final form of the electric enthalpy is given by The constitutive equations (4) for orthotropic materials can be rewritten into a matrix form as 11 11 11, 1 2 33 33, 3 13 13, Using the variational principle of least action, it is possible to derive the governing equations for the considered constitutive equations [31] ,, ( Essential and natural boundary conditions (b.c.) follow from the variational formulation of boundary value problems: where : and the traction vector, and the electric charge are defined as with : The jump at a corner (x c ) on the oriented boundary contour Γ is defined as

The MLPG formulation
The presence of gradients of strains and electric intensity vector in the electric enthalpy requires C 1 continuous interpolations of primary fields, i.e. displacements and electric potential. Recently, the mixed FEM has been developed for considered electro-elasticity problem [32]. The Meshless Local Petrov-Galerkin method (MLPG) with the Moving Least-square (MLS) approximation is convenient approach for problems with higher order derivatives, since the order of continuity can be tuned to a desired value [22][23][24].
The MLPG method is based on the local weak-form with local fictitious subdomains q  constructed for each node q x inside the analysed domain. The geometry of this subdomain can be arbitrary. However, it is appropriate to select a circular shape for simple numerical evaluation of integrals. One can write the local weak-form of the first governing equation (15) as where * () im u x is a test function.
Applying the Gauss divergence theorem to domain integrals in (23) one can write ** , , , where q  is the boundary of the local subdomain which consists of three parts q q q q tu L  =     , in general (see Fig. 1).

If a Heaviside step function is chosen for the test function
the local weak-form (24) is transformed into the local integral equation where Similarly, we get local integral equation for the second governing equation (15) ( ) For numerical solution of the above integral equations (26) and (27), the MLS approximation of trial functions is applied. The primary fields (mechanical displacements and electric potential) are given by [23] 11( where aa d =− xx , and a r is the size of the support domain.
The strains and electric intensity vector is approximated by where the matrices , , , T kk C Λ Φ F are defined in eq. (11) and (13) x . The essential boundary conditions are satisfied in the strong-form at nodal points . If the approximation formulas (28) and (30) Substituting the MLS-approximation (32) and (35) into the local boundary-domain integral equations (26) and (27), we obtain the system of algebraic equations for unknown nodal quantities  Fig. 2 is analyzed by the FEM [31] and the MLPG. The piezoelectric material PZT-5H is chosen for the study. To investigate influence of the strain gradient and electric intensity gradient parameters various integer numbers  and  are selected in numerical analyses.

Gradient theory in thermoelectric materials
The thermoelectric conversion efficiency is high if the thermal conductivity is low. It can be reduced significantly in nano-sized structures. It is due to comparable sizes of phonon mean free path and the structure. Phonons are scattered on interfaces and thermal conductivity is reduced. For this purpose it is needed to develop a theory of heat conduction, where size effect is considered. It is well known that there is no size effect considered in the classical local theory of Fourier heat conduction. Similarly to the elasticity problems in nanostructures it is possible to consider the size effect here through the nonlocal heat transport [19]. The heat flux vector in nonlocal theory is given by where the temperature differences are denoted by where w is the volume density of heat source. By this way it is possible to replace the integro-differential form of the constitutive law in (40) by a more convenient differential form given in (43). Then, higher order derivatives in the governing equation appear in this non-local theory of heat conduction than in the classical local The constitutive equations for thermoelectric materials with higher order heat conduction theory can be written as       Then, the governing equations for stationary thermoelectric problem are given by conservation of energy and electric charge as The weak form of these equations can be written as where P , Q and  are independent boundary densities conjugated with / p  =  n ,  and  , respectively, and given as with  being the heat flux, i n and i  are the Cartesian component of the unit tangent vector on  , and the jump at a corner on the oriented boundary contour  is defined as The rate of work of the external "forces" ( ) ,, PQ  and body source is given by If only the Joule heating plays the role of heat sources, and the governing equations become ,, Furthermore, from the weak formulation, one can deduce the following boundary conditions for coupled thermoelectric problem considered within higher-grade theory essential b.c.: Substituting the constitutive relationships into the governing equations, we obtain the non-linear system of the PDEs for primary field variables  and  ( ) Recall that owing to the Joule heat, the problem is non-linear even if the temperature dependence of material coefficients were neglected. Finally, making use the proportionality relationship ij ij s  = , the system of governing equations become

The MLPG formulation in thermoelectricity
One can see in the previous chapter that MLPG method is based on the local weak-form with local fictitious subdomains q  . The local weak-form of the first governing equation (57) is given as where * () u x is a test function.
Applying the Gauss divergence theorem to two domain integrals in (58), one can write ** , , , where q  is the boundary of the local subdomain q  .
The test function can be arbitrary and we have selected a Heaviside step function stand for the final and starting points on q   with prescribed heat flux.
The local integral equation for the second governing equation (57) is given as The MLS approximation of trial functions is applied for numerical solution of the above local integral equations (61) and (62). The temperature and electric potential are approximated by [23] 11( The essential boundary conditions are satisfied in the strong-form at nodal points

Numerical examples
An axially symmetric thermoelectric problems, as shown in Fig. 6, is analysed in the example.
The thermoelectric material Bi2Te3, is considered. It has the following material constants (Yang In classical thermoelectricity, it is possible to find the analytical solution.

Conclusions
The meshless Petrov-Galerkin (MLPG) method has been successfully applied to multiphysical problems described by advanced continuum models with microstructural effects. Strain-and electric intensity vectorgradients are considered in constitutive equations for electric displacement and stresses in flexoelectricity, respectively. Similarly, the constitutive equations for thermoelectric materials contain higher order derivatives of temperature in the higher-grade heat conduction theory. It allows to describe the heat transfer more realistic in nanostructures. The governing equations are derived for both multiphysical problems, where size effects are considered. These equations contain higher order derivatives of physical fields than in the classical continuum models. Application of classical domain discretization methods to corresponding boundary value problems brings some difficulties with continuity of approximated fields.
The proposed MLPG computational method with the MLS approximation of fields is very convenient to solve governing equations of gradient theory with high-order derivatives. The order of continuity of the MLS approximation can be tuned to a desired value very easily. Therefore, the present computational method is promising to be applied to multiphysical problems described by gradient theories.