- To simulate the behavior of thin dielectric elastomer structures a finite solid shell element formulation is presented. Dielectric elastomers belong to the group of electroactive polymers and their use as actuators is caused by the efficient coupling between electrical energy input and mechanical energy output. Also the large elongation strain of 120-380% of the dielectric elastomer actuators and their light weight are advantages that make the material very attractive. Regarding the electro-mechanical coupling a constitutive model is expounded. For the definition of an electric stress tensor and a total stress tensor the electrical body force and couple are considered in the balance of linear momentum and angular momentum, respectively. The governing constitutive equations are derived and incorporated in a solid shell element formulation based on a Hu-Washizu mixed variational principle considering six fields: displacements, electric potential, strains, electric field, mechanical stresses, and dielectric displacements. This formulation allows large deformations and accounts for physical nonlinearities to capture the main characteristics of dielectric elastomers.
In recent years dielectric elastomers (DE) have become popular for the usage in actuators. The efficient transformation of electrical energy in mechanical energy and the ability to maintain large strains makes them very attractive. A Constitutive model, which describes the specific material behavior, is introduced by Dorfmann and Ogden [1], Steigmann [2] and the references therein. The numerical treatment in the context of the finite element method is presented by Vu et al. [3], who provide a brick element formulation, which incorporates the nonlinear constitutive model. In the present work a shell finite element formulation for DE is presented. It is motivated by the fact that the most DE devices are thin structures, which have a very high length to thickness ratio. The electromechanical coupling is considered in the body force and the couple density, see e.g. Eringen and Maugin [4] or Müller et al. [5]. The angular momentum equation is fulfilled by assuming a Maxwell stress tensor. The linear momentum equation is approximately fulfilled by the finite element method. A mixed solid shell element formulation is introduced. It incorporates specific interpolations for the displacements, strains, stresses as well as for the electric potential, electric field and dielectric displacements. As nodal degrees of freedom three displacements and the electric potential are assumed. The mixed formulation allows for a consistent finite element approximation to avoid electromechanical locking effects. The element formulation is able to simulate large deformations. Some numerical examples show the applicability of the proposed solid shell element.
Let Φ be a deformation that point maps ⃗⃗ of the reference configuration to ⃗ of the current configuration at time , see Figure 1. The deformation gradient is declared as the tangent to Φ and given by
With this definition for the deformation gradient the right
⃗⃗ . (3)
To describe the shell formulation convective coordinates are introduced, where and are the in-plane coordinates and is the coordinate in thickness direction. The covariant contravariant basis vectors are defined in a standard manner by ⃗⃗ ⃗⃗ and ⃗⃗ ⃗⃗ , where is the Kroneckerdelta. With this convective description the deformation gradient also reads ⃗⃗ ⃗⃗ .
A potential is used to describe the electric field ⃗ . With satisfying the Laplace’s equation ⃗ reads ⃗ (⃗ ) ⃗ (⃗ ) order terms results in by neglecting higher ⃗ (∑ ) ⃗ (∑ ) ⃗ , ⃗ - (∑ ⃗ ) . ⃗ by
Applying the pull-back operation the physical electric field ⃗ is observed. The displacement vector ⃗ is defined
, respectively. are given on
To get a macroscopic representation of force and couple density a close look on the microscopic level is presented. Therefore the physical model will start with a particle of the current configuration to deduce electric force in the presence of electromagnetic matter and in absence of a magnetic field. ⃗ ⃗ ⃗ ⃗
⃗ . (4) (5) and ∫ ⃗ ∫ , ⃗ ⃗ ⃗
⃗⃗ . ⃗ ⃗⃗ .
After defining the resultant mass ∑ , charge ∑ , and polarization ⃗ ∑ ⃗ of the particle and performing a simple space averaging ∑ , ∑ , ⃗ ∑ ⃗ leads to the macroscopic body force density ⃗ ⃗ ⃗ , ⃗ ⃗ ⃗ , where ⃗ combines the electric contribution to the body force density. Deriving the corresponding couple density with the same arguments the resultant couple of the particle with respect to the mass center reads ∑ ⃗ ⃗ ∑ ⃗ ⃗ .
Here it is assumed that the mass dipole moment ∑ ⃗ ⃗ . With the resultant polarization ⃗ ∑ ⃗ and a space averaging it follows the macroscopic couple density
Relative to the mass center eccentric by ⃗ not to be confused with the convective coordinates are the point charges and the point masses , as shown in Figure 2. In an electric field a force is acting on each point charge called Lorentz force and it is defined by ⃗ . The gravitation field ⃗ is acting on each point mass producing the Newton force given by ⃗ . Summing over leads to the resultant force on the particle: ∑ ⃗
The position of a point within the particle can be described by ⃗ ⃗ ⃗ . Assuming the gravitation to be constant it follows: ⃗ (⃗ ) ⃗ (⃗ ). This assumption and expanding the external field in a Taylor series (7) (8) (9) (10) (11) (12) (13)
The balance laws incorporating the electric force, couple, and power densities are summarized. After that the global integral forms and the local field equations are presented. The Cauchy stress tensor and an electric stress tensor are introduced. Conservation of mass ∫ , with is assumed. Localization of the material description results in . In quasi static processes the integral form of the balance of linear momentum is given as
where ⃗ is the body force density (8) and the traction vector on . With the Cauchy stress tensor the traction is determined by a linear map of the normal vector ⃗ , according to Cauchy’s stress theorem. Considering conservation of mass, applying the divergence theorem and the localization theorem, the field equation along with the jump condition are observed as ⃗ ⃗ ⃗ ⃗ on in
With as a prescribed traction on the boundary . The boundaries with a given traction and a given displacement satisfy ⋃ and ⋂ . The global form of the balance of angular momentum reads ∫ ⃗ ⃗ ∫ ⃗ where ⃗ is the couple density Using the integral theorem, considering conservation of mass and the linear momentum balance along with the localization theorem results in the field equation ⃗
For dielectric materials the conservation of charge in integral form with the surface charge density on reads
The dielectric displacement vector is denoted by ⃗ and determined by Gauss’ law ⃗ ⃗ . Applying the divergence theorem results in the field equation along with the jump condition ⃗ .
,
Introducing the energy function (⃗ ⃗ ) [(⃗ ⃗ ) where is a function of to fulfill material objectivity, here an Ogden-type material is chosen, and the susceptibility of the material is denoted by , the total stress and the dielectric displacements are derived as ⃗ ⃗ ⃗ ⃗ ( ⃗ ⃗ ⃗ , ( ⃗⃗ ⃗⃗
⃗ ⃗⃗ ⃗⃗
Here, is a prescribed surface charge. For the total boundaries with a given electric potential and a given surface charge it holds ⋃ and ⋂ . With the constitutive equations in matter the dielectric displacements are given as ⃗ ⃗ ⃗ , where is the permittivity in vacuum and the polarization depends on the considered material. An electric stress tensor is introduced as such that ⃗ ⃗ ⃗ and
⃗ . It is remarked that is also known as Maxwell stress tensor.
With it follows that the total stress tensor has to be symmetric. The remaining field equations and boundary conditions are ⃗
⃗ ⃗ ⃗ ) ⃗ ⃗ ⃗ ⃗ ⃗ ⃗ , (19) Considering the 2nd Piola-Kirchhoff stress tensor and , the pull-back of the dielectric displacement ⃗⃗ ⃗ and the transformation of the densities by leads to the material description )⃗ ⃗ on
in on on in in on on in in
Where ⃗ is the traction with respect to the reference configuration.
,
In this section a mixed variational formulation is introduced. Let the set { ⃗ , ( )- ⃗ ⃗ } be the space of admissible displacement variations and { , ( )} be the space of admissible electric potential variations. Further let * , ( )- + , * , ( )- + the spaces of admissible variations of the total stresses and strains and { ⃗⃗ , ( )- } , { ⃗ , ( )- } the spaces of admissible variations of the dielectric displacement and the electric field. Since the variations are arbitrary the field equations (24)-(27), the constitutive equations (29), (20), and the kinetic field equations (3), (4) are rewritten as Applying integration by part, using the divergence theorem and considering the boundary conditions the weak formulation reads )⃗ )
⃗ ) )) ⃗ ⃗⃗ )) ( , , ⃗ (
∫ ⃗ ⃗⃗ ⃗ ) ( ⃗ ⃗ ) ( ⃗⃗
) ) ⃗ ⃗⃗ with ⃗ ,
In this section a solid shell element is introduced. The finite element approximation is constructed in the sense that the whole domain is divided in element domains with ⋃ , where is the total number of elements. The geometry, displacements, and electric potential are approximated as ⃗⃗ ∑ ⃗⃗ , ⃗ ∑ ⃗ , and ∑ with the same interpolation functions ( )( )( ) , at nodes . The vectors ⃗⃗ and ⃗ contain the nodal coordinates and displacements, respectively. Arranging in the matrix , - with , - , the virtual quantities are interpolated as [ ⃗ ]
⃗ , where ⃗ [⃗ ⃗ ⃗ ⃗ ] is the vector of nodal degrees of freedom ⃗ ,⃗ ⃗ ⃗ - .
The gradient fields are defined with respect to the curvilinear coordinates . The constitutive equations and will be given with respect to a local orthonormal coordinate system ⃗ . This necessitates a transformation of the strains and the electric field. Introducing the transformation matrix
The approximation on element level of the virtual gradient field ⃗ reads defined with ⃗ ⃗ ‖⃗ ⃗ ‖ where
⃗ and ⃗ ‖⃗ ‖ , ⃗ , ⃗ ⃗ ⃗ . The Jacobian matrix is denoted by ⃗⃗ ⃗ . For the sake of a compact notation the contravariant components of the virtual strain tensor and the virtual electric field vector are arranged in a generalized vector , and ⃗ Quantities, which are evaluated at the element center are denoted with the index and is a identity matrix. The transformation matrix is obtained by considering and . The matrices and ⃗ are defined as - and [ is defined with some ANS interpolations in [6] ⃗ / ⃗ / ( ( ⃗ ) ) . ) .
⃗ ⃗ ⃗ - . where denotes the partial derivative of the shape function with respect to the curvilinear coordinates. The matrix at the node is determined as The physical stresses and dielectric displacements are derived from the potential function and are arranged in the vector 0 ⃗ ⃗ ⃗ ⃗ 1 . Here, ⃗ are independently assumed quantities for strain and electric field components and are approximated with the following interpolations, see Klinkel and Wagner [7]: ⃗
⃗⃗ , with ⃗⃗
⃗ . The matrices ] , and [
⃗ ⃗ are given as ] , ⃗⃗
, with The matrix and reads
( [ ( ) , ⃗⃗ ( ) , ⃗ /
] . ] The approximation of the independent field is defined as with The weak formulation of will be approximated on element level as following: , with
. In Eq. (59) the body and surface loads are determined by ⃗̃ , ⃗ - and ̃ 0⃗
Having in mind that (50) is solved iteratively the following approximation on element level is obtained , ⃗ ⃗ ⃗ ⃗
⃗ ) ) ⃗ ⃗ , ⃗ ⃗ ⃗ ⃗ ⃗ ⃗ ⃗ ⃗
⃗⃗ , ∫ ∫ ∫ ∫ ∫ +∫ ∫ ( ∫ ( ∫ ∫ ∫ ( ∫ ⃗ ⃗ , ⃗ ⃗⃗ ⃗⃗ ⃗⃗ ∫ ∫ and ⃗ ⃗ ⃗ ⃗
⃗ ⃗
⃗ , ⃗ ⃗ ⃗ ⃗ ] [ ⃗ ]) Taking into account that the finite element interpolations for the fields ⃗ and ⃗ are discontinuous across the element boundaries a condensation on element level yields the element stiffness matrix and the right hand side vector with the unknown incremental nodal displacements and the electric potential. The update of the internal degrees of freedom reads , ⃗ ⃗ , with ⃗ , ⃗ - and ⃗ contains the components of , , and ⃗⃗ accordingly to the vector notation. The weak form is solved iteratively by employing Newton-Raphson’s method. This requires the linearization with
( ⋃ obtains
( ) and ⃗ ) ⃗ . After assembly over all elements , ⃗ ⋃ ⃗ and ⃗ ⋃ ⃗ one Considering the above interpolations and one obtains the following matrices: where
is the matrix of [7] and vectors ⃗ ⃗ ) . (⃗ ⃗
⃗ / ⃗ ) , and .
The relative permittivity is given as ( ) .
In Table 1 it is shown that there is only eigenvalue number 18 much greater than zero. This eigenvalue accounts for the volumetric locking mode which arises by incompressible materials like DE. At this point no other locking modes are observed. The result is also the same when taking an irregular non-cube element. B. Bending Actuator
The second example is presented to demonstrate the valid reflection of the electro-mechanical coupling phenomenon of DEs. Therefore a square plate with the dimension of is investigated. It is clamped on one side and consists of two layers with the thickness of each, see Figure 3. The dataset of the material is the same as given in the first example.
To get a bending answer of the thin structure either the upper or lower layer is loaded by an electric field applied in thickness direction of the structure. Due to the electromechanical coupling the loaded layer responds with an inplane extension. Because of the eccentricity the bending effect of the whole structure occurs. Figure 4 shows several deformed configurations.
Fig. 4 Bending actuator - several deformed configurations In the simulation the surface charge is increased from to and decreases back to again. Then the other layer is loaded in the same way. The tip deflection versus the applied voltage is shown in Figure 5.
80 60 -60 -80 1000 2000 3000 4000 5000 6000
7000
Load j [V]
The presented element formulation is based on a mixed variational approach. It results in an independent interpolation of the displacement, electric potential, strains, electric field, stresses and dielectric displacements. The element possesses only four nodal degrees of freedom, displacements and the electric potential. It allows a consistent approximation of the electromechanical coupled problem. The governing field equations and boundary conditions are presented. The formulation accounts for geometric and material nonlinear behavior. Further tasks will be to embed the material formulation into an improved stabilized element formulation. For more accurate simulations other material models for the mechanical part should be discussed.
We want to thank the Landesforschungszentrum CM² for financially supporting our project Mikro- & Strukturmechanik zur Analyse des nichtlinearen Deformationsverhaltens von dielektrischen und porösen Elastomeren (MSADEL).