On nonlinear dilatational strain gradient elasticity

We call nonlinear dilatational strain gradient elasticity the theory in which the specific class of dilatational second gradient continua is considered: those whose deformation energy depends, in an objective way, on the gradient of placement and on the gradient of the determinant of the gradient of placement. It is an interesting particular case of complete Toupin–Mindlin nonlinear strain gradient elasticity: indeed, in it, the only second gradient effects are due to the inhomogeneous dilatation state of the considered deformable body. The dilatational second gradient continua are strictly related to other generalized models with scalar (one-dimensional) microstructure as those considered in poroelasticity. They could be also regarded to be the result of a kind of “solidification” of the strain gradient fluids known as Korteweg or Cahn–Hilliard fluids. Using the variational approach we derive, for dilatational second gradient continua the Euler–Lagrange equilibrium conditions in both Lagrangian and Eulerian descriptions. In particular, we show that the considered continua can support contact forces concentrated on edges but also on surface curves in the faces of piecewise orientable contact surfaces. The conditions characterizing the possible externally applicable double forces and curve forces are found and examined in detail. As a result of linearization the case of small deformations is also presented. The peculiarities of the model is illustrated through axial deformations of a thick-walled elastic tube and the propagation of dilatational waves.


Introduction
Albeit the simplified version of continuum mechanics (where deformation energy depends on the first gradient of placement, as established by Cauchy and Navier), has been the most studied, exploited and applied in Engineering and Physical sciences, already in the pioneering works by Gabrio Piola (active from 1822 to 1850) [24], which preceded Cauchy's works, it is clearly stated that it is logically possible to conceive more general continuum models. Piola, by means of a careful analysis including also a micro-macro identification procedure, did also find some solid physical motivations urging the need of considering more general continua. Piola's argument clearly establishes the properties to be verified, at micro-level, that impose the adoption, at macro-level, of deformation energies depending on higher gradients of displacement field. Among these, the most relevant seems to be: (i) the complexity of the geometry of the involved micro-structure and (ii) the sharp space variation, at the micro-scale, of the different relevant mechanical material properties, as, for instance, elastic stiffness.
It is to be yet clarified why the intuitions and the results by Piola were not continued until very recently [26]: we do not believe that the reasons reside in some lack of mathematical tools available or in the exceptional and visionary capacities shown by Piola in his works. The second gradient continuum models were more recently revived in the effort of establishing more general models in continuum mechanics (see discussions in [8,14,27,[82][83][84]) and for giving a solid mathematical basis to the effort aimed to design novel metamaterials, see, e.g. [1,9,33,34,102]. The basic idea on which Piola based the theory of second or higher gradient of displacement field can be rather simply explained: when formulating the Principle of Virtual Work, or its particular case represented by the principle of minimum of energy in stable equilibrium configurations, one has to consider the expression of internal work (or of deformation energy) on a virtual displacement field as depending on all its gradients up to order N .
In case of elasticity, in order to generalize Cauchy models, one may assume that the deformation energy may depend on all the gradients of displacement up to order N . As a consequence all other related constitutive equations, in such generalized theories, must depend not only on the first gradient of placement, but also on all its necessary higher-order gradients. The complete study of the objectivity of deformation energy dependence on higher gradient of placement was performed already by Piola in his fundamental work [25].
After works by Korteweg [71], Cahn and Hilliard [15,16] dealing with so-called capillary fluids, by Mindlin [87][88][89] and Toupin [128,129] dealing with more general continua, including elastically deformable solids, and after the more general approaches proposed by Germain [62,63] and Sedov [113], who considered microstructured generalized continua, see also [39], the second gradient elasticity was discussed and enlarged in many works, e.g. [7,12,32,54,134,138]. It has seen a remarkable growth and seems to have been accepted also by those groups of scholars who had initially doubted about its importance.
The growing interest paid to the second gradient elasticity relates to its success in modelling several emerging phenomena. A partial list of some among them can be outlined here: (i) the mechanical behaviour of solids and fluids at the nanoscale, see, e.g. [5,22,58,65,70,79,116] and the references therein, albeit it has to be yet understood in which sense the concepts of quantum mechanics or lattice dynamics are related to the continuum modelling of these specific phenomena; (ii) capillarity phenomena in fluids [18][19][20]31,43]; (iii) material behaviour of multi scale media endowed with complex internal micro-architecture such as beamlattice architectured materials and novel metamaterials [2,6,29,[33][34][35]85,86,91,99,105,107,120,135]. (iv) deformation localization phenomena and dissipation in plasticity, see, e.g. [14,56,57]; (v) damage and fracture phenomena [100,101,103,111,124]; (vi) growth of bone tissues and its effects on bone mechanical properties [60,61,76,118]; (vii) coupled electromechanical phenomena as flexoelectricity and flexomagneticity, see, e.g. [37,46,77,78,80,132,136] and the reference therein; Albeit the previous list seems rather limited and definitely incomplete, it relates to the research fields that are closer to the results presented in this paper. As in the case of simple (Cauchy) continuum models, the crucial problem in material modelling remains, also for second gradient continua, the determination of material parameters to be used when describing a specific mechanical system. One should not forget, while presenting the results concerning second gradient materials, the long debate, mainly occurring in Paris between Navier, Cauchy, Poisson and Lamé, concerning the number of elastic constants to be used in linearized elasticity for isotropic materials. This debate was concluded when George Green introduced, in the debate, the powerful tools supplied by Lagrangian Variational Principles, see [66,73].
As clearly reported by Benvenuto et al. [10] the controversy about isotropic elastic constants was not so easy to settle down, see also [125][126][127]. It is therefore not surprising if a similar debate was started about second gradient isotropic elasticity: an account of available results and approaches is given in [32].
In fact, the formulation of the most general possible constitutive equations for linearly elastic isotropic material in the framework of second gradient continua, as already established by Toupin and Mindlin, requires, in addition to two Lamé classical elastic moduli, the introduction of five more elastic moduli, that introduce constitutively five characteristic lengths.
Seen the initial difficulties in determining the most general form of deformation energy for second gradient continua and the difficulties in finding meaningful analytical solutions in which all the five possible characteristic length appearing in linearized models are simultaneously relevant, in the literature one finds various simplifications of second gradient continuum models, in which less independent material parameters are introduced, see, e.g. [4,123,134,138]. It has to be explicitly stated here that some authors tried to justify the study of a particular class of second gradient continua as determined by a smaller number of constitutive stiffness parameters on the basis of an application of the lex parsimoniae (law of parsimony) or Occam's razor, see [36,Chap. 8].
Of course this claim is not justified, in general. Of course we remember Einstein's quote "Everything should be made as simple as possible, but not simpler". On the other hand, one cannot hope to model a complex physical system with an oversimplified mathematical model. Lex parsimoniae is a logical tool to be used in choosing a postulation scheme for a mathematical model or a set of axioms for a mathematical theory. It cannot be evoked for justifying the introduction of a too simple model for a very complex physical system. If one is not able to handle a complex model, they cannot hope to forecast the behaviour of a correspondingly complex physical system. Therefore, the objections sometimes raised against second gradient models or other generalized models of continua (the claim that these theories require too (!) many material parameters) are totally unjustified from the point of view of a modeller who wants to describe complexity. These objections do not seem to be based on a clear Archimedean argument: indeed one has to choose the logical structure of a mathematical model, then they develops it mathematically and only finally finds the experimental evidence that they is capable to describe. And complex systems obviously need complex models.
As an example of such complex physical system one can consider capillary fluids [59,68]. It is now universally accepted that the capillary fluids can be described by means of second gradient fluid models, since they allow for the prediction of all main capillarity phenomena observed in classical experiments. Being rich enough, second gradient fluid models allow for the introduction of some extra essential or natural boundary conditions, are characterized by some more constitutive coefficients and therefore can predict the behaviour of capillary fluids. All the modelling efforts spent in inventing the second gradient fluid models were based on a very basic assumption: the principle of minimum of energy. A logically less important assumption was to assume that deformation energy depends only on the current mass density and on the norm of the Eulerian gradient of Eulerian density. The suitable boundary conditions to be used for second gradient fluids were determined by means of a mathematical argument (variational approach), and only at the end of the full development of the mathematical theory, it was necessary to interpret physically the introduced boundary conditions. Never in history of science one can observe that a really successful model was formulated by recurring continuously to "physical evidence". Instead successful models are based on a set of postulates (chosen with Occam razor criterion), developed by deducing from the postulates a series of consequences and finally checking some of these consequences with experimental evidence. Only when a model is fully developed, the interpretation of all mathematical results becomes possible by using physical intuition. In [31,115], the previously described modelling process is fully described and all mathematical entities introduced in the process are fully justified when comparing the obtained prediction with the experimental evidences and with the physical phenomena occurring in capillary fluids. More specifically, in [114], the novel boundary conditions involving the normal gradient of mass density at a wall in contact with a capillary fluid is interpreted ex novo as a powerful predictive capacity of the introduced novel model. In fact, wettability of the wall is interpreted as the specific physical property that the novel model can describe and that was not accounted for in previous models. The boundary conditions of the normal gradient of mass density are introduced first on logical grounds and only subsequently they are physically interpreted. More recent discussion of novel boundary conditions that can be introduced in second continuum theories can be found in [21,[41][42][43]47,50,109]. In particular, the interrelations between Toupin-Mindlin strain gradient elasticity with the surface elasticity by Gurtin-Murdoch were analysed in [50], whereas similar comparison of surface phenomena within the lattice dynamics was provided in [51,117].
In this paper, we discuss a particular class of second gradient models that we want to be the object of dilatational second gradient elasticity. The model is based on the introduction of a second energy density U as a function of the placement gradient and of the gradient of its determinant.
Characterizing the place of the discussed model among others let us briefly recall constitutive relations for (hyper)elastic continua: -Cauchy continuum (so-called simple material) -Strain gradient elasticity (Toupin-Mindlin material) -Korteweg or Cahn-Hilliard fluid (gradient fluid) -Dilatational strain gradient elasticity Here, F is the deformation gradient, ρ is a current mass density. The other notations have an obvious meaning: they are introduced in the next sections.
So we can state that the dilatational second gradient elasticity is a reduction in the general second gradient elasticity to the case where simple materials are behaving, only for what concerns second gradient energy dependence, as second gradient fluids. Equivalently one can regard it as an extension of the second gradient fluid model obtained by "solidifying" the first gradient part of deformation energy. Indeed, we can see that the discussed deformation energy density can be treated as a combination of energies for simple materials and second gradient fluids.
It has also to be remarked that the proposed model has also a strong relation to the so-called generalized models of continuum with scalar (one-dimensional) microstructure by Capriz [17] or Eringen [52,53], see also [23,64,94,104,106]. So it can be used for various applications to material modelling of pressure sensitive materials as those studied in [95,131] or to study some phenomena of phase transitions occurring in solid materials [69,110].
In particular, the model could be used for modelling pressure induced phase transitions as it was done for second gradient fluids [30,31,110] but taking into account shear deformations.
The paper is organized as follows. In Sect. 2, we briefly introduce the basics of the second gradient elasticity for media undergoing finite deformations. In the general case for second gradient materials, a deformation energy density depends on two independent deformation measures that are the so called Cauchy-Green strain tensor C and the third-order tensor K related to the gradient of F. As a particular class of nonlinear second gradient media, we focus here on the dilatational second gradient elasticity model. Within this model the deformation energy depends only on C and on the gradient of the determinant of the placement gradient k. The general constitutive equation in this class is presented and few classes of anisotropic materials are discussed. In Sect. 3, the comparison with strain gradient fluids is given. In Sect. 4, the linear dilatational second gradient elasticity is introduced as in [48]. The comparison of the dilatational strain gradient elasticity with other models of continua with scalar microstructure as the Cowin-Nunziato poroelasticity [23,94] is then given. We show that the dilatational second gradient elasticity can be derived using the Lagrange multipliers technique. Using the Hamilton-Ostrogradsky variational principle, in Sect. 6, we extend the models for dynamics. Finally, we discuss linear dilatational waves and plane axial deformations of an elastic thick-walled tube in Sect. 7.
In the following, we utilize an index-free tensor calculus notation as in [44,74,93,119,133]. So vectors, second-and higher-order tensors will be denoted all by boldface symbols.

