Electromechanical Coupling in Electroactive Polymers – a Visual Analysis of a Third-Order Tensor Field

Electroactive polymers are frequently used in engineering applications due to their ability to change their shape and properties under the influence of an electric field. This process also works vice versa, such that mechanical deformation of the material induces an electric field in the EAP device. This specific behavior makes such materials highly attractive for the construction of actuators and sensors in various application areas. The electromechanical behaviour of electroactive polymers can be described by a third-order coupling tensor, which represents the sensitivity of mechanical stresses concerning the electric field, i.e., it establishes a relation between a second-order and a first-order tensor field. Due to this coupling tensor's complexity and the lack of meaningful visualization methods for third-order tensors in general, an interpretation of the tensor is rather difficult. Thus, the central engineering research question that this contribution deals with is a deeper understanding of electromechanical coupling by analyzing the third-order coupling tensor with the help of specific visualization methods. Starting with a deviatoric decomposition of the tensor, the multipoles of each deviator are visualized, which allows a first insight into this highly complex third-order tensor. In the present contribution, four examples, including electromechanical coupling, are simulated within a finite element framework and subsequently analyzed using the tensor visualization method.


INTRODUCTION
E LECTROACTIVE polymers (EAP) feature remarkable proper- ties.They change their shape and other mechanical properties like stiffness under the influence of an electric field.If the electric field vanishes, EAP such as dielectric polymers (DE) resume their original shape and properties.This electromechanical coupling also works vice versa, so changing an EAP's shape influences the electric field.Therefore, EAP facilitate flexible and lightweight smart devices such as actuators and sensors, artificial muscles, energy harvesting systems, lens cleaning, and acoustic devices [1], [2], [3], [4], [5], [6], [7].EAP muscles may be much stronger and faster than biological muscles relative to their weight.Leading companies in the USA, Germany, France, Japan, Saudi-Arabia, and other countries have produced EAPs worth more than $4:5 billion in 2019, according to Market Reports World [8].
Dielectric elastomers (DE) are a particular class of EAPin contrast to ionic polymers, which show only little electromechanical coupling.Electrically activated DE maintain the induced deformation which, however, requires high voltage levels even close to electrical breakdown.Common DE devices are made of silicone or are acrylic-based and allow a fast activation.A prominent example of a dielectric EAP material is VHB produced by 3M which may undergo large strains of more than 500 %.Typical VHB based EAP devices include a thin elastomer layer placed between two thin electrode layers.Upon activation, the electrodes contract so that the elastomer increases its in-plane areal stretch.The material behavior of VHB shows rate-dependent effects, such as electro-viscoelastic response, as well as temperature and moisture dependencies.As this work proceeds, finite deformation electroelasticity shall be considered, whereby a Gent-type model together with material parameters identified for VHB49 are adopted from [9].For overviews on the mechanics of EAP, the reader is referred to, e.g., [10], [11], whereas experimental investigations of VHB are addressed in, e.g., [12], [13], [14].Further model extensions towards the simulation of electro-viscoelastic behaviour are discussed in [9], [15], amongst others.
EAP design is a most challenging and interdisciplinary research field in mechanics, material sciences, biotechnology, and robotics [16].Predictive modeling of EAPs and advanced simulations of EAP based devices contribute to a systematic design and optimization of related smart systems [17].Such design and optimization require the interpretation of the simulation results at material and structural level.The present paper enables such interpretation by advanced visualization techniques of higher-order tensor fields, which represent the (local) material properties.In particular, a third-order tensor H shall be elaborated, which represents the sensitivity of mechanical stresses with respect to the electric field, respectively the sensitivity of dielectric displacements with respect to strains -in other words, the coupling between mechanical and electrical material properties.As we shall elaborate in the following, this field H is inhomogeneously distributed for the boundary value problems, respectively applications considered in this work, which is of major interest in order to understand and illustrate electromechanical coupling as well as its distribution present in EAP based devices.Although it can be computed by suitable finite element formulations, there is, to the best of our knowledge, basically no visualization method available in the literature for electromechanical coupling.
This article intends to change this situation by giving some insight into the use of a glyph visualization for this third-order tensor field.Even though the represented decomposition and the visualization can be used for all three-dimensional tensors of any order, the focus of this work will be the third-order coupling tensor H and its application.To design the glyph, the deviatoric decomposition is used.This decomposition is an irreducible decomposition.Thus, a rotation of the tensor equals the result of the sum of all irreducible parts after the same rotation, and there is no further decomposition of these parts with this property.Accordingly, glyphs that represent tensors equaling each other in a local coordinate system are the same.By analyzing the deviators separately, the number of coefficients to be considered simultaneously increases linearly and not exponentially in comparison to the order of the tensor.
We analyze four different problems to gain a first idea for an interpretation of the third-order coupling tensor H.A radial symmetric cylinder is an illustrative first example because its symmetry simplifies the problem and leads to well interpretable results.The second example is closely related to applications in the field of dielectric elastomer actuators (DEAs).Specifically speaking, a bending beam, which is a component frequently used for gripper-type applications, is analyzed by using the proposed visualization approach.In the third example, the focus is on imperfections, which are naturally present in every material, and their influence on the electromechanical behavior.To this end, a spherical defect in the centre of a cubic region is presented.Since the knowledge on the inhomogeneous distribution of electromechanical coupling response is rather limited, the newly gained knowledge through the visualization approach shall allow for a better interpretation of this example.The last example is from the field of robotics: a bioinspired lens with an electrically tunable curvature is analyzed.Within the scope of this work, a perfectly circular lens as well as a lens with a geometric imperfection, a slightly elliptic shape, are both simulated in order to compare the visualization results.

