Symmetry Classes and Matrix Representations of the 2D Flexoelectric Law

We determine the different symmetry classes of bi-dimensional flexoelectric tensors. Using the harmonic decomposition method, we show that there are six symmetry classes. We also provide the matrix representations of the flexoelectric tensor and of the complete flexoelectric law, for each symmetry class.


Introduction
In recent years, the study of electromechanical effects has attracted the attention of many researchers. The best known of these electromechanical effects are piezoelectricity [1][2][3] and flexoelectricity [4,5]. Piezoelectricity is the appearance of a polarization in a dielectric material when it is subjected to an uniform mechanical strain, while flexoelectricity is the linear response of a polarization to a non-uniform strain or a strain-gradient. Flexoelectricity is a non classical electromechanical coupling and it can be viewed as a higher-order electromechanical effect with respect to piezoelectricity [5].
This higher-order electromechanical coupling has attracted a lot of interest in recent years as it allows electromechanical couplings in materials for which the piezoelectric tensor is zero [6,7]. For instance, this flexoelectricity coupling has been identified in polymers [8], biological membranes [9], and, recently its role in the bone remodeling process has been evidenced [10,11]. The design of mesostructure producing this effect is a new trend in topological optimization [12] since the mastery of this effect makes it possible to consider the development of new sensors, energy harvesting devices, and even the development of polymer actuators. This last application arouses a strong interest in the field of soft robotics, where the optimization of this coupling would allow the control of the geometry of polymers via an electric field [13]. Recently, the flexoelectricity phenomenon is highlighted in 2D materials like 2D crystalline membranes or graphene sheets, and also in thin films [14,15].
Our work is part of a project in line with this perspective of optimal design of higher-order electromechanical mesostructures. Following the approach opened in elasticity to work with tensor invariants [16], we aim at developing an invariant-based topological optimization algorithm. But before considering algorithmic questions related to topological optimization, the first step is to fully characterize the different tensors that contribute to the considered constitutive law. The present contribution will solely concern this aspect of the project. Furthermore, this task has to be undertaken by considering the final application, hence by trying to give a physical content to invariants.
Piezoelectric and flexoelectric effects can be introduced via the following constitutive equation [1,5]: where p is the polarization, q the electric field, ε ∼ the strain tensor, and η the strain-gradient tensor.
The second-order tensor S ∼ , the third-order tensor P , and the fourth-order tensor F ≈ describe the dielectric effect, the direct piezoelectric and flexoelectric effects [5], respectively. It is important to note that in centrosymmetric materials, the third-order tensor P vanishes. This means that piezoelectricity exists only in non-centrosymmetric materials. The vector space of flexoelectric tensors is denoted by Flex. For the optimal design of flexoelectric structures, the understanding of the algebraic structure of Flex is mandatory. Important issues concerning Flex include the following questions: • What types of anisotropies can flexoelectric tensors model? • How does one characterize flexoelectric materials independently of their spatial orientation? • How does one parametrize the different characteristics (anisotropy, coupling, etc.) of flexoelectric materials in an intrinsic way?
In mathematical terms, these questions consist of studying the linear representation of the d-dimensional orthogonal group O(d) on Flex (d = 2 in 2D and d = 3 in 3D) and characterizing the geometry of the orbit space Flex/O(d). Answering these questions in the 3D setting is a rather difficult task. For instance, for the well-known situation of the elasticity tensor, the invariant basis comprises 294 invariants in 3D, for only 5 invariants in 2D. Hence, to develop a mesostructural optimization algorithm, we aim to first formulate the problem in a 2D setting. In addition to the desire to reduce the complexity of the problem, the 2D case has an intrinsic interest for the study of 2D materials or structures with flexoelectric effects.
The objective of the present study is more modest since only the 2D situation is considered. Furthermore, only the first of the three questions will be considered in this article. Our goal is to set the framework to answer two other questions in future work.
The problem of classifying flexoelectric materials by determining the number and types of symmetry classes of flexoelectric tensors has been solved in 3D [17]. In this work, we will solve the same problem in the two-dimensional context by using a completely different approach than the one used in [17]. Compared to the article [17], and in order to move towards the study of invariant algebras, the present contribution states in a more detailed way the structure of the space Flex. Furthermore, we introduce a harmonic decomposition obtained by the Clebsch-Gordan algorithm, which differs from the one of [17] obtained by the Spencer algorithm [18]. Another significant difference is that our decomposition is orthogonal and based on the physics of strain-gradient elasticity while the one of [17] has no clear physical interpretation. The physical content of our decomposition is detailed in the present contribution. As such, the present paper is not a reduced version of the reference [17] but a complete study on its own.
From the above presentation, it can be seen that the flexoelectricity theory is closely related to the strain-gradient elasticity theory [19,20]. Some studies [21,22] provide the complete description of bi-dimensional anisotropic strain-gradient elasticity by determining the different symmetry classes of the constitutive tensors and their associated matrix forms. We will complete these studies in the case where strain-gradient effects are coupled with flexoelectric ones.
Unlike elasticity tensors, it is not always possible to use for the flexoelectric tensors the Voigt notation as pointed out by Yudin et al. in [5]. Therefore, to provide the matrix representation of the complete flexoelectric law we need to construct appropriate orthonormal bases in order to express the different constitutive tensors in a matrix form.
This article is organized as follows. In Section 2, we introduce the notation used throughout the text. Section 3 is devoted to the description of the flexoelectricity law. We present the constitutive equations for a flexoelectric material. Then we describe how we can decompose the space of flexoelectric tensors into a direct sum of irreducible spaces. Finally, we determine the different symmetry classes of the flexoelectric tensor and of the complete flexoelectricity law. In Section 4, we provide the matrix representations of the flexoelectric tensor and of the complete flexoelectricity law. In Section 5, we give the explicit harmonic decomposition of the flexoelectric tensor. We rewrite the matrix representations of the flexoelectric tensor obtained in Section 4 in terms of the components of the harmonic decomposition and discuss the degenerated situations. We wish to emphasize that the explicit harmonic decomposition provided in Section 5 is a new result of practical importance for the computation of physically-based invariants of the flexoelectricity tensor.