Basic equations of nonlinear strain gradient elasticity
Let us consider a finite deformations of an elastic body B within the strain gradient elasticity. As in nonlinear elasticity, a deformation of B is described as an invertible smooth enough mapping from a given reference placement = (B) into a current (deformed) placement χ = χ(B)(t), where t is time. For a given material particle x ∈ B, we denote its position vectors in and χ as X and x, respectively. So a deformation → χ is given by We introduce the corresponding deformation gradient F and its gradient G by formulae where ∇ is the three-dimensional (3D) nabla-operator. For example, in Cartesian coordinates X M with corresponding unit base vectors i M , M = 1, 2, 3, we have where Einstein's summation rule is adopted and ⊗ denotes the dyadic product. For a hyperelastic material there exists a strain energy density function U. Let us introduce it in a general form assuming its dependence on all possible arguments Here, we consider dependence on X for inhomogeneous materials. The principle of material frame indifference [130] states that U must be invariant under changes where a and Q are a constant vector and a constant orthogonal tensor, respectively, and the centred dot stands for the scalar (inner) product of two vectors. Indeed, a and Q may depend on t but this does not change further derivations. So we get the identity U(x, F, G; X) = U(x + a, F · Q, G · Q; X) ∀ a, Q.
From (5), it follows the standard conclusion that U does not depend on x, whereas dependencies on F and G transform into where for brevity we hide dependence on X and T denotes a transposed tensor. Indeed, let us consider the polar decomposition of the deformation gradient [44] F = U · R, where U = (F · F T ) 1/2 is a symmetric second-order tensor and R is an orthogonal tensor. Substituting into (5) Q = R T , we get see [12,44] for more details. This equation can be transformed into where C is the Cauchy-Green strain tensor, and K is the third-order tensor which is the second Lagrangian strain measure in the nonlinear strain gradient elasticity. Note that we keep the same notation U for different functions. Another form of the constitutive equation as was discussed in [32]. Reiher and Bertram [108] also used another Lagrangian strain measure defined as follows: so the constitutive equation takes the form In the framework of the strain gradient elasticity, the equilibrium equation can be derived using the Lagrange variational principle as in [8,43,44], see also "Appendix". The Lagrangian equilibrium equations takes the form [44] ∇ · T + ρ f = 0, (8) where ρ is the mass density in , f is a vector of mass forces, and T is the total stress tensor of the first Piola-Kirchhoff-type. It is given by where P is the first Piola-Kirchhoff stress tensor and M is the first Piola-Kirchhoff hyperstress tensor, which is a third-order tensor. Note that by its definition M has the same symmetry as G = ∇ ∇ x. Introducing the transposition with respect the first two indices of a third-order tensor, we denote this symmetry as So as an example of such third-order tensor, we can consider In what follows, we assume that A consists of a finite number of smooth faces A j separated by edges L m , where normal N undergoes a discontinuity jump or where external line forces are applied, see Fig. 1. The corresponding natural (static) boundary conditions on A = ∪ j A j have the form where N is the vector of unit outward normal to the boundary ∂ V , H = − 1 2 ∇ S · N is the mean curvature of ∂ V , ∇ S is the surface nabla-operator defined as 1 is the 3D unit tensor, t is a vector of the surface traction, and c is a vector of external double forces, see [44, pp. 288-293] for details.
Using (12), we can simplify (11). Indeed, substituting (12) into (11), we get This means that for a curved boundary (H = 0), a surface double force contributes to an external surface traction. Unlike nonlinear elasticity of simple materials, we get also non-trivial boundary conditions along edges. Let L m be an edge between faces A + and A − . Then, we have where t m is a vector of external edge forces defined along edges, N ± are normal vectors to A ± , and ν ± are vectors of unit outward normals to L m lying in the tangent plane to A ± , see Fig. 2.
Let us notice that in the strain gradient elasticity external line forces are admissible, so we can also assign (14) to any surface curve L on ∂ V . In this case, we have that N ± = N and ν + = −ν − = ν, so we transform (14) into a jump compatibility condition where q is a vector of line forces given along L.