QUESTIONS
Electromechanically coupled problems are highly complex due to the, in general, highly non-linear behavior under finite deformations.However, understanding electromechanical coupling is very important from an engineering point of view since electromechanical actuators and sensors are the basis of a huge number of smart devices and other applications.Through a better understanding of the material behavior in such coupled settings, the improvement of the design of EAP-based actuators and sensors is the general goal.
Thus, an investigation of the third-order tensor H, which describes the coupling in electromechanical problems and consequently provides information on the (local) material behavior under load, is of considerable interest.For non-linear problems, such material tensors are not constant but typically change under loading.From an engineering and material science perspective, little information regarding these tensors is currently available.There are a few examples where the properties of the coupling tensor are understood in more detail, but once the examples become more complex, further analysis methods are required.A reason is the complexity of higher-order tensors in general and the missing techniques for visualization and analysis of higherorder, and especially odd-order tensors.
Besides the field of electromechanical coupling, which is the main focus of this work, there are many other coupled problems such as magneto-and thermomechanics, which give rise to similar questions.
The goal of our visualization is to support the understanding of the material behavior within coupled problems by analyzing the respective coupling tensor in detail.We also focus on finding a connection of the higher-order tensors to related tensors of lower order in order to simplify the tensor representations.
All in all, the main goal is to better understand and visualize material tensors for coupled problems, in particular electromechanics, as this work proceeds.The questions can be summarized as follows I. How does the material behave?II.How can a third-order tensor be analyzed?III.How can a third-order tensor be simplified?IV.How can known knowledge from simple examples be transferred to more difficult ones? V. How can the design of EAP-based actuators and sensors be improved?

MODEL
We used four examples to analyze the electromechanical coupling in electroactive polymers, which are introduced in the following.The first one is well established, for the second one, some information, as well as applications are also available, and the third one is more complex and less information is available so far.An analysis of the first and the second example should facilitate the analysis of the third one.The fourth example indicates how the visualization can be an indicator of the functionality of electromechanical devices.

Electrically Activated Cylinder
As a first example, a cylinder is analyzed subject to electric loading, cf.[18].In contrast to the VHB material used for the following examples, a compressible material behavior represented by m ¼ 1:1233 Â 10 4 Pa and k ¼ 2 m is achieved together with plane strain conditions so that the related Poisson ratio corresponds to n ¼ 0:2857.At the outer radius r o ¼ 2 mm of the cylinder, an electric potential of f o ¼ À5 kV has been applied, whereas the potential at the inner radius r i ¼ 1 mm is f i ¼ 5 kV, cf.Fig. 1.The thickness of the cylinder is t ¼ 0:1 mm.Due to the symmetry properties of the structure, it is sufficient to model one quarter of the cylinder together with additional symmetry boundary conditions.In particular, the displacement normal to the symmetry axes is suppressed on both boundaries.In addition, the inner radius of the cylinder is kept fixed, whereas the outer radius is allowed to change.The mesh shown in Fig. 1 consists of n el ¼ 105 elements.
As a consequence of the potential difference Df ¼ f o À f i , the cylinder is compressed in the radial direction.
Due to the particular boundary conditions, this leads to a reduction of the outer radius.The deformation in the last step of the analysis is shown by the colored region in Fig. 1.In addition, in the supplementary material, it is demonstrated how the outer cylinder radius decreases with increasing difference in electric potential, which highlights the overall non-linear response.

Electrically Activated Bending of a Beam
The next example deals with the bending deformation of a beam due to electric loading.Beam-like structures made from electroactive polymers are frequently used in applications with dielectric actuators.Fig. 2 shows the working principle of such an actuator under electric loading.Besides the frequent use in the context of gripper-type devices, cf.[19], dielectric elastomer actuators (DEAs) are also an important component in robotics, see Fig. 3.
For the beam model, a VHB material is employed.It's (quasi) incompressible material behavior shall be represented by m ¼ 1:1233 Â 10 4 Pa and k ¼ 1:1196 Â 10 6 Pa so that the related Poisson ratio corresponds to n ¼ 0:495, cf.[9].In the present contribution, a beam with a length of l ¼ 50 mm, a height of h ¼ 1 mm, and a width of t ¼ 4 mm is considered.The discretization of the model contains n el ¼ 1472 elements.As shown in Fig. 4, only the upper half of the beam is loaded by applying a difference in electric potential.The difference is set to Df ¼ 2 kV.Regarding the mechanical boundary conditions, the right boundary is fixed in the direction of the beam's longitudinal axis.In the other two spatial directions, the boundary nodes are free to move, except for the middle node.The electric load leads to a compression of the elastic material in the upper half of the beam in the direction of the electric field.In addition, an extension in the longitudinal direction of the beam and in the direction of the width of the beam follows because a nearly incompressible material is considered.In the lower half of the beam, no electric field is active.As a consequence of this particular loading condition, together with the mechanical boundary conditions, the structure undergoes a bending deformation because the upper half of the structure is extended, cf.Fig. 4. The vertical displacement u y as well as the horizontal displacement u z of the tip of the beam are shown in dependence on the electric potential differencies displayed in the supplementary material.

Electrically Activated Cube with Spherical Hole
In this subsection, a material that contains defects or inclusions is represented, whereby identical material parameters as in Section 3.2 are chosen.In particular, an EAP cube with one single defect in the shape of a spherical hole is modeled.Fig. 5 shows the back half part of the cube, whereas the front half part can be generated by a reflection of the sketched domain.Due to its symmetry properties, only one eighth of the cube with an edge length of 2 h ¼ 2 mm is considered in the simulation.The discretization of this part is obtained with 783 elements.The radius of the spherical hole in the center of the cube is r ¼ 0:25 mm.Similar to the example in Section 3.1, a potential difference of Df ¼ 10 kV is employed by the application of an electric potential of f l ¼ 5 kV on the left surface of the cube and of f r ¼ À5 kV on the opposite surface.In the vertical symmetry plane, f ¼ 0 kV holds accordingly.In order to capture the symmetry in the simulation, linear constraints are employed for the displacements on the three symmetry surfaces.To be specific, no displacement in the normal direction of these boundary surfaces is allowed to take place.
The deformation pattern of the cube is shown by the colored region in Fig. 5.The cube is compressed in the direction of the electric field.As a consequence of the (quasi) incompressibility of the material, an extension can be observed in the other spatial directions.The respective displacement values u y and u z are presented exemplarily for the upper left corner of the simulated domain in the supplementary material.In accordance with the deformation of the cube, the defect changes its shape from a sphere to an ellipsoid.

