Optimal Design of Electrode Topology of Dielectric Elastomer Actuators Based on the Parameterized Level Set Method

Abstract

Dielectric elastomer actuators (DEAs) are widely used to drive soft machines, and optimal design is required in more advanced designs of soft robots. In this research article, a computational approach is originally proposed and validated for the topology optimization of electrodes of DEAs. The nonlinear finite elements of absolute nodal coordinate formulation are applied to model the dynamic electromechanical behaviors of DEAs. The parameterized level set method is employed to optimize the electrode topology of DEAs. Then, the method of system-wise equivalent static load is employed to convert the optimization problem of dynamic responses into the static one. Based on the sensitivity analysis, the normal velocity field for optimizing the electrode topology is derived. Finally, two case studies are presented to validate the proposed optimization approach. The experimental test is also performed, and the measured results are compared with the numerical ones to further validate the proposed methodology.

Introduction

Dielectric elastomers (DEs) are widely used in a variety of scientific and technical applications, such as soft robots, electronics, actuators, and sensors.^1–6 In recent years, DEs have been increasingly attracting researchers' attention due to their advantages such as the flexibility and adaptability over rigid materials.^7,8 A soft mechanical system actuated by dielectric elastomer actuators (DEAs) usually undergoes both large overall motions and large deformations. A detailed review of the existing studies on the dynamics of soft mechanical systems actuated by DEAs has been given in our previous work.⁹

The absolute nodal coordinate formulation (ANCF) originally proposed by Shabana^10–12 allows for capturing arbitrarily large and coupled displacements, and it is suitable for modeling the flexible bodies undergoing both large overall motions and large deformations. Based on the integration of computer-aided design, robotic algorithms, and ANCF elements, Shabana¹³ proposed a new methodology for capturing the complex geometries and large deformations of soft robotic systems. In our previous work,⁹ a new viscoelastic solid element of ANCF was proposed for meshing the components of DEs, and it was successfully applied to dynamic simulations of soft mechanical systems actuated by DEs. In this study, a new plate element of ANCF is first developed to model the dynamic electromechanical behaviors of the DEAs.

The electrode, especially its topology, is vital to the dynamic actuation performance of DEAs.^14,15 Previous research has mainly focused on the topology optimization of the electrode to improve the static response of the DEAs. Using pairs of Bezier curves, Wang et al.^16,17 optimized the electrode topology to improve the static response of the rod linear and rotary DEAs. Chen et al.¹⁸ used a level-set based optimization method to optimize the electrode topology of planar and three-dimensional (3D) DEAs. However, the dynamic topology optimization of the electrode for soft multibody systems actuated by DEAs has been rarely reported. Therefore, the main objective of this article is to establish an effective topology optimization method of the electrode for soft multibody systems actuated by DEAs.

A level set model can be adopted to implicitly capture the topology of the electrode covered on DEAs. Benefiting from the implicit representation, the level set method can easily and naturally handle complex shape and topology changes such as boundary splitting and merging.¹⁹ The parameterized level set method is an effective approach for shape and topology optimization by introducing the radial basis functions (RBFs) into the level set method.^20,21 Therein, a level set function is decoupled by a linear combination of a set of RBFs and coefficients,²² and this parameterized modeling with RBFs can maintain a relatively smooth level set function without reinitialization operation during the optimization process. In this research, the parameterized level set function is used to describe the topology of the electrode of DEAs.

The method of system-wise equivalent static load (ESL) is employed to perform the dynamic response optimization of a flexible multibody system. The ESL method was initially proposed by Kang et al.²³ for the structural optimization under dynamic loads, and then it was extended to the dynamic response optimization of components in flexible multibody systems.^24,25 The ESL sets are static loads that generate the same displacement fields as those from a dynamic load at a certain time step. Converting dynamic loads into ESL sets at selected time steps can approximately reproduce the evolution of the dynamic response for a continuous-time problem.

And the deviation from the original dynamic problem depends on the size of the selected time step. However, it should be noted that the effect of inertia is neglected or only partially involved, which can lead to inaccuracy in the dynamic optimization process. Via the ESL method, the dynamic response optimization of a flexible multibody system can be converted into a static one under multiple loading conditions. The latter is then efficiently solved by using standard techniques of structural optimization.

Sun et al.²⁶ studied the dynamic response optimization of linear elastic components in flexible multibody systems based on the ESL method. Tromme et al.²⁷ extended the ESL method to the system level, and the system-wise ESL is the static load defined at the system level. Consequently, the shortcomings of formulating an optimization problem based only on local responses can be remedied.

In this research, we propose an effective computational method for optimizing the electrode topology of DE components in a soft multibody system with large displacements. The parameterized level set method is adopted to formulate the dynamic topology optimization problem. The optimization problem of dynamic responses is converted into a static one under multiple loading conditions based on the system-wise ESL method.²⁷ The normal velocity field is derived to determine the evolution of the interface of the electrode topology. With case studies, the proposed computational method is demonstrated to be robust for the optimal design of DE electrodes in a dynamic and systematic manner.

Dynamic Analysis of Soft Systems Actuated by DEs

DEs plate elements of ANCF