Fig. 2 Two faces A ± with common boundary L m
The boundary-value problem (8), (11), (12), and (14) can be also transformed using the Eulerian description. First, we introduce the Cauchy-type total stress tensor T by the relation Using Piola's identities one can easily prove that for any tensor-valued field T, we have where ∇ is the nabla-operator in χ related to ∇ by For example, in Cartesian coordinates x k with corresponding base vectors i k it has the form Using (18), we transform (8) into where ρ = J −1 ρ is a mass density in the current placement χ . Similarly to (16), we introduce the Cauchy stress tensor As M is a third-order tensor, its transformation into Cauchy form is more complex. First, we introduce the Rayleigh product * for dyads and triads as follows: Obviously, it can be easily extended to tensors of any order. Then, the Cauchy hyperstress tensor is given by the relation For the derivation of the Eulerian form of the boundary-value problem, one can use a transformed Lagrange variational principle as described in "Appendix". As a result, one gets x ∈ a f : In (24) Similarly, a f , m , and are images of A f , L m , and L under transformation → χ, a f = χ(A f ), m = χ(L m ), and = χ(L), respectively. In addition, : denotes the double dot product. For two dyads, it is defined as and can be easily extended for higher-order tensors [44,93]. For example, we get

Dilatational gradient elasticity
The peculiarity of the dilatational strain gradient elasticity is a reduced dependence on second gradient of deformation G. Here, we assume that U depends on G only though the gradient of determinant of F, i.e. on ∇ (det F). In this case, we have Applying again the principle of material frame indifference to this dependence, one arrives at where U is an even function with respect to ∇ J : As a result, we consider the strain energy density function within the dilatational strain gradient elasticity in the following form: where U is an even function of k, U(C, k) = U(C, −k) for all C and k in the domain of U. Let us consider an isotropic material. The material symmetry analysis for strain gradient media was developed in [38,92,108]. In case of isotropic media, the strain energy density has the following invariance property In other words, U is an isotropic function of its arguments. Using the theory of invariants [121,137], we get the following representation of U where tr is the trace operator. Instead of traces of C p , p = 1, 2, 3, its principal invariants I 1 = tr C, can be also used.
As an example of a constitutive relation, we can consider the following extension of the Saint Venant-Kirchhoff model with or the harmonic (John) material model where λ and μ are Lamé elastic moduli and α is an additional elastic modulus related to the gradient of dilatation.
Note that for isotropic media, we consider in (32) orthogonal tensors. If we restrict ourselves to proper orthogonal tensors, i.e. rotation tensors, we come to the model of a hemitropic material. Here, we have an additional invariant [137] and the representation of the strain energy density of a hemitropic material reads where × stands for the cross product.
Considering other symmetry groups we can find the related reduced representations of the constitutive equations for anisotropic media. For example, let us consider a transversely isotropic solid. Let e be a unit vector defining an axis of transverse isotropy. Then, a strain energy density should be an isotropic function of C, k and e or e ⊗ e if two directions e and −e are undistinguishable. For the last case, we get the representation in the form 2 , e · C · e, e · C 2 · e, (k · e)(e · C · k) .
This representation corresponds to the one of five possible transverse isotropy groups in the three-dimensional space, i.e. to the group D ∞h for which the structural tensor is e ⊗ e, see [137, p. 563] for details. Orthotropic materials constitute another often observed class of anisotropic materials such as composites or rocks. Here, we have three orthotropic symmetry groups [137]. In particular, the group O 3 includes reflection transformations in planes normal to orthogonal unit vectors e k and the mirror reflection. For this group, the structural tensor is given by S = e 1 ⊗ e 1 − e 2 ⊗ e 2 and the representation of the strain energy density has the form [137, p. 566] In a similar way, we can consider the representations of constitutive relations for other anisotropic media. Let us note that as one faces here with function of a second-order tensor and vector, invariants-based representations of a strain energy density is similar to fibre-reinforced composites [122] of simple materials and simpler than in the general Toupin-Mindlin strain gradient elasticity or even in the micropolar elasticity [49,55].
Equation (31) is a reduced form of (6). Indeed, in (31) a dependence on K appears trough k. Let us express k through K explicitly. First, one compares the identities and get that For brevity we denote the derivative of J with respect to F as J ,F . Here, we used the formula for differentiation Using K = ∇ F · F −1 introduced by (7), one gets another representation of k: With (40) we can transform the general constitutive equations of the strain gradient elasticity (10) for the case of dilatational strain gradient elasticity. To this end, let us consider dU. From (10), it follows that where . . . stands for the triple dot product. For example, for triads, it is given by On the other hand, one gets where m is the first Piola-Kirchhoff-type double force vector. Using the identities one arrives at Comparing (41) and (42), one has Since dF can be treated as an arbitrary second-order tensor from (44), it follows As it was already mentioned, G has the symmetry property G T (1,2) = G, so the same symmetry is inherited by dG. Thus, from (43), it follows that M coincides with the symmetric part of J m ⊗ F −T :