Bioinspired Tunable Lens
Using the same EAP material as in Sections 3.2 and 3.3, a bioinspired tunable lens is modeled.In an experimental setting, as shown in Fig. 6, a prestreched EAP material is placed into a frame and filled with a fluid to obtain a symmetric, biconvex lens, cf.[22].The outer region of the lens is coated with black electrodes.An electric activation causes the material to be compressed in the direction of the thickness.Due to the (nearly) incompressibility of the EAP material, the material accordingly deforms radially towards the center, and the curvature of the lens changes, see Fig. 6.For the numerical simulation, one eighth of the lens is modeled by exploiting its symmetry properties.Consequently, homogeneous Dirichlet boundary conditions are applied to the symmetry planes in order to prevent any displacement into the normal direction of the respective surfaces.In the first part of the simulation, a prestretch is applied in terms of inhomogeneous Dirichlet boundary conditions to the outer edge of the domain.By analogy with the experiment, this prestretch has the advantage of preventing buckling modes that might occur in response to the material compression in the simulation, and it reduces the required electric load since the electrode distance decreases.The  deformed mesh after application of the prestretch is presented in Fig. 7 and consists of n el ¼ 5628 finite elements.Similar to the experimental setting in [22], the outer radius after the prestretch is 10:5 mm and the radius of the lens in the xz-plane is % 3:8 mm.The thickness of the material is approximately 0:1 mm in the prestretched state.Since only the electromechanical behavior is of interest in the context of the present contribution, the fluid filling is not included in the simulation so that the lens remains rather flat before the electric load is applied.In the second part of the simulation, the electromechanical analysis is performed while the displacement of the outer edge is kept fixed at the prestretched state so that the solid frame from the experiment is replicated.Fig. 8 presents the deformed configuration of the lens after electric activation by an electric potential difference of magnitude Df % 1:3 kV being applied between the top and bottom surface of the modeled region.For reasons of numerical stability, the potential on the top surface is increased linearly from zero to the maximum value over a small distance.Similarly to the behavior of the lens in the experiment, the electric activation in the simulation leads to a significant change in curvature when comparing the two configurations in Figs.7 and 8.In addition to the described simulation, in which a perfectly circular lens is investigated, a subsequent simulation has been performed considering a rather elliptic shape that may result from production inaccuracies, see Fig. 9.In this regard, the visualization shall help to identify acceptable deviations from the perfect shape and to distinguish them from deviations that would cause a loss of functionality.Besides a deviation from the perfectly circular shape, also different kinds of imperfections could be analyzed, such as changes in material properties, imperfect attachment of the electrodes, or defects within the material itself, as discussed in detail in Section 3.3.

Coupling Tensor
Before introducing the particular model considered, some basics of nonlinear continuum mechanics in the context of electromechanical coupling and mathematical fundamentals are briefly summarized.
We denote the n-dimensional euclidean vector space as V n .Its scalar product is a bilinear mapping of two vectors v 1 ; v 2 to a real number and is denoted as v 1 Á v 2 .We assume that we are given an orthonormal basis fe 1 ; . . .; e n g, i.e. e i Á e j ¼ d ij with the Kronecker delta d ij .Because no distinction between co-and contravariant tensors is required in this work, we define an n-dimensional tensor of order q as multilinear map of q vectors to the real numbers, As a multilinear map, the tensor can also be described by its coefficients with respect to a fixed orthonormal basis of V n , i.e.
. .e i q : Therefore, a zeroth-order tensor can be represented as a scalar, the coefficients of a first-order tensor as a vector of dimension n and the coefficients of a second-order tensor as an n Â n matrix.The coefficients of higher-order tensors can be represented as arrays of order q.For example, the coefficients of a three-dimensional third-order tensor can be represented as a 3 Â 3 Â 3 array.With respect to notation, the single contraction is denoted by a single Á, whereas represents the standard dyadic product.The main tensor of this work is the third-order coupling tensor It represents the sensitivity of the Piola-Kirchhoff-type stress tensor S with respect to the electric field vector E. In other words, it shows how stresses change with changing electric field.This work will use the coupling tensor H defined by  with C a (positive definite) deformation related tensor.The subsequently elaborated electromechanically coupled and nonlinear boundary value problems are solved within the framework of the Finite Element Method.For further details on the finite element implementation the reader is referred to the supplemental material and [23], [24] as well as the references cited therein.

RELATED WORK
Now that the problem has been presented, it can be analyzed from the visualizer's point of view.As mentioned before, tensors appear in many application areas.Thus, there are also many visualization techniques to represent this type of data.Most tensor visualizations are known from diffusion tensor imaging (DTI) and mechanical engineering.Glyphs are widely used to represent tensors.Most glyphs represent symmetric second-order tensors.A popular example is the superquadric glyph by Schultz and Kindlemann [25].It represents such tensors with both positive and negative eigenvalues.Therefore, the tensor norm, the eigenvectors, the rotational invariant, and the uniform scaling invariant part of the tensor are displayed.
However, there are also glyphs that describe asymmetric second-order tensors.Gerrits et al. [26] designed a glyph representing a general two-dimensional, second-order tensor.The glyph preserves the invariance under isometric domain transformations, is scaling invariant, encodes the real eigenvalues and eigenvectors, is unique, and is continuous.A glyph representing three-dimensional tensors was analogously designed.
Glyph-based methods are local methods and, therefore, focus on a detailed representation of the data points.Most glyph-based visualization methods are not used to give an overview of a whole field because they include too much information.Therefore, geometric methods are used, like, for example line-based methods.The best known line-based methods are vector and tensor lines.Tensor lines were introduced by Dickinson [27], [28].To complement this work, Delmarcelle and Hesselink [29] presented hyperstreamlines.These line-based methods are also often used to show the tensor field topology [30], [31].
All the methods mentioned above represent tensor fields of second-order tensors.Some works also exist with respect to higher-order tensor visualization.Most of them deal with applications from the DTI area.For DTI, practically all tensors are symmetric tensors of even-order.For these tensors, Schultz and Kindlemann [32] generalized the tensor ellipsoid and Florack et al. [33] analyzed these tensors by a deviatoric decomposition.A generalization of a line-based method to higher-order, totally symmetric even-order tensors, so-called HOT-lines, were given by Hlawitschka and Scheuermann [34].They also used the deviatoric decomposition and calculated the extreme values on the sphere.Some other works deal with the tensor visualization in mechanical engineering.The most visualized tensor in this application area is the stiffness tensor.There are some glyph approaches in the literature which represent this tensor [35], [36], [37].
To the best of our knowledge, Zobel et al. [38] have so far presented the sole visualization-based work regarding third-order tensors.They have presented a simplification and a visualization of the stress gradient.In one of their visualizations, e.g., two super-positioned distinct-colored ellipsoids present the gradient.One of the ellipsoids represents the average change of all stress vectors, the other one the actual stress tensor.As the considered third-order tensor H in our work represents electromechanical coupling, the work by Zobel et al. [38] is not applicable because their glyph strongly relates to the idea of a spatial stress gradient.
Choosing the right colormap for each visualization is an important point.Thus, there are many works and surveys on this topic.A discussion on colormap guidelines is given by Bujack et al. [39].As the problems discussed in this work include norm variations in different orders of magnitude, the color maps must address this special problem.Nardini et al. [40] designed colormaps for these types of data.They use combinations of different continuous colormaps in a single colormap.