Notation
In this section we provide some notation and conventions that will be used throughout the article. Let R m be the m-dimensional real Euclidean space. We equip this space with a rectangular Cartesian coordinate system with origin O and an orthonormal basis P = {e 1 , . . . , e m }. We associate with each point in R m a unique m-tuple of real numbers (x 1 , . . . , x m ). We can then view R m as a m-dimensional vector space with elements x = (x 1 , . . . , x m ) = x 1 e 1 + · · · + x m e m = x i e i , where Einstein summation convention is used, i.e., when an index appears twice in an expression, it implies the summation of that term over all the values of the index.
Below are provided some specific notations and conventions used in this article. Groups: The following matrix groups are considered in this article: • GL(2) the group of invertible 2 × 2 matrices over R; • O(2) the group of 2 × 2 matrices Q satisfying Q ∈ GL(2) and Q −1 = Q T , where Q −1 and Q T stand for the inverse and the transpose of the matrix Q. This group is called the orthogonal group; • SO(2) the subgroup of O(2) of matrices Q of determinant 1, called the special orthogonal group.
As a matrix group, O(2) is generated by: where R(θ) is a rotation by an angle θ and P(n) is the reflection across the line normal to n. SO (2) corresponds to the group of plane rotations generated by R(θ). The following finite subgroups of O(2) will be used: • 1 the identity group; • Z k (k ≥ 2) the cyclic group with k elements, generated by R(2π/k). For k = 1, we have Z k = 1; • D n k (k ≥ 2) the dihedral group with 2k elements, generated by R(2π/k) and P(n). For k = 1, the group D n 1 will be denoted by Z n,π 2 .
It has to be noted that, up to conjugacy by an element of SO(2), any D n k is conjugate to D e 2 k , will simply be denoted by D k .
Subgroups of O(2) can be divided into four subsets according to the nature of their generators. A group will be said to be: • Centrosymmetric (denoted by I) if its contains the inversion R(π), and non-centrosymmetric (I) otherwise; • Chiral (denoted by C) if its does not contain reflection or mirror π := P(n), and achiral (C) otherwise.
The different situations are reported in Table 1.
Tensor products: • ⊗ denotes the usual tensor product of two tensors or vector spaces; • ⊗ n denotes the n-th power of the tensor product, e.g., ⊗ n V := V ⊗ · · · ⊗ V (n copies); • ⊗ s denotes the symmetrized tensor product; • ⊗ indicates the twisted tensor product defined by: Tensor spaces: • T n := ⊗ n R 2 is the space of n-th order tensors with no index symmetry on R 2 ; • S n := S n (R 2 ) is the space of completely symmetric (by completely symmetric we mean symmetric with respect to all permutations of indices) n-th order tensors on R 2 ; • K n is the space of n-th order harmonic tensors (i.e., completely symmetric and traceless tensors), with dim(K n ) = 2 if n > 0, 1 if n = 0, −1; • K 0 is the space of scalars; • K −1 is the space of pseudo-scalars (i.e. scalars which change sign under improper transformations).
We will use tensors which are of different orders. Tensors of order −1, 0, 1, 2, 3, 4, 5, and 6 are denoted by β, α, v, and a contractions are written ·, . ., . . ., and . . . ., respectively. In components with respect to P, for general tensors A and B, these notations correspond to: When needed, index symmetries of both spaces and their elements are expressed as follows: (..) indicates invariance under permutations of the indices in parentheses and .. .. indicates symmetry with respect to permutations of the underlined blocks. For example, a ∼ ∈ T (ij) means that a ij = a ji .