As shown in Figure 1a, let X be the position vector of an arbitrary material point of the DE body in the reference configuration, and r is defined as the corresponding position vector in the current configuration. In the current configuration, the elastomer is subject to an electric potential $ϕ$ , and the induced electric field is defined as $E$ . The induced electric field is mapped back to the reference configuration and defined as the nominal electric field $Ẽ$ .

FIG. 1.

(a) Configuration description of the motion of a DE body; (b) schematic of a DEA with an electrode pattern $Ω_{0}$ . DE, dielectric elastomer; DEA, dielectric elastomer actuator.

In the framework of ANCF,^28,29 the generalized coordinate vector of one element is composed of nodal coordinates, namely, three position coordinates and nine position gradient coordinates of each node. Figure 1a illustrates a spatial plate element of ANCF with four nodes. The plate element in the reference configuration has dimensions of length l, width b, and thickness h. The position vector of an arbitrary point on the plate element can be written as^28,29 $r = S e .$ (1)

In Eq. (1), $S = [\begin{matrix} s_{1} I & s_{2} I & \dots & s_{16} I \end{matrix}]$ is the shape function matrix with $I$ denoting the 3 × 3 identity matrix, and the detailed definitions of $s_{i} (i = 1, 2, \dots, 16)$ can be found in the work by Mikkola and Shabana.²⁸ $e$ is the generalized coordinate vector of the element formulated by $e = {[\begin{matrix} e_{1}^{T} & e_{2}^{T} & e_{3}^{T} & e_{4}^{T} \end{matrix}]}^{T},$ (2)

The constitutive model of DEs

As shown in Figure 1a, the matrix of position gradients F of an arbitrary point of the element can be written as²⁸

where $X$ and $X^{l} = {[\begin{matrix} x & y & z \end{matrix}]}^{T}$ are, respectively, the global and local position vectors in the initial reference configuration, and r is the global position vector at time t in the current configuration.

As shown in Figure 1b, the reference configuration of the DE body is denoted as D, and it is selected as the design domain for the electrode pattern. $Ω_{0}$ denotes the domain distributed with the electrode, and $D ∖ Ω_{0}$ denotes the domain of the passive materials. Assuming that the nominal electric field is only along the thickness direction and the electric field is constant along the thickness direction,³⁰ the nominal electric field can be given as follows: $Ẽ = {[\begin{matrix} 0 & 0 & Ẽ_{3} \end{matrix}]}^{T},$ (4)

where $Ẽ_{3}$ is stated as follows:

According to thermodynamics, the principle of increase of entropy can be fulfilled by minimizing the free energy of DEs.^6,9 The free energy density is defined as follows: $ψ = ψ^{e q} + ψ^{n e q} + ψ^{e},$ (6)

where $ψ^{e q}$ denotes the hyperelastic strain energy density of the elastomer, $ψ^{n e q}$ denotes the non-equilibrium strain energy density accounting for viscoelasticity, and $ψ^{e}$ denotes the electrical free energy density. To simplify computation, we choose the compressible neo-Hookean material modeloi to describe the nonlinear elastic behaviors of elastomers.

Generally, the neo-Hookean model is more appropriate for cases of small strains (<40%), and more advanced models such as the Mooney-Rivlin model, the Yeoh model, or the Gent model should be considered for cases of larger strains or with strain-stiffening effects.^31,32 Nevertheless, the neo-Hookean model can be well fitted into the data of tensile tests of elastomers,³³ and the simulation results have shown good consistency with the experimental ones for the investigated range in our previous studies.^9,33,34 Specifically, the strain energy density $ψ^{e q}$ of the neo-Hookean model can be expressed as $ψ^{e q} = \frac{μ}{2} (I_{1} - 3) - μ ln J + \frac{K}{2} {(ln J)}^{2},$ (7)

where $μ$ denotes the shear modulus, K denotes the bulk modulus, $I_{1} = t r (C)$ is the first invariant of the right Cauchy-Green strain $C = F^{T} F$ , and $J = det (F)$ is the determinant of $F$ . The non-equilibrium strain energy density $ψ^{n e q}$ can be calculated by $ψ^{n e q} = \frac{μ^{n e q}}{2} (C_{11} ξ_{1}^{- 2} + C_{22} ξ_{2}^{- 2} + C_{11}^{- 1} C_{22}^{- 1} ξ_{1}^{2} ξ_{2}^{2} - 3),$ (8)

where $μ^{n e q}$ denotes the non-equilibrium shear modulus, $C_{11}$ and $C_{22}$ are the components of $C$ representing the in-plane stretches, and $ξ_{1}$ and $ξ_{2}$ represent the viscous components of these stretches. Here, $ξ_{1}$ and $ξ_{2}$ can be obtained by solving a set of nonlinear ordinary differential equations in iteration, the details of which can be found in our previous work.⁹ And the electrical free energy density $ψ^{e}$ can be calculated by $ψ^{e} = - \frac{ε}{2} J D_{33} Ẽ_{3}^{2},$ (9)

where $ε$ is the permittivity of the material, and $D_{33}$ is the component of the inverse of the right Cauchy-Green tensor $D = C^{- 1}$ .