TENSOR DECOMPOSITION
A three-dimensional vector can be easily interpreted in most cases.The interpretation of tensors of second-order, e.g. in terms of matrices representing their coefficients, is quite a challenge.A well known method for the analysis of a second-order tensor is the decomposition into its eigenvalues and -vectors.These eigenvalues and -vectors are then interpreted instead of the tensor itself, which works straightforwardly in the case of symmetric second-order tensors where both, eigenvalues and -vectors are real.Higher-order tensors can (in most cases) not straightforwardly be decomposed by an eigendecomposition (or spectral decomposition).However, another method has been established in order to represent a tensor of arbitrary order up to dimension three.Each of these tensors can be described by a unique set of vectors, so-called multipoles, and scalars.In order to compute these multipoles, the tensor must be decomposed into deviators which can be decomposed into multipoles and scalars.Even though a tensor of arbitrary order in dimension three can be decomposed into these multipoles, this work will focus on three-dimensional third-order tensors.This is done because the EAP application focuses on the understanding of the electromechanical coupling as described by H.The deviatoric decomposition is also known as irreducible decomposition [41].The rotation of these irreducible parts represents the rotation of the overall tensor.The multipoles, on the other hand, clearly represent the deviators as well as the rotation of the deviators.Thus, tensors can be represented by the multipoles uniquely and independent of the local coordinate system.The number of independent coefficients of the deviator increases linearly and not exponentially in terms of order.The analysis of each deviator using the multipoles reduces the amount of values to be considered at the same time, compared to the analysis of the original tensor without loss of information.

Deviatoric Decomposition
Hergl et al. [42] described the deviatoric decomposition and the multipole representation based on this decomposition.Each tensor of any order up to dimension three with or without any index symmetry can be decomposed into traceless and totally symmetric tensors, called deviators.The deviatoric decomposition of a third-order three-dimensional tensor is an orthogonal decomposition, in relation to the inner product, of the 3 3 -dimensional space.
Every tensor of order three can be split into a totally symmetric and an asymmetric part.Through the well known relation between totally symmetric tensors and spherical harmonics, the totally symmetric part can be decomposed into deviators.For the asymmetric part, an isomorphism from the totally symmetric tensors into the asymmetric tensors of order three exists.These totally symmetric tensors can now, through the above named relation, also be decomposed into deviators.This highly non-trivial procedure is described by Backus [43] who also showed that it is a unique decomposition.The deviatoric decomposition of a second-order tensor is common and given by the symmetric part and a vector, representing the antisymmetric part.Auffray [44] presented the decomposition of a general thirdorder tensor H with the left symmetry, i.e.H ijk ¼ H jik , into a third-order deviator D, a second-order deviator D and two first-order deviators d sym and d asym by with the permutation tensor " " " " " " " and the Kronecker delta d ij .Each of the four summands is an irreducible part of the tensor.The first two summands describe the symmetric part, the other two the asymmetric part.An explicit description of all deviator coefficients can be found in the supplementary material.

Multipoles
Higher-order tensors, including third-order tensors, generally cannot be decomposed into eigenvectors and eigenvalues.A generalization of the concept can lead to eigentensors and eigenvalues in some cases, but the eigentensors are of an order higher than one.Even though there are some cases where a mapping can be used to obtain a second-order tensor and to perform an eigendecomposition, a generalization of tensor lines by using this decomposition has so far not been found.Especially a generalization of eigenvalues and eigenvectors for third-order tensors is so far unknown.There is also no other well known third-order tensor decomposition.
Here, we use a finding by Sylvester [45] who proved that each nth-order d-dimensional deviator ID with d > 2 can be represented by the symmetrization sðÁÞ of a, except for a scalar, unique set of d-dimensional vectors m i with cardinality n by (To be precise, the m i are unique up to an even number of sign changes.)These vectors m i are called multipoles and can be interpreted in a physical way.The calculation of the multipoles is described in more detail by Zou et al. [46] and summarized in the supplementary material.
There is a connection between these mathematical multipoles and electrical multipoles.A single charge is called a (electrical) monopole.An (electrical) dipole refers to two opposite charges, where a Coulomb force F C results in a shift u 1 .When taking two (electrical) dipoles and shifting these by u 2 the resulting field is called (electrical) quadropol.This can be done analogously for other (electrical) poles, e.g.octopoles and so on.The (mathematical) multipoles describe exactly this phenomenon.Except for the scalar a from Equation (5), the (mathematical) multipoles describe the shift of the charges.The term multipole generalizes all of these poles.
The multipoles can also be interpreted in terms of material symmetries.The symmetries of the tensor equals the intersection of the symmetries of the deviators.Furthermore, the symmetries of the deviators can be determined from the multipoles.The material represented by the tensor has at a point a symmetry plane, if all material properties fulfill this symmetry.A plane is a symmetry plane of a deviator, if all multipoles fulfill this symmetry.The intersection of the symmetry planes of the different deviators equal the symmetry planes of the tensor.
Further information regarding the group theoretical background of this decomposition is, for example, given by Hergl et al. [42] and in more detail by Hamermesh [47].