Special tensors:
• 1 ∼ is the second-order identity tensor with components δ ij , where δ ij is the Kronecker delta; is the fourth-order identity tensor on S 2 ; • ∼ denotes the 2D Levi-Civita tensor defined by: Matrix spaces: • M n is the space of n × n dimensional square matrices; • M s n is the space of n × n dimensional symmetric square matrices; • M p,q is the space of p × q rectangular matrices.
Miscellaneous notation: denotes hereafter an isomorphism.

Flexoelectricity Law
In this section we describe the flexoelectric law in the static regime. In [6,23] the flexoelectric effect was introduced via an energy density function. Here, we introduce this effect via constitutive equations.
For a linear flexoelectric material, the constitutive equations define the stress tensor σ ∼ , the hyperstress tensor τ , and the polarization p as linear functions of the strain tensor ε ∼ , the strain-gradient tensor η , and the electric field q. These constitutive equations are given by: is the dielectric susceptibility tensor; • P ∈ Piez := T ∈ ⊗ 3 R 2 | T ∈ T (ij)k is the piezoelectric tensor; is the transpose of the piezoelectric tensor P .
The transposition is defined as (P ) ijk = P jki ; The transposition is defined as (F The first two right-hand side terms in Equation (5) describe the direct piezoelectric and flexoelectric effects (we use the convention to represent the direct piezoelectric and flexoelectric effects by the transposed tensors P and F ≈ , respectively). Expressing constitutive law, we distinguish two types of tensors: State tensors describe point-wisely the different physical fields (primal and dual) considered in the model. Constitutive tensors model the influence of the matter on these state tensor fields, more precisely they describe how primal and dual fields are connected by the matter. The behaviour of a linear flexoelectric material is described by the constitutive tensors. In what follows, we will denote by F the sextuplet S ∼ , P , F ≈ , C ≈ , M , A ∼ ∼ ∼ and by F lex the complete space of the flexoelectricity law, i.e., Our goal is to classify flexoelectric materials with respect to their spatial symmetry properties. To do this, we need to introduce the notion of symmetry class.

Notions of Symmetry Group and Symmetry Class
We give here some abstract definitions which contain all the information we need to classify flexoelectric materials. We begin by defining the action of the orthogonal group O(2) on the space T n . The group O(2) acts on the space T n by the Rayleigh product defined, with respect to a basis, by Here and throughout the rest of the article, we use the convention to treat elements of O(2) as orthogonal second-order tensors, this justifies the notation Q ∼ .
We will say that a given tensor T ∈ T n is fixed under an orthogonal transformation Q We define the symmetry group of T, denoted G T , as the set of all orthogonal transformations that leave T fixed: It is known [24,25] that each symmetry group is a closed subgroup of O(2). We say that two subgroups is conjugate to exactly one element of the following collection: We will say that two tensors T 1 and T 2 are symmetrically equivalent (we write T 1 ∼ T 2 ) if their symmetry groups G T 1 and G T 2 are conjugate. We define the symmetry class of a tensor T as the set [G T ] of all subgroups of O(2) conjugate to G T , i.e., In other words, two tensors are symmetrically equivalent if and only if their symmetry groups belong to the same symmetry class. One can show [26] that the relation ∼ is an equivalence relation on T n . An equivalence class for this equivalence relation is called a stratum [27]. For each element H ∈ S, we distinguish two types of strata: The open stratum denoted Σ [H] which is the stratum of tensors whose symmetry class is exactly [H] and the closed stratum denoted Σ [H] which is the stratum of tensors whose symmetry class is at least [H] (in the set of all symmetry classes we introduce a partial ordering relation as follows:

Symmetry Classes of the Complete Flexoelectric Law
We now apply the previous notions to the complete space F lex defined in Equation (6). The action of O(2) on the space F lex is given by: where is the action defined in Equation (7). The symmetry group of F is defined as the intersection of the symmetry groups of the constitutive tensors, i.e., and its symmetry class is defined as: In 2D, the number and types of symmetry classes of each tensor space in Equation (6) are present in the literature except for Flex up to authors best knowledge. The calculations of the types of symmetry classes and their number are based on the technique of harmonic decomposition [26,28]. The theorems leading to the results summed-up in the following tables can be found in [28]. In the following tables, #indep(T) indicates the number of independent components of the tensor T in each symmetry class and the in-parenthesis number indicates the minimal number of components in an appropriate basis. It is important to note that for tensors belonging to a stratum of type Σ [Z k ] , and conversely to elements of Σ [D k ] , bases adapted to symmetry elements are determined up to a rotation angle. This free rotation parameter can hence be used to decrease the number of minimal parameters in these optimized bases. In general, this rotation is not uniquely defined.

•
Con: • Piez [29,30]: Name Digonal Orthotropic Tetrachiral Tetragonal Hexachiral Hexagonal Hemitropic Isotropic #indep(A ∼ ∼ ∼ ) 21 12 9(8) 6 7(6) 5 5 4 To complete the picture, it remains to determine the symmetry classes of flexoelectric tensors. The tool we use to do this is the isotypic decomposition, which consists in decomposing a finite-dimensional vector (2) and T ∈ K) and its only invariant subspaces are itself and the null space. It is known that each O(2)-irreducible space is isomorphic to a direct sum of spaces of harmonic tensors K n [28,32]. Such a decomposition is interesting since the O(2)-action on K n is elementary and given by ρ n [33], with for n ≥ 1: The O(2)-action on K 0 is the identity and the O(2)-action on K −1 is given by the determinant of the transformation.
We will now determine the isotypic decomposition of the space Flex.

Isotypic Decomposition of the Space of Flexoelectric Tensors
Due to the linear relation τ = F ≈ · q, we can identify F ≈ with a linear map from T l into T (ij)k . Therefore, we have Flex T (ij)k ⊗ T l . To determine the isotypic decomposition of Flex, we use the Clebsch-Gordan formula which allows us to decompose a tensor product of two irreducible spaces into a direct sum of irreducible spaces. For more details, we refer to [28]. We remark that T l is nothing but K 1 . It then remains to determine the isotypic decomposition of the space T (ij)k . For this, we use the following result, the proof of which is found in [28].

Lemma 1.
For every integers p > 0 and q > 0, we have the following isotypic decompositions, where the meaningless products are indicated by ×: The space Flex admits the following isotypic decomposition: where the notation m k K k := m k l=1 K k is used (m k is called the multiplicity of K k ).
Proof. We have Flex T (ij)k ⊗ T l . Let us decompose the space T (ij)k : Using again Lemma 1, we deduce that: We obtain the decomposition of Equation (9) by grouping spaces of the same order.
To determine the symmetry classes of the flexoelectric tensor, we will use the following result.
Lemma 2 (Proposition 4.10 of [28]). Let T 2p (p > 0) be a space of even-order tensors. If T 2p contains the irreducible space K −1 in its isotypic decomposition, then where I(T 2p ) denotes the set of symmetry classes of the elements of T 2p .
We remark that the case of the space Flex is an immediate corollary of this result. Indeed, Flex is a space of fourth-order tensors containing the irreducible space K −1 in its isotypic decomposition of Equation (9). Thus, the symmetry classes of Flex are: The space Flex is hence divided into 6 strata: The partition of Equation (12) is called the stratification of Flex. It remains to calculate the dimension of each symmetry class in Equation (11). There are formulas in [28] to do this. We recall these formulas in the following lemma.
Lemma 3 (see Lemma 3.10 of [28]). Let T n a space of n-th order tensors having the isotypic decomposition: where the notation p|k means p divides k and (T n ) H is the fixed-point space associated to H ⊂ O(2) defined as: Now applying these formulas to the space Flex, remembering that: we obtain: We can summarize the previous results in Table 2.
Remark 1. Notice that the use of these formulas does not allow the determination of the stratification of the space [28].
The associated closed strata satisfy the following inclusion relations: The connection between these closed strata is provided by the diagram given in Figure 1, where an arrow from O O Figure 1. Bifurcation diagram of Flex.

Remark 2.
Although the elasticity tensors and the flexoelectric ones are both fourth-order tensors, the structure of flexoelectric tensors is richer than that of the elasticity tensors. Indeed, it has been shown [26,31] (2)]. These classes characterize the so-called chiral-sensitivity.
Now by collecting the results given in the literature and those we obtained for Flex, we get all the symmetry classes of the complete space F lex: and then F lex is divided into 14 strata. Our goal in the following section is to determine the matrix representations of the complete flexoelectricity law, for each symmetry class.

Matrix Representations
It is common practice in continuum mechanics [34] to rewrite a constitutive tensorial law in matrix form for being able to implement it in a FEM code. This procedure requires the construction of appropriate bases which allows us to express state tensors as vectors and constitutive tensors as matrices. The matrix representations associated to tensors S ∼ , P , C ≈ , M , and A ∼ ∼ ∼ are given in the literature [22]. We will now determine the matrix representations associated to the tensor F ≈ .

Matrix Representations of the Flexoelectric Tensor
Since F ≈ is a fourth-order tensor in R 2 , it possesses generally 2 4 = 16 components. However, the symmetry relation F ijkl = F jikl reduces the number of independent components from 16 to 12. To rewrite the tensorial relation τ = F ≈ · q into a matrix form, we need to construct appropriate orthonormal bases for the vector spaces T (ij)k and T (ij)kl . Let B i be a basis for the space of i-th order tensors constructed from the canonical basis B := {e 1 , e 2 } of R 2 . For T (ij) we chose: For T (ij)k , we chose: The ordering for basis elements in B 3 is chosen so that the matrix forms of A ∼ ∼ ∼ are block diagonals for the dihedral classes [21]. Since this ordering is the one considered in the references [21,22], the results from these different papers can directly be combined. In the basis B 3 , the hyperstress tensor τ and the strain-gradient tensor η are represented by vectors [ τ] ∈ M 6,1 and η ∈ M 6,1 : For T (ij)kl we chose: In the basis B 4 , the flexoelectric tensor F ≈ is represented by a matrix F ∼ ∈ M 6,2 : We will now determine the matrix forms of the flexoelectric tensor in each stratum.

Z 2 -Class
The generator of the group Z 2 is R(π). The condition that (Q ∼ F ≈ ) ijkl = F ijkl is trivially satisfied because the minimal symmetry class of an even-order tensor is Z 2 . Indeed, the orthogonal transformation in the Z 2 situation is R(π) = −I. Moreover, in the transformation of an even-order tensor (cf. Equation (7)), the orthogonal transformation acts an even number of times. As a consequence R(π) is in the symmetry group of any even-order tensor or, said differently, any even-order tensor is centrosymmetric. We have 12 independent coefficients, as indicated below:

Z 4 -Class
The generator of the group Z 4 is R(π/2). The condition Q ∼ F ≈ ijkl = F ijkl yields: Thus, we have 6 independent coefficients, as indicated below:

SO(2)-Class
Consider a rotation of order n. For n > 4 the order of the rotation exceed the order of the tensor F Thus, we have 4 independent coefficients, as indicated below:

Remark 3.
As in [22], one can show that there exist rotations which allow us to reduce the number of parameters of F ∼ Z 4 and F ∼ SO(2) from 6 to 5 and from 4 to 3, respectively.

D 2 -Class
The generators of the group D 2 are R(π) and P(e 2 ). The condition (Q ∼ F ≈ ) ijkl = F ijkl yields: We have 6 independent coefficients, as indicated below:

D 4 -Class
The generators of the group D 4 are R(π/2) and P(e 2 ). The condition Q ∼ F ≈ ijkl = F ijkl adds to Equation (15) the new conditions: We have 3 independent coefficients, as indicated below:

O(2)-Class
Exploiting again the Hermann's theorem, O(2)-invariant tensors are obtained by imposing D 8 -invariance, with generators R(π/4) and P(e 2 ), to F We have 2 independent coefficients, as indicated below: We recover the number of independent components for each symmetry class of F ≈ given in Table 2.

Matrix Representations of the Complete Flexoelectric Law
The determination of the matrix representations of the complete flexoelectric law is done in two stages: First to express the sate tensors as vectors and the constitutive tensors as matrices in appropriate orthonormal bases and then to assemble these matrices. The orthonormal bases associated to the state tensors and the constitutive tensors S [ τ] [p] [ η] [q] For each symmetry class, the assembly of the sub-matrices is done taking into account two rules: 1. Any odd-order tensor vanishes for the centrosymmetric symmetry classes : 2. An even-order tensor can not see an odd-order material invariance (in 2D), the invariance seen will be twice the order of the former (see [36] (Theorem 5)).
Using the above rules, we obtain the following matrix representations of the complete flexoelectricity law in each stratum: It is important to notice that the matrix representations of the complete flexoelectricity law as given above do not provide a clear picture of their physical content (symmetry, coupling, etc.). For example, if we only focus on the relation p = F ≈ . . . η describing the direct flexoelectric effect, we can not predict the coupling between the different mechanisms: Stretch-gradient and rotation-gradient, that occur in strain-gradient elasticity [19]. The tool that will allow us to analyse such couplings is the harmonic decomposition, i.e., the decomposition of a tensor into a sum of harmonic tensors, via an explicit isomorphism [26]. We will rewrite the matrix representations of F ≈ in terms of the components of the harmonic tensors of the harmonic decomposition. This procedure will be called the harmonic parametrization of the flexoelectric tensor.

Harmonic Parametrization of the Flexoelectric Tensor
In Section 3.4, we have shown that the space Flex has the isotypic decomposition: In this subsection, we give an explicit isomorphism that realizes this decomposition. This decomposition, as provided here, is not the restriction of the one proposed in [17] to the 2D situation, but another, albeit isomorphic, decomposition. This harmonic decomposition is referred to as the Clebsch-Gordan harmonic decomposition, and aims at giving a precise meaning to its elements. This procedure will be detailed in a forthcoming paper [37]. We just provide here the main results.
Due to the multiplicities of the spaces K 2 , K 0 , and K −1 in the isotypic decomposition, the isomorphism that realizes this decomposition is not uniquely defined [28]. Among the different possibilities, some decompositions have more physical content than others. The decomposition we consider for the constitutive tensor is chosen in order to be compatible with the explicit harmonic decomposition of the state tensor space T (ij)k as proposed in the following subsection.

Decomposition of the State Tensor Space T (ij)k
We recall that Flex is isomorphic to T (ij)k ⊗ T l . The decomposition we choose here consists first of decomposing the state space T (ij)k into a sum of O(2)-irreducible spaces. We recall from Section 3.4 that the space T (ij)k has the following isotypic decomposition: An explicit decomposition is carried out as follows. We split any tensor T ∈ T (ij)k into a complete symmetric tensor S ∈ T (ijk) and a remainder V r as: with When applied to the strain-gradient tensor η , the decomposition of Equation (19) has a well-known mechanical interpretation in the strain-gradient elasticity theory [19,20]: The tensor S is the stretch-gradient tensor and the tensor V r is the rotation-gradient tensor.
It is easy to check that the tensor V r belongs to the space: and that the two subspaces T (ijk) and H r(3,1) are in direct sum. Since dim(H (3,1) ) = 2 = dim(K 1 ), we can associate a vector v r ∈ K 1 to the tensor V r via an embedding Φ ≈ r(3,1) : K 1 − → H r (3,1) . This embedding is given by: which can be rewritten in an intrinsic form as where i 1 ≈ , i 2 ≈ , and i 3 ≈ are the fourth-order isotropic tensors with components: In order to construct an O(2)-irreducible decomposition of the tensor T , we decompose the tensor S into the sum of a harmonic tensor K ∈ K 3 and an element V s belonging to the subspace: where P ∼ ∼ ∼ K 3 is the projector from T (ijk) onto K 3 defined by: the identity tensor on T (ijk) with components: Since dim(T (ijk) ) = 4 and dim(K 3 ) = 2, we deduce that dim(H s(3,1) ) = 2. As a consequence, we can associate a vector v s to the tensor V s via an embedding Φ ≈ s(3,1) : K 1 − → H s(3,1) , given by: Finally, an explicit O(2)-irreducible decomposition of the tensor T ∈ T (ij)k is given by: When applied to the strain-gradient tensor η , this decomposition has a clear interpretation: The vectors v s and v r represent the vector parts of the stretch-gradient tensor S and the rotation-gradient tensor V r , respectively.

The Harmonic Basis
With those results at hands, the harmonic parametrization of an element of T (ij)k can be obtained. Since dim(K n>0 ) = 2 in 2D, each of the harmonic tensors K , v 1s and v 1r has two independent components. Therefore, we can represent them as vectors in appropriate orthonormal bases for K 3 and K 1 . Let us denote by (K 1 , K 2 ), (v 1s 1 , v 1s 2 ), and (v 1r 1 , v 1r 2 ) the independent components of the tensors H , v 1s , v 1r in K 3 and K 1 , respectively. Since the space K 1 coincides with R 2 , elements v 1s and v 1r are standard vectors: For the space K 3 , we consider the basis K 3 := E 1 , E 2 with: (e 1 ⊗ e 1 ⊗ e 1 − e 2 ⊗ e 2 ⊗ e 1 − e 2 ⊗ e 1 ⊗ e 2 − e 1 ⊗ e 2 ⊗ e 2 ), In this basis, each tensor K ∈ K 3 is represented as a vector: Using the relations given in Section 5.1 and a formal calculation software, we obtain: Using this parameterization, any element T ∈ T (ij)k can be expressed in the basis B 3 under the following form: For the case of the strain-gradient tensor η , the splitting between the stretch-gradient and rotation-gradient contributions is not apparent when the stretch-gradient tensor S and the rotation-gradient tensor V r are expressed in the basis B 3 . To have a better insight, we have to define an appropriate orthonormal basis H 3 , called the harmonic basis associated to the basis B 3 . This harmonic basis H 3 is constructed from the explicit harmonic decomposition of T (ij)k . By applying the projectors defined by Equations (28) In this new basis, any tensor T ∈ T (ij)k has the following representation: In the above representation, for the strain-gradient tensor η , the first four components are associated with the stretch-gradient part of the strain-gradient while the last two components are associated with the rotation-gradient part of the strain gradient. Hence elementary mechanisms that work in the strain-gradient tensor appear clearly separated once expressed in H 3 . Now let us go to the decomposition of the constitutive tensor F ≈ .

Explicit Harmonic Decomposition of the Flexoelectric Tensor
From the tensors Φ ≈ s(3,1) and Φ ≈ r(3,1) defined in Equations (25) and (22), we construct a family of orthogonal projectors P ∼ ∼ ∼ 1s , P ∼ ∼ ∼ 1r , P ∼ ∼ ∼ 3 on H s(3,1) , H r(3,1) , and K 3 , respectively, as follows: where I ∼ ∼ ∼ Piez is the identity tensor on T (ij)k with components: An explicit harmonic decomposition of the flexoelectric tensor is given in the following result. The proof of this result will be given in [37].

Proposition 2.
A tensor F ≈ ∈ Flex admits the uniquely defined Clebsch-Gordan decomposition associated to the family of projectors P ∼ ∼ ∼ 3 , P ∼ ∼ ∼ 1s , P ∼ ∼ ∼ 1r : Those elements are determined from F ≈ as follows: Above, P ≈ 2 is the tensor defined by P To denote explicitly the components of F ≈ ∈ Flex in the harmonic decomposition in Equation (29), we use the compact notation: In this case, we define the O(2)-action of F ≈ as:

Associated Matrix Representations
Since dim(K n>0 ) = 2 in 2D, each of the harmonic tensors H components. Therefore, we can represent them as vectors in appropriate orthonormal bases for K 2 and K 4 . As in [38,39], for K 2 we choose the basis K 2 := {E 1 ∼ , E 2 ∼ } with: (e 1 ⊗ e 2 + e 2 ⊗ e 1 ), and for K 4 we choose the basis K 4 : the scalars α 1s , α 1r , and the pseudo-scalars β 1s , β 1r . Using the relations given in Proposition 2 and a formal calculation software, we obtain: By resolving the above system with respect to the components of the tensor F ≈ , we obtain the following matrix: According to the mechanical interpretations given in Section 5.1 and following the procedure introduced in Section 5.2, it would be interesting to highlight the coupling encoded by the flexoelectric tensor between the stretch-gradient and the rotation-gradient tensors. When the matrix representations of the flexoelectric tensor are expressed in B 4 , these couplings are not apparent. As in Section 5.2, to clear things up we have to define an appropriate orthonormal basis H 4 , called the harmonic basis associated to the basis B 4 . The vectors of the basis H 4 are constructed from the vectors of H 3 and B as: where we used the slight abuse of notation to label the vectors of the basis H 3 by the subscript I according to the order given in Equation (27). In the basis H 4 , the matrix of the flexoelectric tensor is written: . Now using the above matrix, we can see clearly the couplings encoded by the flexoelectric tensor. Indeed, the structure of the former matrix is as follows: Returning to the partition of a third-order tensor (Equation (19)), the first two blocks correspond to the coupling between q and S while the third block corresponds to the coupling between q and V r . In more physical terms, the basis H 4 makes apparent the fact that the matrix of the flexoelectric tensor can be split into two parts: 1. The first part F ∼ S relates the generation of double-forces [19,40] to the electric field. And, on its side, its transpose describes how the stretch-gradient generates an electric polarization; 2. The second part F ∼ R relates the generation of a couple-stress [41,42] to the electric field. And, on its side, its transpose describes how the rotation-gradient generates an electric polarization.
In a formula, we have: F For a Koiter continuum, also known as a constrained couple-stress continuum, the flexoelectric tensor reduces to the sole contribution of F ≈ R [43]. As a consequence, and conversely to a complete strain-gradient continuum, the flexoelectric tensor behaves like a second-order one with respect to orthogonal transformations. At the other end, a continuum for which the contribution of higher-order effects to the mechanical energy is restricted to the symmetric part of the strain-gradient tensor is given by the sole contribution of F ≈ S .
With this interpretation in mind, let us consider the flexoelectric matrices for the different symmetry classes: Using this decomposition, it can be observed, for instance, that the restriction of a strain-gradient model to a Koiter continuum, can not describe the tetragonal anisotropy.
A second interesting point is that the explicit harmonic decomposition allows us to construct a graphical representation of a flexoelectric tensor in R 2 by identifying its harmonic tensors with their associated vectors, as illustrated in Figure 2.
In this graphical representation, the vectors associated with the different harmonic components are transformed by an orthogonal transformation according to the irreducible space K n they belong to (cf. Equation (8)).

Discussion
In this work we solved the problem of the determination of the number and types of symmetry classes of the complete flexoelectricity law. Our work completed previous ones [21,22] where attention focused only on the determination of the number and types of symmetry classes of the strain-gradient elasticity law.
Furthermore, we introduced an explicit harmonic decomposition of the flexoelectric tensor. This decomposition helps to understand the nature of the different electromechanical couplings that took place in a flexoelectric material, such as the contribution of the stretch-gradient or the couple-stress to the electric polarization. This decomposition will be of prime importance for the determination of an explicit basis of polynomial invariants for the flexoelectric tensor.
Using the introduced parametrization, the next step, object of a future contribution will consist in the construction of a finite set of such invariants. This set of invariants will then be used to formulate objective functions for the optimization algorithm.