Lagrange variational principle
Let us find the Lagrangian equilibrium conditions within the dilatational strain gradient elasticity using directly (31). To this end, the Lagrange variational principle is applied by using the same technique as in [40,43]. Let V be a domain occupied by B in a reference placement, V = (B). The total energy stored in V is given by Equilibrium conditions follow from the variational equation where δ is the symbol of variation and δA is a work of external actions which will be specified later. Let us note that δA does not constitute the first variation of any functional, in general. Calculating the first variation of E, one gets where v = δx and P = ∂U ∂F , m = ∂U ∂k (50) are the first Piola-Kirchhoff stress tensor and the first Piola-Kirchhoff double force vector, respectively. Using the identity and integrating by parts in (49), one gets Integrating again by parts as follows: In order to transform the last integral in (52), the decomposition of ∇ v into normal and tangent parts is used as follows: In order to integrate by parts in (53), the surface divergence theorem is used, which consists of the relation for any differentiable tensor-valued field Z and any smooth enough surface A with a contour L. Here, H is the mean curvature and ν is the outward unit normal vector to L, such that N · ν = 0. Using (54), one has the following formula of integration by parts on a surface With (55), the last integral in (53) is transformed as follows: where [[(. . .)]] m means a jump across L m . For example, if L m is a common boundary between two faces A − and A + as shown in Fig. 2, one has where m and J F −T are also calculated as unilateral limits at L m . As a result, one gets the formula for δE Obviously, δ E includes integrals over domain V , its surface A = ∂ V and edges L m , where the surface integrals are linearly dependent on v and its normal derivative. So one can assume the same form for δA In (58) the natural (static) boundary conditions on A f and the static conditions along L m and L, respectively, Introducing the total stress tensor of the first Piola-Kirchhoff-type T as Eq. (59) is transformed into the classic form Using (64) and (61), one can also simplify (60) as follows: whereas (62) and (63) become Obviously, from (66) it follows that the scalar surface double force gives a contribution to the apparent surface traction t. Conditions (67) and (68) constitute a certain compatibility condition for applied loadings g, t m and q.
Let us reformulate the boundary-value problem using Eulerian description. To this end, one can transform δE and δA using Eulerian coordinates. Changing coordinates X → x in (49) and using (19) and (21), one has where v = χ(B) is the volume of body B in the current placement χ . As δ J = J ∇ · v, one gets also where m is the Cauchy double force vector introduced as follows: So δE takes the form Now the derivation performed above for Lagrangian coordinates can be mimicked. Integrating by parts in (72), one gets where n is the vector of unit outward normal to ∂v. Using formula 5 from [74, p. 19] and (18) For further derivations, the decomposition of ∇ on the surface is introduced ∇ = n ∂ ∂n + ∇ s and the surface integral in (73) is transformed as follows: Using the surface divergence theorem formulated in Eulerian coordinates, one rewrites the integration by parts formula (55) as Here, h is the mean curvature of a, h = −1/2∇ s · n, η is the outward normal vector to ∂a such that η · n = 0. With (77), one obtains the formula Here, m = χ(L m ) are images of edges under mapping → χ . Finally, one gets the Eulerian form of δE given by In a similar way, we δA is transformed where da d A and ds dS are the elementary changes in surface area and in curve length, respectively. Now the Eulerian boundary-value problem within the dilatational strain gradient elasticity takes the following form x ∈ a f : x ∈ : [[J (n · m)η]] = q dS ds . From the last equation, it follows that for g = 0 an action of a line force may produce an edge in a current configuration. From the Eulerian statement, it is clear that the dilatational strain gradient phenomena result in appearance of normal forces on a free surface and "hydrostatic" pressure in the equations for the bulk. In Sect. 5 our model is compared with models of materials with voids where pressure plays an important role.