VISUAL ANALYSIS
Arrows, lines, and streamlines are well-known visualizations for (steady) first-order tensors.Second-order (symmetric) tensors are commonly visualized by using various glyphs and tensor lines based on the eigenvector (spectral) representation.In contrast, not much research on the visualization of higher-order tensors has so far been made.
The analysis of third-order tensors is often limited to the investigation of the tensor norm.In this section, we propose an approach for the visualization and analysis of third-order tensors in three dimensions (with left index symmetry).For this purpose, we base our visualization technique on a unique set of lines, the multipoles, and a few scalars at every point.
To provide an adequate structural perception, we render the set of lines representing the multipoles as tubes with spherical caps with a user-defined constant width and a constant length.We achieve this by utilizing real-time GPU raycasting of the tube primitives and spheres.In order to reduce clutter and consistently display overlapping multipoles, we split the deviators into separate linked multiple views [48].In this process, we render the first-order symmetric, the first-order asymmetric, the second-order symmetric, and the third-order asymmetric multipoles in four separate views.As a result, comparing deviators in terms of direction is more straightforward.Our approach allows the user to freely assign the glyphs to one or multiple views.This allows the user to recognize subtle differences in terms of the direction of the multipoles.To visualize scalar quantities associated with the deviators, we apply color mapping of the norm to the coloring of the glyphs.The norms of the different irreducible parts have considerable differences in their ranges.Using a linear colormap would lead to a significant loss of information.Thus, we use a color mapping that enables the representation of small value ranges in different orders of magnitude.We design suitable color maps using the CCC-Tool introduced by Nardini et al. [40].
Each value range is mapped using a linear color map.For the different ranges, we used base colors that can be distinguished clearly.The combination of these color maps forms the basis of the color map representing the irreducible part norm.The numbers between the value ranges occurring in the data set are represented by a neutral color.The mapping is linked with all views.This allows a quantification of the scalars and an intuitive comparison with scalars in other views.
The necessary steps for our glyph at each position are to 1) Compute the four deviators.
2) Compute the multipoles and the norm of the corresponding irreducible part for each deviator.3) Render the glyph at the given position in corresponding linked views.We tested the visualization for the four examples introduced in Section 3 and interpreted the results from an engineering perspective as well.Electric loading is increased in time for all examples considered in this work.In the first time step, all multipoles corresponding to all deviators nearly vanish for all three examples, which also reflects that H ¼ 0 for E ¼ 0, cf.Equation (3).With increasing time, respectively increasing electrical loading, the multipoles elongate, and the orientation may undergo small changesthe general orientation, however, remains.