By differentiating the free energy function formulated by Eq. (6) with respect to the position gradients, the first Piola-Kirchhoff stress tensor can be obtained as $P = \frac{\partial ψ}{\partial F} = \frac{\partial (ψ^{e q} + ψ^{n e q})}{\partial F} + \frac{\partial ψ^{e}}{\partial F} = P^{m} + P^{e},$ (10)

where $P^{m}$ and $P^{e}$ denote, respectively, the mechanical and electric parts of P and their detailed expressions are given in Section S1 in Supplementary Data S1.

The equations of motion

By using the principle of virtual work,³⁵ the weak form of the dynamic equations of a soft mechanical system actuated by DEAs can be expressed as $\{\begin{matrix} \int_{D} P^{m} : δ F d D + \int_{Ω_{0}} P^{e} : δ F d Ω_{0} = \int_{D} {(\tilde{B} - ρ \ddot{r})}^{T} δ u d D \\ + \int_{S_{t}} {\tilde{T}}^{T} δ u d S, \\ u = \bar{u}, o n S_{u}, \forall δ u \in U_{a d} . \end{matrix}$ (11)

In Eq. (11), $δ F$ is the variation of the position gradients $F$ , $\tilde{B}$ denotes the body force, $ρ$ denotes the mass density, and $\tilde{T}$ denotes the surface force. $δ u$ denotes the virtual displacement, and $U_{a d}$ is the kinematically admissible displacement field. $\bar{u}$ is the prescribed displacement on Dirichlet boundary S_u.

In the dynamic analysis of soft multibody systems actuated by DEs, Eq. (11) is a set of differential algebraic equations, and the detailed derivation of the dynamic equation has been shown in our previous work.⁹

Formulation of Topology Optimization

RBFs-based parameterized level set method

In this subsection, a brief introduction to the RBF-based parameterized level set method is presented.

As shown in Figure 2, based on the level set method, the two-dimensional (2D) electrode topology is represented implicitly through a higher-order dimensional scalar function, 3D level set function, which is usually defined as follows:

FIG. 2.

An example of the level set model: (a) 3D level set function; (b) 2D electrode topology. 2D, two-dimensional; 3D, three-dimensional.

where $x = {[\begin{matrix} x & y \end{matrix}]}^{T}$ defined in the 2D coordinate system is the position vector of an arbitrary point in $\bar{D}$ , and $\bar{D}$ is the projection of design domain D. ${\bar{Ω}}_{0}$ denotes the projection of the domain distributed with the electrode $Ω_{0}$ , $\partial {\bar{Ω}}_{0}$ is the interface of ${\bar{Ω}}_{0}$ , and $\bar{D} ∖ {\bar{Ω}}_{0}$ denotes the projection of the domain of the passive materials $D ∖ Ω_{0}$ .

In the level set-based method, moving the interface of the shape is equivalent to transporting the scalar implicit function by solving the evolution equation.¹⁹ Assuming that the shape ${\bar{Ω}}_{0}$ evolves in a pseudo time $τ \in ℛ^{+}$ , where $ℛ^{+}$ denotes the positive real number field, then $Φ (x (τ), τ) = 0, i f x (τ) \in \partial {\bar{Ω}}_{0} (τ) .$ (13)

By differentiating Eq. (13) with respect to $τ$ , the following evolution equation, also known as the Hamilton-Jacobi equation,³⁶ is obtained as $\frac{\partial Φ}{\partial τ} + {(\nabla Φ)}^{T} \frac{d x}{d τ} = \frac{\partial Φ}{\partial τ} + {(\nabla Φ)}^{T} V = 0,$ (14)

where $\nabla Φ$ represents the gradient vector of the scalar function $Φ$ , and $V = \frac{d x}{d τ}$ is the velocity of the point x with respect to $τ$ . Eq. (14) can be rewritten as follows: $\frac{\partial Φ}{\partial τ} - V_{n} |\nabla Φ| = 0,$ (15)

where $V_{n} (x, τ) = - \frac{{(\nabla Φ)}^{T}}{|\nabla Φ|} V$ is the normal velocity field.

The MultiQuadric (MQ) spline, a type of RBFs, is used here to construct the level set function, and can be written as³⁶ $g_{i} (x) = \sqrt{{|| x - x_{i}^{k n o t}| |}^{2} + c^{2}}, x_{i}^{k n o t} \in \bar{D},$ (16)

where $|| \cdot| |$ denotes the norm of a vector, $x_{i}^{k n o t} = {[\begin{matrix} x_{i}^{k n o t} & y_{i}^{k n o t} \end{matrix}]}^{T}$ defined in the 2D coordinate system is the position vector of the coordinates of the ith knot, and c is the corresponding shape parameter, which is commonly assumed to be a constant with a small value for all knots. In this study, c is set to be $1 0^{- 4}$ .³⁶

Then, the level set function at point x is interpolated by n MQ splines, which can be written as

By substituting Eq. (17) in Eq. (15), a set of simpler ODEs can be obtained as follows:

where $G = [\begin{matrix} A & L \\ L^{T} & 0 \end{matrix}],$ (19)

A = [\begin{matrix} g_{1} (x_{1}^{k n o t}) & \dots & g_{m} (x_{1}^{k n o t}) \\ ⋮ & ⋱ & ⋮ \\ g_{1} (x_{m}^{k n o t}) & \dots & g_{m} (x_{m}^{k n o t}) \end{matrix}],

(20)

L = [\begin{matrix} 1 & x_{1}^{k n o t} & y_{1}^{k n o t} \\ ⋮ & ⋮ & ⋮ \\ 1 & x_{m}^{k n o t} & y_{m}^{k n o t} \end{matrix}],

(21)

\frac{\partial g (x)}{\partial x} = {[\begin{matrix} \frac{\partial g_{1} (x)}{\partial x} & \dots & \frac{\partial g_{m} (x)}{\partial x} & 0 & 1 & 0 \end{matrix}]}^{T},

(24)

\frac{\partial g (x)}{\partial y} = {[\begin{matrix} \frac{\partial g_{1} (x)}{\partial y} & \dots & \frac{\partial g_{m} (x)}{\partial y} & 0 & 0 & 1 \end{matrix}]}^{T} .

(25)

The detailed derivation of Eq. (18) can be found in Section S2 in Supplementary Data S1. It should be noted that the original Hamilton-Jacobi equation has been reduced to a simpler coupled ODEs.

Topology optimization formulation

The purpose of the topology optimization of the electrode is to improve the driving performance of DEAs. Based on the parameterized level set method, the formulation for the topology optimization of the electrode of a mechanical system actuated by DEs can be mathematically expressed as

where denotes the objective function, V denotes the volume of the soft body distributed with the electrode, f is the prescribed allowable volume fraction, and $|D|$ is the volume of the design domain D.

ESL Based Optimization

System-wise ESL method

To solve the dynamic topology optimization problem for a mechanical system as shown in Eq. (26), the system-wise ESL method is utilized here to convert the dynamic response optimization into a static one. This enables the optimization process of complicated dynamic problems to be computationally feasible but at the cost of accuracy loss.

Based on ANCF, the ESL sets are defined as the static loads that generate the same deformations in static analysis as those from dynamic analysis at an arbitrary time step at the system level. According to the work by Tromme et al.,²⁷ the system-wise ESL sets are computed by $F_{i}^{e q} = \int_{D} {(\tilde{B} - ρ {\ddot{r}}_{i})}^{T} d D + \int_{S_{t}} {\tilde{T}}^{T} d S,$ (27)

where subscript i denotes the selected ith time step, $i = 1, \dots, N_{t}$ , and N_t is the number of selected time steps. $F_{i}^{e q}$ is used as the external force vector for the static problem at the selected ith time step. Here, the continuous inertia forces are transformed into the discrete static loads at selected time steps, so the inertia effect is neglected or only partially involved when the time step is small enough. It should be noted that this approach may not be appropriate if the inertia effect becomes dominant or high precision in control of soft robots is required. However, for this study, the precision loss from this approximation is acceptable.

The static problem at the selected ith time step is defined as follows: $\{\begin{matrix} \int_{D} P_{i}^{m} : δ F_{i} d D + \int_{Ω_{0}} P_{i}^{e} : δ F_{i} d Ω_{0} = F_{i}^{e q} δ u_{i}, \\ u = \bar{u}, o n S_{u}, \forall δ u_{i} \in U_{a d} . \end{matrix}$ (28)

In Eq. (28), $i = 1, \dots, N_{t}$ , and N_t is the number of time steps.

In this work, the expert evaluation method of weight based on the Grey theory³⁷ is utilized to define the objective function as follows:

where N_t is the number of time steps and $⨂ w_{i}$ is the grey number to denote the weight for the objective function at the ith time step, which is defined based on the expert evaluation method of weights.³⁸

Based on the system-wise ESL method described earlier, the dynamic response optimization can be converted into the following static one:

where $i = 1, \dots, N_{t}$ , and N_t denotes the number of time steps.

Sensitivity analysis

In the optimization process, topology changes of the electrode are achieved based on the evolution of the electrode boundary, and the sensitivity analysis of the objective function needs to be calculated. In this subsection, a sensitivity analysis is carried out by conducting the material derivatives combined with the adjoint sensitivity method.³⁹ The Lagrangian of the optimization problem in Eq. (30) is constructed by

According to the adjoint sensitivity method,³⁹ the first-order variation of the objective function (with respect to a pseudo time variation $δ τ$ ) can be calculated by

where ${(\cdot)}^{'} = \frac{\partial (\cdot)}{\partial τ}$ denotes the partial derivative of $(\cdot)$ with respect to $τ$ .

In Eq. (32), $\begin{matrix} δ L_{i}^{a} = \int_{D} (H_{i}^{m} : F_{i}^{'} : {\bar{F}}_{i} (v_{i})) δ τ d D \\ + \int_{Ω_{0}} (H_{i}^{e} : F_{i}^{'} : {\bar{F}}_{i} (v_{i})) δ τ d Ω_{0} + ⨂ w_{i} {(\frac{\partial J_{i}}{\partial u_{i}})}^{T} u_{i}^{'} δ τ, \end{matrix}$ (33)

δ L_{i}^{s} = \int_{D} (P_{i}^{m} : \nabla v_{i}^{'}) δ τ d D + \int_{Ω_{0}} (P_{i}^{e} : \nabla v_{i}^{'}) δ τ d Ω_{0} - F_{i}^{e q} v_{i}^{'} δ τ,

(34)

where $i = 1, \dots, N_{t}$ . $H^{m} = \frac{\partial P^{m}}{\partial F}$ and $H^{e} = \frac{\partial P^{e}}{\partial F}$ are the fourth-order material elasticity tensors corresponding to P ^m and P ^e , respectively. The terms $H^{m}$ and $H^{e}$ can be viewed in Section S3 in Supplementary Data S1.

It is obvious that the equation $δ L_{i}^{s} (u_{i}, v_{i}^{'}) = 0, (i = 1, \dots, N_{t}),$ (36)

denotes the virtual work principle at the ith time step, and it is exactly the same as Eq. (28).

The adjoint displacement v _i at the ith time step is obtained by solving the adjoint equation at the ith time step,³⁹ that is, $δ L_{i}^{a} (u_{i}, u_{i}^{'}, v_{i}) = 0, (i = 1, \dots, N_{t}),$ (37)

with the same displacement boundary conditions as that of the static problem in Eq. (28).

Finally, the first order variation of the objective function can be calculated by $δ J = δ L^{r} = \sum_{i = 1}^{N_{t}} \int_{\partial Ω_{0}} (P_{i}^{e} : {\bar{F}}_{i}) V_{n} δ τ d S + λ \int_{\partial Ω_{0}} V_{n} δ τ d S,$ (38)

Further, the steepest ascent direction is directly obtained by setting $V_{n} = \sum_{i = 1}^{N_{t}} (P_{i}^{e} : {\bar{F}}_{i}) + λ .$ (39)

Topology optimization computational procedure

The overall flowchart of the topology optimization of the electrode for soft multibody systems actuated by DEs using the parameterized level set method is established in Figure 3, where the steps of the algorithm are briefly listed as follows.

FIG. 3.

The flowchart of topology optimization of the electrode for soft multibody systems actuated by DEs.

Step 1. Initialize the design variables and set the iteration counter K to 0.

Step 2. Establish the dynamic equation of soft multibody systems actuated by DEs based on ANCF.

Step 3. Solve the dynamic equation and compute the system-wise ESL sets, then the dynamic response optimization is converted into a static one.

Step 4. The static response optimization is conducted based on the parameterized level set method and the adjoint sensitivity method. The number of iterations to solve the static response optimization is denoted as k.

Step 5. Estimate the convergence criterion. If the convergence criterion is satisfied, terminate the optimization process; otherwise, update the design variables, $K = K + 1$ and go to step 3.

Case Studies

This section presents two case studies to illustrate the effectiveness and validity of the proposed method for the topology optimization of the electrode of DE components in soft multibody systems, that is, the static and dynamic response optimization of a planar DEA and the dynamic response optimization of a revolute-spherical-spherical-revolute (RSSR) mechanism.

Static and dynamic response optimization of a planar DEA

First, the static response optimization of a planar DEA is considered as shown in Figure 4a. The static response optimization of the planar DEA has also been studied by Chen et al.¹⁸ The design domain is the planar DEA, whose geometric parameters are set to be $L = 40 m m$ and $H = 1 m m$ . The design domain is discretized by 20 × 20 DE plate elements of ANCF.

FIG. 4.

(a) Geometric schematic of the planar DEA; (b) initial topology; (c) optimized topology.

J = u_{O} \cdot T_{O},

(40)

where $u_{O}$ denotes the nodal displacement of the point O, and the transformation vector is $T_{O} = {[\begin{matrix} 1 & 0 \end{matrix}]}^{T}$ .

A symmetric initial topology is given in Figure 4b, where many circles are contained to allow for the sufficient geometric flexibility for potential topological evolutions. The optimized topology is shown in Figure 4c, which is obtained at step 150, and the corresponding objective function is 0.3215 mm.

The iteration histories of the volume fraction and the objective function in the optimization process are shown in Figure 5a and b, respectively. It can be seen that the objective function (the X-displacement of the point O) is improved; meanwhile, the volume constraint is satisfied.

FIG. 5.

Iteration histories of the optimization process: (a) volume fraction; (b) objective function.

To further validate the proposed optimization methodology, the dynamic response optimization of the planar DEA described earlier is performed. A sinusoidal voltage is applied to the planar DEA, which is set to be $ϕ = 2.5 sin (\frac{2 π}{0.1} t) k V$ . The analysis time is set from 0 to 0.1 s, and the size of the time step is 0.01 s. Thus, there are 10 time steps altogether, that is $N_{t} = 10$ . It means that the number of the loading conditions is 10, too. The regularized whitening values of the weights in Eq. (29), that is, $⨂ w_{i}$ , are defined as follows²⁶ $\begin{matrix} w = [\begin{matrix} 0.06208 & 0.08719 & 0.10488 & 0.13951 & 0.13270 \end{matrix} \\ {\begin{matrix} 0.15044 & 0.10419 & 0.09658 & 0.06814 & 0.05429 \end{matrix}]}^{T} . \end{matrix}$ (41)

The optimized topology for the dynamic response optimization is shown in Figure 6a. As shown in Figure 6a, the optimized topology for the dynamic response optimization of the planar DEA is almost the same as the static response optimization result shown in Figure 4c. The reason is that the dynamic behavior of the planar DEA takes place near the equilibrium position, which is very close to the static behavior.

FIG. 6.

(a) Optimized topology for the dynamic response optimization; (b) iteration history of the objective function; (c) experiment setup; (d) planar DEA; (e) comparison between the numerical and experimental results of a planar DEA with 40% electrode volume.

Figure 6b presents the iteration history of the objective function. As shown in Figure 6b, the topology optimization procedure converges fast within a few optimization iterations, and the optimized topology relies on that in the first iteration very much. It should be pointed out that the number of iterations refers to that based on the system-wise ESL method, that is, the number of outer iterations in Figure 3.

Further, an experimental test rig is also established to validate the results of the dynamic simulation. As shown in Figure 6c, the sinusoidal voltage is generated by a signal generator (Tektronix AFG31000 series) and amplified by a voltage amplifier (Aigtek ATA-7025). The deformation of the DEA is recorded by a high-speed camera (Protron FASTCAM SA5). The DEA used in the experiment is VHB 4910 produced by the 3M company.

As shown in Figure 6d, a piece of VHB 4910 sheet with equal-biaxial pre-stretches of two is supported by a stiff acrylic frame, and the compliant electrodes of carbon grease are coated on the surfaces of the elastomer sheet. From Figure 6e, it is shown that the numerical results match well with the experimental results, which indicates a good prediction of the DEA dynamics.

Finally, the effectiveness of the proposed optimization approach is further demonstrated in terms of the actuation performance of DEAs. The comparison is made on dynamic responses of the planar DEAs with the initial topology, an intermediate topology during iteration, and the optimized topology of electrodes. Here, our main concern is to evaluate the actuation performance of the DEA and the effects caused by different topologies need to be well distinguished by our experimental measurement. Thus, the volume constraint of electrodes is released.

And the optimization objective is to maximize the dynamic distance between points A and B at equal initial distances of 8 mm from the origin point O on Y axis, as shown in Figure 4a. A sinusoidal voltage of $ϕ = 3 sin (2 π t) k V$ is applied to the planar DEAs. This more proper voltage signal with larger amplitude and smaller frequency is used to recognize the small actuation by the initial topology and the intermediate topology in the experiment.

The simulation time is set from 0 to 1 s, and the size of time step is 0.1 s. And the corresponding experiments are also performed. The experimental setup is the same as shown in Figure 6c. Geometric parameters of the square DEA are set to be lengths $L = 32 m m$ and thickness $H = 1 m m$ . The boundary of the DEA is fixed under equal-biaxial pre-stretches of 2.5. Then, the lengths L of the DEA and the initial distance between points A and B become 80 and 40 mm, respectively. As shown in Figure 7, though deviations between the numerical and experimental results exist for different topologies, the evidence is clear that the DEA with the optimized electrode topology provides the best actuation on the relative displacement of the two points.

FIG. 7.

Comparison of dynamic responses of planar DEAs without volume constraint but with three different topologies of electrodes: (a) numerical results; (b) experimental results.

Optimization of dynamic responses of a RSSR mechanism

The RSSR mechanism is considered to validate the proposed topology optimization methodology. As shown in Figure 8, the RSSR mechanism consists of four links⁴⁰: the fixed link, the crank shaft, the coupler, and the rocker. The crank is connected to the ground (fixed link) by a revolute joint at O and to the coupler by a spherical joint at A. The rocker is connected to the coupler by a spherical joint at B and to the fixed link by a revolute joint at E. As shown in Figure 8, the point P is the center point of the coupler.

FIG. 8.

Geometric schematic of the RSSR mechanism. RSSR, revolute-spherical-spherical-revolute.

The dynamic response optimization of the RSSR mechanism actuated by DEs is performed. The model parameters for the RSSR mechanism used in the optimization are presented in Table 1. The design domain is the coupler, which is discretized by 36 × 18 DE plate elements of ANCF. The volume fraction is set to be $f = 0.4$ . The objective function is set to maximize the elastic strain energy in the design domain, that is,

Table 1.

Model Parameters for the Revolute-Spherical-Spherical-Revolute Mechanism Used in the Optimization

Links	Crank shaft	Coupler	Rocker
Dimensions (mm)	21.33 × 4.19 × 2.10	80 × 40 × 2	53.33 × 0.42 × 0.21
Locations of body reference (mm)
x	5.333	37.333	64.000
y	9.239	25.806	16.567
z	0.0000	28.892	36.892
Orientations of body reference in Euler angles (rad)
φ	1.0470	0.2680	1.5708
θ	0.0000	−1.5708	1.5708
ψ	0.0000	−0.8074	0.9007
Material properties	Rigid	μ = 40 kPa, K = 3 GPa, μ^neq = 90 kPa, τ^neq = 200 s, ɛ = 3.6 × 10⁻¹¹ F/m, ρ = 1000 kg/m³	Rigid

J = \int_{D} ψ^{e q} d D,

(42)

The angular velocity of the crank shaft is defined as follows: $ω = \{\begin{matrix} 10 π t, i f 0 \leq t \leq 0 . 1; \\ π, i f 0 . 1 < t \leq 0 . 3 . \end{matrix}$ (43)

The analysis time is set from 0.1 to 0.3 s, and the size of the time step is 0.02 s. Thus, there are 10 time steps altogether, that is $N_{t} = 10$ . Figure 9 illustrates the optimized topologies for the 10 time steps in the first iteration. It is easy to identify that the optimized topologies of the electrode for the coupler corresponding to different time steps are different.

FIG. 9.

Optimized topologies of the electrode from (a) to (j) at 10 time steps in the first iteration.

To study the influence of regularized whitening values of the weights on the optimized topology of the electrode, two sets of regularized whitening values are selected in the optimization process. One set of regularized whitening values is w in Eq. (41), and the other set of regularized whitening values is given by $w_{1} = {[\begin{matrix} 0.1 & 0.1 & 0.1 & 0.1 & 0.1 & 0.1 & 0.1 & 0.1 & 0.1 & 0.1 \end{matrix}]}^{T} .$ (44)

The optimized topologies of the electrode corresponding to two different weights, that is, w and $w_{1}$ , are presented in Figure 10a and b, respectively. As shown in Figure 10a and b, there is some difference between the optimized topologies corresponding to different weights. In other words, the choices of the weights affect the optimized topology of the electrode. Dynamic configurations of the RSSR mechanism with applied DEA at t = 0.2 s are presented in Figure 10c and d.

FIG. 10.

The final optimized topologies corresponding to different weights: (a) w; (b) $w_{1}$ . And dynamic configurations of the RSSR mechanism with optimized DEA at t = 0.2 s with: (c) w; (d) $w_{1}$ . (e) Iteration history of the objective function. (f) Z-displacement of the point P.

The iteration history of the objective function is shown in Figure 10e. The objective function values corresponding to different weights are different. The Z-displacement of the point P corresponding to the optimized topologies is shown in Figure 10f. From Figure 10f, the values of the weights have little impact on the Z-displacement of the point P.

Conclusions

Based on the ANCF nonlinear finite elements, the parameterized level set method, and the system-wise ESL method, a computational methodology for the topology optimization of electrode in a mechanical system actuated by DEs is proposed. A plate element of ANCF based on the constitutive model of DEs is established to perform the dynamic analysis of soft multibody systems actuated by DEs.

The parameterized level set method is used to describe the electrode topology and formulate the optimization problem. Based on the system-wise ESL method, the dynamic response topology optimization problem is converted into the static one. The adjoint sensitivity analysis is performed to derive the normal velocity field of the electrode topology. Two numerical examples are presented to validate the capability and effectiveness of the proposed dynamic optimization method for the electrode topology. First, the static and dynamic response optimization of a planar DEA is studied, and numerical results of the dynamic behavior for the optimized topology agree well with the experimental results.

Subsequently, the dynamic response optimization of an RSSR mechanism is conducted to validate the proposed method, and the optimization results indicate that the choices of the weights affect the optimized topology of the electrode. The application of the proposed method to the design of soft robots will be further investigated in our future studies.

Footnotes

Author Disclosure Statement

No competing financial interests exist.

Funding Information

This work is supported by the National Natural Science Foundation of China under grant nos. 12125201, 11832005, and 11902028.

Supplementary Material

References

O'Halloran

, O'Malley

, Mchugh

. A review on dielectric elastomer actuators, technology, applications, and challenges. J Appl Phys, 2008; 104:071101.

Liu

, Urban

. Recent advances and challenges in designing stimuli-responsive polymers. Prog Polym Sci, 2010; 35:3–23.

Romasanta

, Lopez-Manchado

, Verdejo

. Increasing the performance of dielectric elastomer actuators: a review from the materials perspective. Prog Polym Sci, 2015; 51:188–211.

Brochu

, Pei

. Advances in dielectric elastomers for actuators and artificial muscles. Macromol Rapid Commun, 2010; 31:10–36.

Feher

, Filipcsei

, Szalma

, et al. Bending deformation of neutral polymer gels induced by electric fields. Colloids Surf A Physicochem Eng Asp, 2001; 183:505–515.

Suo

. Theory of dielectric elastomers. Acta Mech Solida Sin, 2010; 23:549–578.

, Zhu

, et al. A survey on dielectric elastomer actuators for soft robots. Bioinspir Biomim, 2017; 12:011003.

Hossain

, Vu

, Steinmann

. Experimental study and numerical modelling of VHB 4910 polymer. Comput Mater Sci, 2012; 59:65–74.

Tian

, Zhang

, Luo

. Dynamics of soft mechanical systems actuated by dielectric elastomers. Mech Syst Signal Process, 2021; 151:107392.

10.

Shabana

. Flexible multibody dynamics review of past and recent development. Multibody Syst Dyn, 1997; 1:189–222.

11.

Shabana

, Zhang

. ANCF curvature continuity: application to soft and fluid materials. Nonlinear Dyn, 2020; 100:1497–1517.

12.

Shabana

, Xu

. Rotation-based finite elements: reference-configuration geometry and motion description. Acta Mech Sin, 2021; 37:105–126.

13.

Shabana

. Continuum-based geometry/analysis approach for flexible and soft robotic systems. Soft Robot, 2018; 5:613–621.

14.

Shintake

, Rosset

, Schubert

, et al. Versatile soft grippers with intrinsic electroadhesion based on multifunctional polymer actuators. Adv Mater, 2016; 28:231–238.

15.

Zou

, Li

, Qu

, et al. Active shape control and phase coexistence of dielectric elastomer membrane with patterned electrodes. J Appl Mech, 2014; 81:031016.

16.

Wang

, Guo

, Chen

, et al. Design of a rotary dielectric elastomer actuator using a topology optimization method based on pairs of curves. Smart Mater Struct, 2018; 27:055011.

17.

Wang

, Guo

, Chen

, et al. Design of dielectric elastomer actuator using topology optimization method based on genetic algorithm. Smart Mater Struct, 2019; 28:065013.

18.

Chen

, Liu

, Wang

, et al. Automatic design of soft dielectric elastomer actuators with optimal spatial electric fields. IEEE Trans Robot, 2019; 35:1150–1165.

19.

Sethian

, Wiegmann

. Structural boundary design via level set and immersed interface methods. J Comput Phys, 2000; 163:489–528.

20.

Wang

, Wang

. Radial basis functions and level set method for structural topology optimization. Int J Numer Methods Eng, 2010; 65:2060–2090.

21.

Wang

, Wang

. Structural shape and topology optimization using an implicit free boundary parameterization method. Comput Model Eng Sci, 2006; 13:119–147.

22.

Wang

, Lim

, Khoo

, et al. An extended level set method for shape and topology optimization. J Comput Phys, 2007; 221:395–421.

23.

Kang

, Choi

, Park

. Structural optimization under equivalent static loads transformed from dynamic loads based on displacement. Comput Struct, 2001; 79:145–154.

24.

Hong

, You

, Kim

, et al. Optimization of flexible components of multibody systems via equivalent static loads. Struct Multidiscipl Optim, 2010; 40:549–562.

25.

Kang

, Park

, Arora

. Optimization of flexible multibody dynamic systems using the equivalent static load method. AIAA J, 2005; 43:46–852.

26.

Sun

, Tian

, Hu

. Topology optimization based on level set for a flexible multibody system modeled via ANCF. Struct Multidiscipl Optim, 2017; 55:1159–1177.

27.

Tromme

, Sonneville

, Guest

, et al. System-wise equivalent static loads for the design of flexible mechanisms. Comput Methods Appl Mech Eng, 2018; 329:312–331.

28.

Mikkola

, Shabana

. A non-incremental finite element procedure for the analysis of large deformation of plates and shells in mechanical system applications. Multibody Syst Dyn, 2003; 9:283–309.

29.

Dufva

, Shabana

. Analysis of thin plate structures using the absolute nodal coordinate formulation. Proc Inst Mech Eng K, 2005; 219:345–355.

30.

Pelrine

, Kornbluh

, Joseph

. Electrostriction of polymer dielectrics with compliant electrodes as a means of actuation. Sens Actuators A Phys, 1998; 64:77–85.

31.

Shabana

. Computational Continuum Mechanics. New York: Cambridge University Press, 2012.

32.

Colonnelli

, Saccomandi

, Zurlo

. The role of material behavior in the performances of electroactive polymer energy harvesters. J Polym Sci Part B Polym Phys, 2015; 53:1303–1314.

33.

Luo

, Rothemund

, Whitesides

, et al. Soft kink valves. J Mech Phys Solids, 2019; 131:230–239.

34.

Luo

, Tian

, Hu

. Dynamic modeling, simulation and design of smart membrane systems driven by soft actuators of multilayer dielectric elastomers. Nonlinear Dyn, 2020; 102:1463–1483.

35.

Bower

. Applied Mechanics of Solids. Boca Raton: CRC Press, 2009.

36.

Wei

, Li

, et al. An 88-line MATLAB code for the parameterized level set method based topology optimization using radial basis functions. Struct Multidiscipl Optim, 2018; 58:831–849.

37.

Deng

. Introduction to grey system theory. J Grey Syst, 1989; 1:1–24.

38.

Luo

, Chen

, Zhang

, et al. A new hybrid fuzzy-goal programming scheme for multi-objective topological optimization of static and dynamic structures under multiple loading conditions. Struct Multidiscipl Optim, 2006; 31:26–39.

39.

Choi

, Kim

. Structural Sensitivity Analysis and Optimization 1: Linear Systems. New York, NY: Springer, 2006.

40.

Chang

, Shabana

. Nonlinear finite element formulation for the large displacement analysis of plates. J Appl Mech, 1990; 57:707–718.

Supplementary Material

Please find the following supplemental material available below.

For Open Access articles published under a Creative Commons License, all supplemental material carries the same license as the article it is associated with.

For non-Open Access articles published, all supplemental material carries a non-exclusive license, and permission requests for re-use of supplemental material or any part of supplemental material shall be sent directly to the copyright owner as specified in the copyright notice associated with the article.

0.24 MB

0.10 MB

0.00 MB