Reduction to gradient fluids
As it was mentioned in Introduction the dilatational strain gradient elasticity is closely related to the model of capillary fluids [12,44,116]. In order to discuss this matter in more details, let us consider another form of the strain energy function. One can use the decomposition of F into its dilatational and distorsional parts Here F describes isochoric (distorsional) deformations. Such decomposition is used for description of nearly incompressible materials [96]. With this decomposition, one can write the constitutive equation (29) in the form U = U (F, J, k).
For isotropic materials instead of (33), one can use the following representation where I 1 = tr C, Since the referential ρ and current ρ mass density relates to each other through the mass balance Here, the following substitutions were performed So one can see that a strain energy density of capillary fluid U K is a dilatational part of U measured per unit volume in a current configuration Using (89), one can see that the dilatational strain gradient elastic material may demonstrate the same behaviour as a capillary fluid. In particular, for dilatational deformations interfacial layers of finite thickness may appear between two phases as observed in [31,110]. For further extensions towards viscous gradient fluids, one is referred to [13,72]. Similarly, the dilatational strain gradient elasticity can be extended towards models with viscosity.

Linear dilatational strain gradient elasticity
For infinitesimal deformations, the vector of displacements is introduced Here, there is no difference between and χ , and one has the following approximations ∇ ≈ ∇, F = 1 + ∇u, C ≈ 1 + 2ε, K ≈ K ≈ G = ∇∇u, where ε and ε are the linear strain tensor and the linear dilatation. For an isotropic material (33) transforms into where λ and μ are Lamé moduli and α is an additional elastic modulus. So unlike the Toupin-Mindlin strain gradient elasticity with five additional moduli here one has only one related to higher gradients. For a general anisotropic material one has a strain energy expressed by a quadratic form where Λ is a forth-order tensor of elastic moduli which has the same symmetries as in linear elasticity, i.e. 21 independent elastic moduli, Γ is a symmetric second-order tensor (three moduli), and Υ is a third-order coupling tensor having 18 independent components. So for general anisotropy the total number of elastic moduli becomes equal to 42. For centrosymmetric materials, Υ = 0, and one has no coupling.
For infinitesimal deformations, one has also that stress tensors and vectors coincide to each other as follows: For isotropic materials, one gets Unlike finite deformations, for linear isotropic solids, stresses and hyperstresses are decoupled, as the stress tensor depends only on strains, whereas the double force vector depends on the gradient of dilatation. Considering the representation of T as a sum of spherical and deviatoric parts, we can see that one can see that the deviatoric part does not depend on the second gradients of displacements, whereas the spherical part depends on the gradient of the divergence of displacements Here, dev means the deviatoric part of a tensor. So the discussed model is a particular class or strain gradient model proposed by [138]. Relation (46) between the third-order hyperstress tensor M and the double force vector m takes now the form For a homogeneous material, the equilibrium equation (65) is The well-posedness of boundary-value problems within the linear dilatational strain gradient models was studied in [48].