INSIGHTS IN USE CASES
This section presents the visual analysis of the four different electromechanically coupled boundary value problems described in Section 3. Fig. 10 depicts the multipole visualization of the simulation results for the cylinder model introduced in Section 3.1.The plots correspond to the state of maximum loading, i.e.Df ¼ 10 kV.It can be observed that the norm of the irreducible parts corresponding to the second-and third-order deviators is much smaller than those of the two first-order deviators.The multipoles of the first-order deviators are located planar in the x-y plane and are radial symmetric.Regarding the second-order deviator on the other hand, one of the corresponding multipoles is oriented tangentially to the curvature of the cylinder, the other one is orthogonal to the horizontal surface.
This cylinder example is of particular interest since the electric field vector E is always aligned with a principal direction (eigenvector) of the deformation as represented here by the tensor C. For the particular example considered, this direction corresponds to the radial direction in which the two first-order deviators are also oriented.As mentioned in Section 5.2, the multipoles can be interpreted in terms of material symmetries.The symmetries of the material under the influence of the electric field are represented by single sets of multipoles.There are two deviators of first order.Thus, there are two sets of one multipole each.A deviator that is represented by one multipole can have two different types of symmetries.The first one is described by a multipole that vanishes.Then, the deviator is isotropic.If the multipole does not vanish, the deviator has one symmetry plane orthogonal to the multipole and all planes that enclose the multipole are also symmetry planes.In this example, there are two deviators of first order.Thus, there are five different cases of how the material can behave under the influence of the electric field.First, both multipoles corresponding to the deviators of first-order vanish.Then, the symmetries are given by the deviators of order two and three.Second, one multipole vanishes.Then, the symmetries are defined by the non-vanishing multipole.The third case is given if the multipoles equal each other.Consequently, the symmetries of the intersection of the deviators are the same and represent the material's symmetry.The fourth case is given if the multipoles are orthogonal.The material can have three orthogonal symmetry planes at most.When both multipoles are different and not orthogonal, as in the fifth case, the material does not behave symmetrically to any plane.An illustration of these five cases can be found in the supplemental material.In this example, the two multipoles corresponding to the deviators of first-order equal each other so that the symmetries of the material under the influence of the electric field are given by a subset of the plane orthogonal to the radial symmetric direction and all planes that enclose the radial symmetric direction.The deviators of order two and three restrict the symmetries to three orthogonal symmetry planes that are given by the x À y-plane, the plane orthogonal to the radial direction and the plane orthogonal to these two planes.
It is also interesting to note that the plane strain conditions considered together with strain in circumferential direction remaining comparatively small, yield the state of deformation in the plane perpendicular to the electrical field direction to be almost isotropic.Conceptually speaking, the principal strains, respectively, stretches in the plane perpendicular to E are almost identical.This yields the norm of the second-order deviator of H to be comparatively small, also in comparison to the norm of the third-order deviator of H.Moreover, in the case of a perfectly isotropic deformation state in that plane, this norm would vanish identically.Such scenarios correspond to extremal states of energy -in other words, given the level of electric loading, the orientation of the electric field results in an extremal state of deformation and stiffness.Such properties are most attractive from an optimal design perspective of smart EAP based devices.Similar problems are established within the design of, e.g., fibre-reinforced composites, where the fibre orientations remain constant.In view of the electromechanical coupling considered in the present work, the electric field can be compared with the (transversely isotropic) fibre contribution -the main difference is that both referential orientation and magnitude of the electric field generally change during the loading process.Since the electric field vector is always aligned with a principal axis (eigenvector) of the deformation for the cylinder model -in fact E is aligned with the eigenvector related to the minimum eigenvalue (compression) of C -referential stresses and deformation measures commute, i.e.C Á S ¼ S Á C. Such states are well-established so as to correspond to extremal states of energy within the theory of anisotropic hyper-elasticity.A key contribution of the visualization is that vanishing contributions of the second-order deviator of H indicate such extremal states of energy which, from an engineering design point of view, are most interesting and relevant for an improved design.Moreover, the norm of the irreducible part corresponding to the first-order deviators shown in Fig. 10 both increase from the outer radius r o to the inner radius r i -in other words, electromechanical coupling increases from r o to r i .From an engineering perspective, such effect is intuitively expected (and highlighted here by the visualization) since the level of deformation also increases from r o to r i .In summary, the particular example of the cylinder exhibits properties such as i) radially oriented electric field, ii) coaxiality of stresses and strains and iii) increasing electromechanical coupling from r o to r i , which are intuitively expected on the one hand and directly reflected by the visualization established in this work on the other.In this regard, the example may be considered as a proof of concept based on which more advanced examples can be investigated where interpretations without visualization are not that intuitive.
The second example of a bending beam, as described in Section 3.2, is more complex but still allows for some intuitive interpretation.The visualization results are shown in Fig. 11.Specifically speaking, a detailed view of the supported end of the beam is presented.For the purpose of illustration, the beam is shown in its undeformed shape.However, all quantities illustrated correspond to the state of maximum loading, i.e.Df ¼ 2 kV.Since only the upper half of the beam (in y-direction) is electrically activated, the lower part shows no electromechanical coupling and all multipoles related to the deviators of H vanish identically.Regarding the upper half of the beam, the norm of the irreducible part corresponding to the first-order deviators is approximately 10 2 times larger than the respective norm corresponding to the third-and second-order deviators.
The multipoles of the first-order deviators are aligned with the y-direction.All multipoles of the third-order deviators lie in the y-z-plane and, in particular, one of those is oriented in y-direction as well.The other two multipoles have the direction of the first one as a bisecting line.In contrast to these characteristics, the multipoles belonging to the second-order deviator are not oriented uniformly, which we attribute to noises resulting from the very small norm of the irreducible part of this deviator.
These orientations of the multipoles again allow an interpretation of the material properties under the influence of the electric field in terms of symmetries.Taking the norm of the irreducible part into consideration, the multipoles corresponding to the deviator of order two, do not affect the symmetries much.Thus, the bending beam behaves under the influence of the electric field which is nearly symmetric to the three planes that are in this undeformed representation given by the coordinate planes.
From an engineering point of view, the chosen electrical loading of the beam induces a mode similar to pure bending (in addition to the compressive deformation in the thickness direction of the upper half of the beam).Such a state of deformation does not include classic beam shear contributions.By analogy with the previous example of the cylinder, this is reflected by the fact that the two first-order deviators are aligned with the direction of the electric field E which, moreover, is a principal direction (eigenvector) of the deformation measure C. Within the finite deformation setting considered, the beam undergoes curvature in both longitudinal (z) direction and transverse (x) direction.One of the Bernoulli hypotheses, namely that cross-sections of the beam remain plane during the deformation, is typically violated -in particular for such finite deformation setting.The second-order deviator of H yields non-vanishing contributions in the upper half of the beam, which indicates that the deformation state is not isotropic (equi-biaxial) in the plane perpendicular to the direction of the electric field.Moreover, it can clearly be seen in Fig. 11 that the norm of the irreducible part related to the deviators increases from the mid-layer (neutral axis) of the beam towards the top surface.While the longitudinal strain is (quasi) zero along the neutral axis, it reaches a maximum at the top surface (tension state), which is clearly sensed by the third-order deviator of H.
The cube, described in Section 3.3, is a more complex example due to the included spherical defects.The general response is, from an engineering point of view, less intuitive compared to the previous examples.Considering the cube as a homogeneous model without the presence of any defects, a compression in y-direction due to an electric field in the same direction would result in the first-order deviators of H being oriented in y-direction as well.This is illustrated in Fig. 12 when focusing on the regions off the spherical hole.In this part of the cube, the response of the corresponding homogeneous setting with the first-order deviators of H oriented in y-direction is indicated, whereas the orientation of the deviators changes in the inhomogeneously deforming region near the spherical defect.Regarding the third-order deviator of H, two of the corresponding multipoles have the first-order multipoles as a bisecting line which, moreover, also corresponds to the orientation of the third multipole of the third-order deviator.
The norm of the irreducible part related to the first-order deviators is approximately 10 15 times larger than this related to the third-and second-order deviators.In direct neighborhood to the hole, the response becomes rather complex and, from an engineering perspective, the degree of deformation inhomogeneity significantly increases.The multipoles corresponding to the two first-order deviators are oriented around the defect -in other words, near the surface of the hole, the first-order deviators are parallel to the surface.This clearly indicates the (general) orientation of the electric field (by analogy with, e.g., a fluid flow around an obstacle).Nonvanishing second-order deviators may indicate, as mentioned in the previous examples, deviation from extremal energy states with isotropic state of deformation in the plane perpendicular to the electric field.From an engineering point of view -the area near to the surface of the hole and along the z-axis is practically unloaded, whereas the area nearby the surface of the hole and along the y-axis undergoes a state of maximum compression.It is interesting to note that this is clearly reflected by the norm of the irreducible parts, respectively the color of the multipoles: small values are present in regions of low electromechanical loading, whereas larger values are present in regions of higher electromechanical loading.From an engineering perspective, the visualization shows that the spherical hole i) lowers electrical and mechanical loading in some regions and increases electrical and mechanical loading in other regions, ii) lowers electrical coupling in some regions and increases electromechanical coupling in other regions, iii) influences the orientation of the electric field in a streamline manner, and iv) induces inhomogeneous states of deformation.In view of the design of smart EAP based devices, such information is most valuable since the level of electromechanical coupling influences the efficiency of the smart device on the one hand, and maximum electrical and mechanical loading levels influence lifetime as well as its failure properties on the other hand.
The lens datasets described in Section 3.4 can be used to analyze slight deviations of the geometrically perfect model.Such simulations can be used to define acceptable deviations of the perfect structure and deviations that would affect the functionality.
The glyph representation of the perfect lens is pictured in Fig. 13a and the imperfect lens, the one with the greater deviation, in Fig. 13b.More figures of the two datasets can be found in the supplemental material.Figs.13a and 13b both use 1500 seed points.Using a laptop with an AMD Ryzen 5 3500u with radeon vega mobile gfx Â 8 it took 30 seconds to calculate and visualize all four multipole glyphs of the coupling tensor.
To provide a visual spatial division, black lines have been added to divide the lens into four equal areas.Analyzing the norm of the irreducible parts shows that the influence of the first-order deviators is much larger than the influence of the second-and the third-order deviator for the perfect as well as for the imperfect lens.For both models, all illustrations of the multipoles show a ring around the inner part of the lens, where the orientation of the multipoles and the norm of the irreducible parts, as well as the norm of the tensor itself, change considerably.This can be attributed to the simulation setting.The electric potential has been prescribed to be linearly increasing in this region, from a value of zero up to the maximum applied value acting on the rest of the flat lens surface.By doing so, numerical instabilities, which might follow from a sharp jump in the electric potential, are avoided.Thus, the change of the multipole orientations is a side effect of the simulation and will not be analyzed any further.
The norm of the irreducible parts describes the influence of the corresponding deviator.For the perfect lens, the irreducible part norms are radial symmetric.This is not the case for the geometrically imperfect lens.For this example, the norms of the irreducible parts change by a rotation around the y-axis.The norm is highest around the x-axis and decreases by a rotation up to the lowest values around the z-axis.
As described before, the multipoles can be analyzed in terms of symmetries.Thus, the perfect lens is symmetric to the x À z-plane.It can be seen in Fig. 13a that the multipole corresponding to the first-order deviator describing the  symmetric part aligns with the y-axis.The same applies to the multipole corresponding to the first-order deviator describing the asymmetric part.Due to the norm of the irreducible part corresponding to the second-order deviator, the influence of this deviator is much smaller than the influence of the two first-order deviators.The multipoles of the firstorder deviator describing the symmetric part have the same orientation as those of the first-order deviator describing the asymmetric part.Therefore, taking into account that the norms of the irreducible parts corresponding to the deviators of order two and three is much smaller, it can be assumed that the material state under these conditions is nearly symmetric with respect to a rotation around the z-axis.
From an engineering point of view, the comparison of the visualization results corresponding to the geometrically perfect and imperfect lens model shows some clear differences which are relevant regarding improved design.For the perfect lens, the irreducible part norms corresponding to the symmetric part of the first-order deviator show a perfect radial symmetry in the flat lens region under consideration, cf.Fig. 13a.Since the geometrically imperfect lens does not yield such a uniform distribution, the irreducible part norm of the symmetric first-order deviator could be used as an indicator for the quality of the shape of the electromechanical device.Its asymmetric counterpart, on the other hand, shows rather similar results with regard to both deviator norms and orientations of the multipoles when comparing the perfect and imperfect lens model.Thus, the geometric imperfection effects mainly the symmetric part.From the visualizations of the higher-order deviators, further conclusions can be drawn.For the second-and third-order deviators, cf.Fig. 13a, belonging to the perfect lens model, we observe a symmetry line at an angle of 45 between the x-and z-axis.In the case of a geometric imperfection, cf.Fig. 13b, this symmetry line is shifted in the circumferential direction toward the x-axis or, in other words, towards the side where the radius of the curved lens region increases.From this knowledge gained from the visualizations, the design of electromechanical devices could, in general, be improved with regard to energy efficiency by aiming for a uniform or at least perfectly symmetric distribution of the respective tensorial quantities.For the particular example of the bioinspired tunable lens, the visualizations could additionally be an indicator of the functionality of the device since production errors such as indicated by the imperfect lens, geometry may lead to an inaccurately adjusted curvature of the lens.