Cowin-Nunziato poroelasticity with constraints
As it was mentioned in Introduction, the above-discussed dilatational strain gradient elasticity has some relations to the so-called generalized models with one-dimensional microstructure [17]. Among such models, it is worth noticing the poroelasticity model proposed by Cowin and Nunziato [23,94]. Let us briefly recall the basic equations of this model. In addition to x, one has an additional kinematic descriptor ϕ which describes the voids fraction. So the strain energy density has the form In the following a particular class of materials is considered with W given by Applying the principle of material frame indifference to (92) one gets Calculating the first variation of the energy density functional one gets where v = δx and ω = δϕ. The consistent form of the virtual work of external actions is given now by where the extrinsic equilibrated body force c and the equilibrated surface stress τ have been introduced. From the variational equation one derives the Lagrangian equilibrium equations and the natural boundary conditions In order to underline the similarity of the model with scalar microstructure and the dilatational strain gradient elasticity, let us introduce the following constraint For small deformations, it takes the form ϕ = tr ε. Now, one can apply again the variational approach using the Lagrange multiplier technique. The modified functional has the form where Λ is a Lagrange multiplier. Repeating the calculation of the first variation of E CN , one gets Here, the identity δ J = J F −T : ∇ v has been used. Now using (102) from (97), one gets (100) and the modified equilibrium conditions Using (104), Λ can be excluded from (103) and (105). As a result, one gets Λ = ∇ · μ + ρ c and (103) transforms as follows whereas (105) takes the form Introducing the new mass force vector f by the relation one can see that (107) takes, up to notations, the same form as (65) where T = P − (∇ · μ )1.
As a result, one can conclude that the relation between the Cowin-Nunziato poroelasticity and the dilatational strain gradient elasticity is the same as the correspondence between the model of micromorphic medium and the strain gradient elasticity. For application of the Lagrange multipliers technique in the case of higherorder and generalized models of continuum, one is referred to [11,28].

Dynamics
Using the least action principle (Hamilton-Ostrogradsky principle) [8,43,44], one can extend the boundaryvalue problems to dynamics. To this end the kinetic energy density is introduced as follows: where the overdot stands for the derivative with respect to t and ς is the microinertia scalar. Within the Cowin-Nunziato model ς is called the equilibrated inertia [23,94]. In general, ς may depend on the body deformations. Following motivations [23,94], it is assumed that ς depends on volumetric strains, i.e. ς = ς(J ). The least action functional H takes the form where t 1 and t 2 are two time instants. Considering the variational equation for kinematically admissible variations v = δx with constraints v t=t 1 = 0, v t=t 2 = 0, the governing equations within the dilatational strain gradient elasticity in dynamics are derived.
Since the first variation of E has been already derived, one has to concentrate on δK. Using the identities and integrating by part one gets In (116) Within the Eulerian description this formula takes the form As a result, from (113) and (116) or (117) one gets the equation of motion and dynamic boundary conditions in both Lagrangian and Eulerian descriptions: or x ∈ a f : n · P − ∇ s · [J (n · m)] − 2h J (n · m)n (120) Unlike classic elasticity here the boundary conditions (119) and (121) for stresses contain inertia terms. Note that other boundary conditions for m and m did not change.
In the case of small deformations the equation of motion for a homogeneous isotropic solid takes the form where ς > 0, i.e. it is assumed to be a positive constant.

Dilatational waves
Let us consider a wave propagation in an elastic space. Within the dilatational strain gradient elasticity, the higher-order terms in the strain energy density have influence on dilatational waves (P-waves), whereas shear waves (S-waves) are the same as in classical linear elasticity [3]. Introducing a dilatational potential Φ as follows u = ∇Φ, (122) is transformed into a scalar equation Here, Δ = ∇ · ∇ and it is assumed that f = 0. Looking for a solution in the form where k is the wave vector, ω is the frequency, i = √ −1, and Φ 0 is a constant wave amplitude, one arrives at the dispersion relation, i.e. to an equations which relates ω and k, where k 2 = k · k. Introducing the phase velocity c = ω/k, one gets the dispersion relation in the form where c p = √ (2μ + λ)/ρ is the speed of longitudinal waves in an elastic space [3]. Eq. (125) contains two length-scale parameters, that are the static characteristic length s and the dynamic characteristic length d defined as With these notations the dispersion relation takes the form Obviously, if s = d a dilatational wave is a dispersive wave. For long waves, i.e. as k → 0, c ≈ c p . For short waves (k → ∞), one has the asymptotic value c = c p d / s . For other gradient-type models, such dispersion relation was analysed in [7] in more detail and for the couple stress theory in [90].
Characteristic dispersion curves are shown in Fig. 3. Here, the dashed horizontal line corresponds to the classic relation c = c p , the dotted horizontal lines show the asymptotes c = c p d / s . In addition dimensionless variablesc = c/c p andk = s k have been introduced. The group velocity c g = dω/dk is shown in Fig. 4.

Axial deformations of an elastic tube
As an example, let us consider the plane axial deformation of an elastic tube. Let a and b be the inner and outer radii of an elastic cylinder in a reference placement, respectively, see Fig. 5. Hydrostatic pressures and constant double forces are assumed on the boundaries. In the following, we use the polar coordinates in both reference and current placements. So an axial deformation has the following form [74,130] where r (R) is an unknown function. The deformation gradient is given by where e R , e Φ , e Z and e r , e φ , e z are unit base vector related to Lagrangian and Eulerian polar coordinates, respectively, and the prime stands for the derivative with respect to R. So one has also In the following attention is restricted to infinitesimal deformations. So one has where u = u(R) is a displacement. For the axial deformations, the equilibrium equations transform into the following scalar equation:

Fig. 5 Cross section of an elastic tube in a reference placement
With α = 0 it has simple solution given by where C 1 and C 2 are integration constants, see [75]. For α = 0 the solution of (131) is more complex. First one gets the solution for ε as follows where I 1 and K 1 are the modified Bessel functions of the first and second kinds, respectively, γ = √ (2μ + λ)/α, and C 3 and C 4 are other integration constants. Note that ξ = γ −1 plays a role of a length-scale parameter which characterizes the influence of gradient of dilatation. Then, integrating this equation using the properties of Bessel functions one gets u(R) For the sake of simplicity, it is assumed that the internal boundary (R = a) is free, whereas at R = b, one has the action of both hydrostatic pressure p and double force g. As a result, for axial deformations, one gets the following boundary conditions R = a : αε = 0, 2μu 2μu Boundary conditions (134)-(137) constitute a system of linear algebraic equations for the integration constants which depend linearly on p and g. Obviously, if g = 0 from (134) and (136) one has C 3 = C 4 = 0 and u has the classical form (132). On the other hand, if p = 0 one has a nonzero surface traction related to g, see (137). In the next figures, the normalized distributions of stresses, double stresses and total stresses for various of thickness and the material length-scale parameter ξ are presented. These stresses are given by In the following dimensionless variables R = R/b, P = P Rr /μ, T = T Rr /μ,m = m R /(bα),p = p/μ, g = bg/μ have been used. In addition, it was assumed a value λ = 1.5μ that corresponds to Poisson's ratio The stress distributions for the particular case − p + 1 b g = 0,ḡ = 0.01, are given in Fig. 7. Here, one has P R R = 0 at R = a, b, but since g = 0 one gets nonzero distribution of stresses.
In Fig. 8, the normalized double stressm = m/(αb) is presented as a function of R for various values of a/b and various values of ξ . Clearly, the presented solutions have the form of a boundary-layer-type solution, i.e. localized near the boundary R = b. This behaviour is more pronounced for thick cylinders, see, e.g. Fig. 8a, b.
Finally, in Fig. 9, the normalized total stress T as a function of R is presented. Other parameters are the same as in Figs. 6 and 7.
For thick-walled cylinders, one can observe strong localization phenomena in stresses, double stresses and total stresses, see Fig. 10, where a/b = 0.01. Here, boundary-layer-type solutions for stresses and total stresses in the vicinity of the free surface R = a are observed, whereas double force is localized near R = b.