LESSONS FOR VISUALIZATION RESEARCH
Overall, the visual analysis of the electromechanical coupling represented by the third-order tensor H was a success.This can especially be seen from the analysis of the third and the fourth example, the cube with a defect and the two lenses.The engineers wrote: "such information is most valuable since the level of electromechanical coupling influences the efficiency of the smart device on the one hand, and maximum electrical and mechanical loading levels influence its lifetime as well as its failure properties" as their evaluation in the last section where, based on our visualization, more insight could be gained.Looking back at the original questions from Section 2 posed by the engineers, it can be seen that we were able to I. describe the behavior of the material by visualizing information on the electrical loading and the directional information, including information on the electrical field and the deformation.II.present a tensor decomposition which is so far not well known, respectively established in the community.III.reduce the third-order tensor of the presented examples to a few scalars and vectors.In the presented examples the normalized vectors regarding the two first-order deviators and one of the vectors regarding the third-order deviator equal each other, which reduces the number of independent vectors.The visualization of the multipoles makes it possible to gain information on the electromechanical coupling.This allows the analysis of the third-order tensor H in a way that was not possible before and which creates new insights for engineers.IV. apply established knowledge to obtain a first idea for an interpretation of the multipole representation and use the gained knowledge to obtain new information.This shows how knowledge from simple examples can be transferred to more difficult ones.In the end, the engineers were able to obtain quite a good understanding of the material behavior, even in the more difficult and complex case of a cube with a spherical defect and the lenses.The information gained from the lens example can be used in the future to improve the design.V. show information that may support and improve component design in the future.Even if the visualization has not yet immediately improved the design of EAP-based actuators and sensors, the engineers indicate that, in view of the design of smart EAPbased devices, such information is most valuable.Thus, the gained information concerning the coupling and the electrical field can be used to support the design of such components.Even though this work may enable just a small step towards a complete interpretation of the third-order tensor H, the visualization helps to answer the questions asked by the engineers and can hopefully be used to improve the design of sensors and actuators in the future.Therefore, we have found answers to all questions from the engineers, even if there is more to be obtained in the future.
There are three main lessons for visualization researchers which can be drawn from this contribution.
1) The first lesson is that the design of electroactive polymer based smart devices and the analysis of different coupling problems is a most challenging and highly interesting interdisciplinary research field.The missing (visual) analysis methods complicate a systematic design and the optimization of related smart systems.2) At the heart of the design, there is the electromechanical coupling described by a third-order tensor.Thus, the analysis of third-order tensors turns out to be a second task that the visualization community should concentrate on.3) A third lesson is based on the well-known fact that large and complex data requires a good overview.In visualization literature, a wide array of tools are at our disposal.By utilizing GPU-based rendering, applying linked multiple views, and using the recently introduced CCC-Tools, we were able to provide a well-structured overview of the data, taking into account different scales of ranges.Later in the research process, the visualization design may be reworked and specialized.Nevertheless, the straightforward visualization using these tools was sufficient for this basic applied research and contained a vast amount of information.