Conclusions and research perspectives
In this paper, it is proposed to call "nonlinear dilatational strain gradient elasticity" the theory concerning second gradient continua where the deformation energy depends on: (i) the gradient of placement (via its objective part) and (ii) the gradient of the determinant of the gradient of placement (which is obviously objective).
The class of deformation energies considered here is singular, in the sense given by the hypotheses demanded by the theorems presented in [67,81,97,98]. It is therefore clear that the study of mathematical well-posedness of the equilibrium and dynamical problems for dilatational strain gradient elasticity require the adaptation and/or the generalization of the standard arguments used in the theory of elasticity, in a similar way to what has been done, for the linear case, in [45,48]. In these last papers, the considered problem was linear, while here one deals with a nonlinear one. Therefore, sophisticated arguments involving lower-semicontinuity seem necessary, and will be the object of future investigations. Instead in the present paper, in order to give a strong mechanical motivation to the needed mathematical efforts demanded by the novel energies considered, an approach having the same style as the one usually used in mechanical literature has been adopted, based on the techniques originally employed by Lagrange himself. Without any formal specification of the functional space where the introduced energy functional are introduced, their first variation have been formally calculated.  Fig. 6 whereas for c and d as in Fig. 7, respectively The used heuristic procedure will lead to a formal expression of the work functional for externally applied loads: it is proven that together with standard volume forces and surface contact forces, dilatational strain gradient continua can support, on edges already present in the reference configuration, externally applied edge forces. Moreover, see [63], one can apply on contact surfaces also external double forces, that it interactions that expend work on normal (to the contact surface) derivative of virtual displacement. Rather interestingly it has been presented a reasoning which indicates that one can apply, on the regular parts of the boundary (i.e on the boundary faces) of the considered dilatational strain gradient continua, external forces concentrated on curves that are not already edges, in the reference configuration. These linear forces may, in the actual configuration, give rise to new edges, edges that did not already exist in the reference configuration. One has to remark that nonlinear dilatational strain gradient elasticity is a very interesting particular case of the complete Toupin-Mindlin nonlinear strain gradient elasticity. Indeed, in the class of the continua considered here the only second gradient dependence of deformation energy involves the gradient of the volume change from the reference to the actual configuration: in other words, it has been considered those phenomena where the inhomogeneity of the dilation state produces an extra contribution to deformation energy, while the other eventual parts of second gradient of placement are "energy indifferent". One has to remark that dilatational second gradient continua must be related to other models introduced in generalized continuum theory. In particular, it is immediate to regard them as the continua obtained by considering the models introduced in the so-called poroelasticity (see, in this context, e.g. [112]). In poroelasticity, a scalar (one-dimensional) porosity field is introduced to describe the porevolume changes, which are distinct from the overall-continuum-volume changes. When a constraint linking the porosity kinematical descriptor to the determinant of placement gradient is introduced, then poroelasticity reduces to the presently discussed theory. A rigorous study of this relationship must consider the introduction of Lagrange multipliers and should parallel the considerations presented in [11].
On the other hand, the set of nonlinear dilatational strain gradient continua not only includes capillary fluids that have been extensively studied already, and that, in practice, are the only second gradient continua being considered in the literature, but also encompass a rather wide class of solids, whose second gradient behaviour is, in the sense specified before, singular. One could say that these dilatational strain gradient solids are a kind of "solidified" second gradient Korteweg or Cahn-Hilliard fluids. By exploiting the techniques developed by Lagrange and Piola, it has been managed to derive the expressions for virtual work equality (including, eventually, inertial contributions) in both Lagrangian and Eulerian description and to obtain the formulas for Piola's transport for the considered class of continua. Some relevant conditions characterizing the set of double forces and edge forces applicable to dilatational strain gradient are derived: in particular, it has been proven that contact edge forces cannot have, in Eulerian description, any component in the edge tangent direction. In order to check the logical consistency and the predictive capacity of considered continuum models, the linearized equations governing two important paradigmatic problems (involving wave propagation and deformation of an elastic tube under external pressure and double force) have been considered; the results were obtained should be regarded as very promising. This for two reasons: first for their intrinsic value, but also because they will represent a valid checking of future numerical codes to be developed for studying nonlinear problems, of interest, for instance, in the theory of metamaterials.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.
Funding Open access funding provided by Università degli Studi di Cagliari within the CRUI-CARE Agreement. This study was funded by the Russian Science Foundation (Grant No. 20-41-04404).

Declarations
Ethical approval This article does not contain any studies with human participants or animals performed by any of the authors.

Conflict of interest
Author VAE declares that he has no conflict of interest. Author AC declares that he has no conflict of interest. Author FdI declares that he has no conflict of interest.

Appendix: Lagrange variational principle for the nonlinear strain gradient elasticity considering edges
Within the Toupin-Mindlin strain gradient elasticity the Lagrange variational principle has the form δE TM − δA TM = 0, (138) where E TM is the functional of total energy and δA TM is the work of external actions given by On a part A f of the boundary A = ∂ V the external actions are prescribed. Here, ρ is the referential mass density, f is a vector of mass forces, t is a vector of the surface traction, c is a vector of surface external double forces, δx is the variation of x, ∂/∂ N is the normal derivative, t m and q are vectors of line forces given along L m and L, respectively. Note that δA TM does not constitute the first variation of any functional, in general. Such functional exists only for conservative loadings.
Introducing the Piola-Kirchhoff stress P and hyperstress M one gets the first variation of E TM in the form where v = δx. Then, integrating by parts and under standard assumptions of the calculus of variations from (138) we derive Eqs. (8), (11), (12), and (14), (15), see [8,43,44] for more details.
The same technique can be applied for deriving Eulerian form of the respective boundary-value problems. To this end, (140) is transformed into the Eulerian form changing coordinates X → x and using the identities dv = J dV and F · ∇ = ∇ δE TM = v (P : ∇v + M . . . ∇∇v) dv, where are the Cauchy stress and hyperstress tensors, respectively, and * stands for the Rayleigh product. Applying the integration by parts one gets The Eulerian form of δA TM is given by In (145) one has used the identities and da d A = J √ N · F −T · F −1 · N, ds dS = τ · F T · F · τ , τ = N × ν.