CONCLUSION AND FUTURE WORK
In this contribution, we present a successful visual analysis of a simulation of an electromechanically coupled problem.The simulation results in tensor field data sets of different (especially third) order, which are difficult to study due to their complexity and the missing analysis methods of thirdorder tensors.We demonstrate that the multipole representation allows the reduction of the third-order tensors of these examples to a set of two to five vectors and two to three scalars and allows a first interpretation of these thirdorder tensors, which have so far only been analyzed by their norm.As can be seen from the insights gained by the engineers, our visualization allows us to answer the original questions posed by the engineers.Our tensor visualization provides illustrative information on the coupling response of the material -here electromechanical coupling -with information on both the level of coupling and the directional coupling properties.So far, there is only very little information regarding this tensor in the literature.In contrast, the described method allows the interpretation of its directional properties and absolute values, which is most valuable for optimal design approaches and lifetime, respectively failure predictions.It can be seen in the examples that nearly identical tensors are difficult to tell apart without color coding.In the future, the visualization of tensors that are similar should be analyzed.The deviatoric decomposition is a powerful algebraic tool that may be helpful for other higher-order tensor visualization problems as well.To give some indication, we mentioned at the beginning of this paper that there are other coupled problems apart from the electromechanical coupling studied here, for example, electromagnetic couplings or diffusion problems.Analyzing these problems raises the same questions regarding higher-order tensors.Thus, the presented visualization may also help domain experts in these applications.Even if the problems may differ in the kind of coupling, the interpretation developed in this paper may be used for these settings as well.Furthermore, the visualization can be used for each tensor up to dimension three with any index symmetry, which could also be analyzed in the future.

Fig. 1 .
Fig. 1.Cylinder model: The initial mesh containing n el ¼ 105 finite elements is shown including the particular boundary conditions.The colored region shows the deformed configuration and indicates the electric potential f.

Fig. 2 .
Fig.2.Working principle of a dielectric elastomer actuator (DEA).The second picture shows the actuator at rest while the activated state is presented in the first and third picture.Reproduced from[20] under the terms of the Creative Commons Attribution License (CC BY).

Fig. 3 .
Fig. 3. Working principle and fluorescence image of a jellyfish-inspired and DEA-based robot.Reproduced from [21] under the terms of the Creative Commons Attribution License (CC BY).

Fig. 4 .
Fig. 4. Beam model: The initial mesh containing n el ¼ 1472 finite elements is shown including the particular boundary conditions.The colored region shows the deformed configuration and indicates the electric potential f.

Fig. 5 .
Fig. 5. Cube model: The initial mesh containing n el ¼ 783 finite elements is shown including the particular boundary conditions.The colored region shows the deformed configuration and indicates the electric potential f.

Fig. 7 .Fig. 8 .
Fig. 7. Lens model: The deformed mesh containing n el ¼ 5628 finite elements is shown after a radial effective prestretch of % 10% has been applied.

Fig. 10 .
Fig. 10.Cylinder model: Multipole representation of the third-order tensor H using our visualization.Two different views of the following deviators are given: first-order deviator d sym representing the symmetric part, first-order deviator d asym representing the antisymmetric part, second-order D , thirdorder deviator D. The arrangement of the deviator fields of the side view (y À z plane) corresponds with these of the top-view (x À y plane).The color map represents the norm of the respective irreducible part.

Fig. 11 .
Fig. 11.Beam model: Multipole representation of the third-order tensor H using our visualization.Top left -first-order deviator d sym representing the symmetric part, Bottom left -first-order deviator d asym representing the antisymmetric part, Top right -second-order D , Bottom right -third-order deviator D. The color map represents the norm of the respective irreducible part.If the norm is smaller than 1 Á 10 À16 , we display white spheres.

Fig. 12 .
Fig. 12. Cube model: Multipole representation of the third-order tensor H using our visualization.Top left -first-order deviator d sym representing the symmetric part, Bottom left -first-order deviator d asym representing the antisymmetric part, Top right -second-order D , Bottom right -third-order deviator D. The color map represents the norm of the respective irreducible part.

Fig. 13 .
Fig. 13.Lens model: Representation of the multipole glyphs for the different deviators.The color map represents the norm of the respective irreducible part.The black lines highlight the angles 1=8p , 1=4p, and 3=8p.

Chiara
Hergl received the BSc degree in computer science, in 2017, and the master degree (diplom) in mathematics from the Leipzig University.She is currently working toward the PhD degree with the Leipzig University, Institute of Computer Science, Image and Signal Processing Group.Her research interests includes the analysis and visualization of higher-order tensors.CarinaWitt received the BSc degree in mechanical engineering and the MSc degree from TU Dortmund University, in 2018 and 2019, respectively.In 2020, she started working as a research assistant with the Institute of Mechanics, TU Dortmund University.Her research interest includes material modelling of fiber-reinforced composites as well as electromechanically coupled problems.Baldwin Nsonga received the BSc degree from the TU Bergakademie Freiberg, in 2013, and the MSc degree from the Leipzig University, in 2016.He is currently working towards the doctorate degree with the Leipzig University, Institute of Computer Science, Image and Signal Processing Group.His research interest includes flow visualization.Andreas Menzel received the diploma degree in civil engineering from Leibniz University Hanover, in 1997.The same year, he moved to the Department of Mechanical and Process Engineering, Technical University Kaiserslautern., and the Dr.-Ing.degree in 2002, he continued as a postdoc and was awarded the Habilitation for Mechanics in 2006.The following year, he held a temporary professorship with the University of Siegen.He joined the faculty with the Department of Mechanical Engineering, TU Dortmund University, in 2007 and holds a double affiliation with the Division of Solid Mechanics with Lund University.His main research interests include the fields of computational mechanics, material modelling, and smart materials.Gerik Scheuermann received the master's degree (diplom) in mathematics and the PhD degree in computer science from the Technical University of Kaiserslautern, in 1995 and 1999, respectively.He is a full professor with the Leipzig University since 2004.He has coauthored more than 250 reviewed book chapters, journal articles, and conference papers.His current research interests include focus on visualization and visual analytics, especially on feature and topology-based methods, flow visualization, tensor visualization, environmental visualization, document visualization, and visualization for life sciences.He has served as paper cochair for Eurovis 2008, IEEE SciVis 2011, IEEE SciVis 2012, and IEEE PacificVis 2015.He has co-organized TopoInVis 2007, AGACSE 2008, EuroVis 2013, and IEEE VIS 2018, as well as three Dagstuhl Seminars on Visualization.