Pseudo-elastic cavitation model: part I—finite element analyses on thin silicone adhesives in façades

This study investigates the structural behavior of adhesive bonds of glass and metal using thin, structural silicones in heavily constrained applications. This special type of connection may lead to triaxial stress conditions under axial loading, which can lead to dilatation failure due to the abrupt growth of cavities (cavitation effect). Cavitation failure leads to significant stress softening and loss of stiffness; however, it increases connection’s ductility. These material deformations should be considered when designing glass-metal connections. Therefore, a constitutive model is developed to account for cavitation in hyperelastic materials. The volumetric component of the model is equipped with a non-linear Helmholtz free energy function that accounts for isotropic void growth under hydrostatic loading. An energy coupling term is then added that numerically explicates strain energy under isochoric deformation, while also guaranteeing physical material behavior. The energy contribution is calculated internally by analysing the geometric evolution of inherent voids. The extended volumetric–isochoric split enables one to numerically calculate heavily constrained silicone joints under arbitrary deformation modes. Three-dimensional finite element calculations on uniaxial tension, bulge, and pancake tests validate the constitutive model. All experiments could be validated with one set of material parameters through numerical simulations. The numerical calculations were robust and efficient without any underlying mesh dependencies.

contribution is calculated internally by analysing the geometric evolution of inherent voids. The extended volumetric-isochoric split enables one to numerically calculate heavily constrained silicone joints under arbitrary deformation modes. Three-dimensional finite element calculations on uniaxial tension, bulge, and pancake tests validate the constitutive model. All experiments could be validated with one set of material parameters through numerical simulations. The numerical calculations were robust and efficient without any underlying mesh dependencies.

Motivation and problem statement
When constructing conventional glass façades, a critical design element is the connection between the glass pane and the supporting secondary load bearing members. To avoid stress peaks in the glass, the joint must be ductile and be able to withstand environmental influences (wind, rain, sun, etc.), while also being architecturally appealing. In modern construction, one typical example of a ductile connection is achieved by directly bonding glass to metal with thin silicone adhesives (see Fig. 1). This is also known as a laminated connection due to a special production process (Bedon and Santarsiero 2018). When an adhesive joint is subjected to an axial load, the disability of lateral strains may cause abrupt cavity growth in the silicone. This behavior is known as cavitation (Drass et al. 2017) and under its effects the connection experiences significant increases in ductility before final dilation failure occurs (see Fig. 2).
At present, such thin adhesive joints may not be employed or studied because they do not meet the geometric requirements of the Guideline for European Technical Approval (ETAG 002 2012). The adhesive joint height to width ratio must be in the range of 0.33 < H/B < 1 to prevent cavitation failure. In contrast to thick silicone adhesive joints-which are typically used in façade constructions-thin adhesive joints offer clear advantages in terms of design, rigidity, struc- tural diversity and architectural aesthetics. So far there are no simple numerical models to describe the material behavior of thin silicone adhesives in façade applications, so analytical models were developed in order to realistically depict the structural behavior (Descamps et al. 2017). Ongoing research in the field of structural silicones and MS polymers lies in creating physically accurate, constitutive models of structural behavior (Dispersyn et al. 2017;Drass et al. 2018d) and limit state analyses of such connections (Staudt et al. 2018;Rosendahl et al. 2018;Santarsiero et al. 2018).
The attempt to approximate the special structural behavior with classical compressible hyperelastic material models using the volumetric-isochoric split in accordance to Flory (1961) leads to non-physical effects (Ehlers and Eipper 1998;Li et al. 2007). This fact was proved by the work of Danielsson et al. (2004) and Drass et al. (2017), who analysed multivoided representative volume elements (RVE) under hydrostatic tension loading. A key component is that the structural behavior exhibits a high initial stiffness followed by strong stress softening caused by excessive void growth. Except the constitutive models proposed by Danielsson et al. (2004), Li et al. (2007) and Drass et al. (2017), none of the classical volumetric Helmholtz free energy functions (Hartmann and Neff 2003) are able to represent this physical phenomenon (see "Appendix B"). Hence, conventional compressible hyperelastic material models are not suitable to represent volumetric stress softening effects caused by void growth. Even the Gurson model, which accounts for void growth in metals, cannot be applied to polymeric adhesive connections (Silvestru et al. 2018).  Focusing on point fixing connections in modern façades, classical finite element models accompanied with a compressible hyperelastic material formulation lead to an overestimation of stresses, thus, underestimating the ductility of this special type of connection. This fact is illustrated in Fig. 2, which presents the experimental and numerical results of an axially loaded point fixing. From Fig. 2, it is obvious that the apparent stress softening effects due to cavitation cannot be simulated utilizing conventional hyperelastic material models. Considering these drawbacks, the design process for thin adhesive glass to metal connections in façades design issues a challenge for structural and façade engineers.

State of the art
Considerable effort has been put forth analysing experimentally thin, structural silicone adhesives as a bonding material for point fixed connections (Hagl et al. 2012a, b;Drass et al. 2018b). Analytical prediction models based on extensive experimental investigations were proposed by Santarsiero et al. (2016) to describe the tensile resistance. Regarding the numerical modelling of such connections, especially the stress softening under constrained axial loading, only Drass et al. (2017Drass et al. ( , 2018b proposed a nonlinear volumetric Helmholtz free energy function and a continuum damage formulation accounting for the degradation of the bulk modulus due to void growth. A modified version of the model of Ogden and Roxburgh (1999) was utilized to calculate isochoric softening due to the Mullins effect (Ioannidou-Kati et al. 2018). Yet still, the proposed models lack either a robust and efficient numerical calculation providing mesh-independent results or are not able to account for stress softening due to cavitation.
Pioneering experimental work focusing on the cavitation effect in rubber-like materials was presented by Busse (1938), Yerzley (1939) and Gent and Lindley (1959). To obtain an overview of ongoing studies on the cavitation phenomenon in an analytical and numerical framework, the authors refer to Fond (2001). Dal et al. (2018) categorized the constitutive description of hyperelastic materials accounting for cavitation into three essential thematic categories. The first subsection splits the Helmholtz free energy function Ψ into an isochoric part Ψ iso and a volumetric part Ψ vol . Regarding volumetric Helmholtz free energy functions, there exist only a few proposals in literature to describe volumetric deformations under a physical viewpoint (Bischoff et al. 2001;Danielsson et al. 2004;Drass et al. 2018b). The main focus in the development of the volumetric Helmholtz free energy functions was on an incompressible and nearly incompressible material behaviour, which ensures the element's stability even under small isochoric stiffness (Bonet and Wood 2008, p. 171). The penalty function method describes a procedure where the strain energy potential is calculated as a sum of an incompressible material plus a penalty term enforcing the incompressibility constraint reading J = 1 (Simo and Taylor 1982). Following (Belytschko et al. 2013, p. 253), the penalty parameter K should be chosen large enough to enforce the incompressibility constraint, but not so large that numerical ill-conditioning occurs. An overview on the volumetric Helmholtz free energy function is given in Hartmann and Neff (2003).
Since the present paper utilizes the volumetricisochoric split to analyse cavitation, the second and third categories in accordance to Dal et al. (2018) are only briefly summarized. The second category analytically analyses the evolution of a single void as it reaches critical stresses in an undamaged medium. Additionally, the kinematics of a macroscopic continuum with an inherent void used as a representative volume element were examined by Ball (1982) and Hou and Abeyaratne (1992). The third category makes use of homogenization schemes accounting for an incompressible medium with a given void fraction. Danielsson et al. (2004) captured the idea of Hou and Abeyaratne (1992) to determine a homogenized Helmholtz free energy function for poro-hyperelastic materials. Li et al. (2007) developed a compressible hyperelastic material law accounting for cavitation and cavity growth using a Voigt-type homogenization.

Scope and outline
The present study utilizes the Flory-type volumetricisochoric split (Flory 1961) to develop a constitutive modelling approach accounting for stress softening in structural silicones under triaxial loading. Stress softening due to the Mullin's effect under isochoric deformation modes is excluded in this paper, but part of the investigations of Drass et al. (2019a). Section 2 provides a short overview on the pseudo-elasticity approach, which is used to account for modelling cavitation. In Sect. 3, the novel pseudo-elastic cavitation model is presented with a detailed derivation of the approach, including single-element tests as benchmark tests and the algorithmic setting. The new model provides an extension of the classical volumetricisochoric split by adding a coupling term Ψ vol,couple . With this approach, cavitation and the corresponding stress softening effect due to void growth in silicone adhesives can be analysed numerically for heavily constrained axially loaded connections. To validate the present constitutive model, Sect. 4 discusses three experiments from literature that are simulated with the novel constitutive model. Here, three-dimensional finite element calculations are performed simulating a simple uniaxial tension test (Santarsiero et al. 2016), a biaxial tension test in terms of a bulge test (Drass et al. 2018c), and a nearly triaxial tension test (which can be accomplished by testing the so-called pancake test) (Hagl et al. 2012a, b). In addition, pancake tests with two different diameters are simulated to prove the accuracy of the model.

General concept of pseudo-elasticity
In the investigations of Fung (1980) on living tissues, it was found that the loading and unloading paths under cyclic loading differ significantly. Each can be described separately by a unique relationship between stresses and strains. Since both load paths were represented independently by two single classical hyperelastic material models or with different material parameters for one hyperelastic material model, but the unloading path is strictly speaking inelastic, this phenomenological concept was termed pseudoelasticity. Lazopoulos and Ogden (1998) proposed the general concept of pseudo-elasticity by introducing a Helmholtz free energy function Ψ (F, η) that is dependent on the deformation gradient F and additionally on an internal scalar variable η. The extension of Ψ by a continuous or discontinuous internal variable η additionally motivates the designation of this theory as pseudo-elasticity.
For the extended functional Ψ (F, η), the first Piola-Kirchhoff stress tensor P can be calculated through assuming an elastic body, which is mapped from the reference to the current configuration by the deformation gradient F. Since the internal variable depends additionally on F, the functional must be differentiated according to η, whereby η must be subsequently differentiated with respect to F. In the absence of body forces the first and at the same time also usual equilibrium equation in the Lagrangian description reads where ∇− → X · (•) represents the divergence operator with respect to the reference configuration. Lazopoulos and Ogden (1998) formulated an additional equilibrium condition derived from a variational principle (stationary energy principle), which reads This equilibrium equation implicitly relates the internal variable η to F and thus continuously modifies the strain energy density function with progressing deformation. Based on these basic equations, Ogden and Roxburgh (1999) defined a pseudo-elastic energy function of type where φ (η) describes a damage function to describe the Mullin's effect. The inclusion of η and φ (η) offers the possibility to change the form of the energy function during deformation. According to Dorfmann and Ogden (2003), the damage parameter was assumed to be inactive during loading conditions and active for unloading with η ∈ [0, 1].
The internal variable was explicitly derived by inserting Eq. (4) in Eq. (3), which reads Since the damage function φ (η) is intended to determine η depending on the state of deformation, φ (η) can be chosen arbitrarily. However, Eq.

Pseudo-elastic cavitation model
A pseudo-elastic Helmholtz free energy function is derived in the following, which takes into account the effective stress softening caused by cavitation in rubber-like materials (Gent and Lindley 1959;Drass et al. 2018a). In contrast to the approach of Ogden and Roxburgh (1999), in which the theory of pseudoelasticity was used to calculate Mullin's damage in filled rubbers, here the contribution of the strain energy due to the anisotropic geometrical evolution of inherent cavities in the bulk material for any deformation is determined. The pseudo-elasticity approach in accordance to Ogden and Roxburgh (1999) was deliberately chosen to avoid multi-scale modelling (Dal et al. 2018), homogenization methods (Lopez-Pamies and Castañeda 2007) or even the FE 2 method (Feyel 2003) for reasons of simplicity and time saving due to the above-mentioned numerically complex calculation methods. The advantages are that the hyperelastic Helmholtz free energy function can be explicitly described at the macro-scale, the Flory-type volumetric-isochoric split can be used, which is ideal for the separate determination of material parameters, damage and healing effects can be modelled without great effort and the renouncement of time-consuming multi-scale modelling methods. Especially the description of damage and healing effects under cyclic loading can be described very clearly by this type of formulation. A detailed description on the extension of this approach towards cyclic loading is presented in part II of this paper.

General concept
To define a macroscopic Helmholtz free energy function that calculates the volumetric Helmholtz free energy of a porous rubber-like material under any deformation phenomenologically, the pseudo-elasticity approach according to Lazopoulos and Ogden (1998) is used. In principle, we assume a Flory-type hyperelastic material law with the special feature that the volumetric part of the functional is able to describe isotropic void growth. Drass et al. (2018b) presented a nonlinear volumetric Helmholtz free energy function (see Eq. (43)), which is capable of representing the cavitation effect and thus effective stress softening due to growing voids caused by triaxial loading (see "Appendix A"). Hence, the strain energy density function reads The advantages of this formulation lie in the application of the Flory-type volumetric-isochoric split, the possibility of representing stress softening due to cavitation and the adaptability of the volumetric formulation. The disadvantage are described in "Appendix B" for reasons of clarity. However, in summary the main disadvantage lies in the fact that under isochoric loads anomalies can be observed in the structural response that only arise from the application of the Flory-type split of the Helmholtz free energy function. Therefore, the pseudo-elastic formulation is chosen in order to avoid the disadvantage mentioned above. Nonetheless, to take the benefits of the Flory-type hyperelastic material and the nonlinear volumetric Helmholtz free energy function Ψ vol,ND , which is also capable of representing stress softening due to cavitation, an additional strain energy couple term Ψ couple is added. This couple term expresses the amount of the strain energy due to anisotropic geometric evolution of a void in relation to the dissipated energy due to isotropic void growth. Keeping this in mind, the modified Helmholtz free energy function reads To be more specific, the strain energy coupling term is defined by a shape function Ω ∈ [0, 1] based on the geometrical evolution of a deforming void. The variable Ω has therefore two main objectives. On the one hand, it governs the amount of dissipated strain energy D cav due to void growth, and on the other hand, it basically decides whether there is a volumetric or isochoric deformation. In the case of a purely volumetric deformation, Ψ vol,ND alone is able to calculate the strain energy as a result of void growth, so that the additional dissipated strain energy D cav due to cavitation must Fig. 3 Illustration of the total energy dissipation D cav due to void growth become zero. Consequently, Ω must be zero for a volumetric load. Looking at isochoric deformations, Ω must assume values of > 0 to compensate for the loss strain energy caused by pore growth, as no significant pore growth is expected under isochoric deformations. The dissipated strain energy D cav is defined as the difference between an unvoided and voided hyperelastic material under hydrostatic loading. For reasons of simplicity, it is assumed that the dissipated strain energy reads which serves as a working hypothesis. In this context, Ψ vol,classic represents the classical formulation of the volumetric part of a Helmholtz free energy function which is explicitly defined in "Appendix A" by Eq. (37). An exemplary illustration of the dissipated strain energy due to void growth is shown in Fig. 3. Since the coupling term depends on the global deformation gradient F, i.e. on isochoric and volumetric deformations, it can be understood as an intermediate configuration between the volumetric and isochoric configurations. To illustrate this, Fig. 4 shows the volumetric, isochoric and intermediate configuration and the associated Helmholtz free energy functions. Since the volumetric deformation gradient is a second order unit tensor weighted with a scalar value, the order of the decomposition of the total deformation gradient is arbitrary. To clarify that the computation of the Helmholtz free energies for the respective configuration is concerned, the undeformed and deformed material bodies (B 0 , B J ,B and B) including a deforming cavity for the respective configuration are shown.
To be in line with the theory of pseudo-elasticity, a function φ depending on Ω must also be added, which serves to implicitly calculate the shape function Ω in terms of the state of deformation. In contrast to the approach of Ogden and Roxburgh (1999), in which φ (Ω) was understood as a kind of damage function to describe the Mullin's effect under cyclic loading of filled rubbers, φ (Ω) represents the counterpart of a damage function in this case and thus serves only to determine the shape function Ω. This approach is consistent with the proposed ansatz of Lazopoulos and Ogden (1998), Ogden and Roxburgh (1999) and Dorfmann and Ogden (2003). Accordingly, the entire pseudo-elastic Helmholtz free energy function can be written as where now the strain energy couple term Ψ couple (J, Ω) is defined as the sum of Ω D cav and φ (Ω).
The calculation of the Cauchy stress tensor for a Flory-type hyperelastic material is generally done by For the stress calculation of the couple term, Ω must also be considered, since it depends on the macroscopic deformation gradient F or on the left Cauchy Green tensor b respectively. Considering only Ψ couple from Eq. (9), the Cauchy stress tensor is calculated as follows: Keeping in mind that Lazopoulos and Ogden (1998) determined an additional equilibrium equation based on equilibrium thermodynamics, which only results from the inclusion of an internal variable in the constitutive model, here Ω, the calculation of the Cauchy stress tensor for the couple term simplifies considerably (Dorfmann and Ogden 2003). Since the additional equilibrium equation must be fulfilled according to Lazopoulos and Ogden (1998), the Cauchy stress tensor is calculated after where p represents in general the hydrostatic pressure. Up to this point the shape function is still unknown, but the constraint equation Eq. (12) implicitly defines the variable Ω. Hence, the pseudo-elastic Helmholtz free energy function is now differentiated with respect to Ω resulting in This equation is used in the following to derive the shape function Ω. However, additional constraint equations for φ must be defined beforehand. An additional constraint equation according to Ogden and Roxburgh (1999), which must be fulfilled to determine the shape function Ω, reads If Ω is inactive, i.e. Ω = 0, then the Helmholtz free energy describing isotropic void growth according to Eq. (6) must come out. A third condition is described by in particular, intended to provide a softer structural response on the unloading path due to the Mullin's effect. Since the present material model initially ignores material softening due to cyclic loading, the condition of Eq. (16) is not further considered. A self-evident requirement, which must also be met, is that the shape function is within the limits of Ω ∈ [0, 1]. Fulfilling the constraint equations, the choice of φ (Ω) is arbitrary. In this case, φ (Ω) was chosen in such a way that φ (Ω) reads where Θ ∈ [0, 1] and Π ∈ [0, 1] are additional variables that distinguish between volumetric and isochoric deformations or describe the aspect ratio between a deformed to an undeformed cavity. The exact derivation of both variables is described in the following sections. Returning to the shape function, it can be explicitly described by inserting Eq. (17) into Eq. (14), which leads to Since Θ and Π are limited in their functional values, but the dissipated energy D cav can take any positive values depending on the applied deformation, the value range of the trigonometric function is additionally given with Looking at the right side of Eq. (18), the first term aims for the value one and is monotonically increasing under the condition of increasing isochoric deformation, whereas the second term becomes zero and monotonically decreases. This shows that Ω lies in the range of Ω ∈ [0, 1]. As the shape function is now fully described, the energy function φ (Ω) can be determined by integrating Eq. (17) with respect to Ω, which leads to On the basis of this procedure, the first condition that ∂Ψ couple /∂Ω must be zero is directly fulfilled, since this condition has just determined the shape function. Verifying the second constraint equation according to Eq. (15) by inserting Ω = 0 gives To meet this condition, i.e. to obtain a constitutive model that represents isotropic void growth when Ω is inactive, a new condition results for Θ, which governs the shape function Ω. The new condition to fulfil Eq. (15) is that the parameter Θ must be zero at purely hydrostatic load for which φ achieves the limit value zero. In contrast, it should take the value one for isochoric loads, so that this parameter can be described using a Heaviside step function, for example. Since the general concept of the Flory-type pseudo-elastic material model was presented, the parameters Ω and Θ of the shape function are defined more precisely in the following.

Isoperimetric inequality
As described in the general concept of the pseudoelastic cavitation model, there is a shape function that governs the strain energy due to the anisotropic geometric development of a deforming void in relation to the dissipated strain energy due to isotropic void growth. For this purpose, a scalar quantity must be found which defines the geometric evolution of a deforming cavity via a shape parameter within the limits zero and one. To do so, one of the fundamental problems in the classical calculus of variations is employed, which is to find the geometric figure with maximum area at a given perimeter (Zwierzynski 2016). The proposed problem can be described by formulating a so-called isoperimetric inequality defined in the Euclidean space R 2 . The isoperimetric inequality holds if a simple closed curve s of length L s , enclosing a region of the area A s , is a circle. Thus, ϑ describes the deviation between a deformed closed curve s and a circle. The classical formulation of the inequality ϑ is defined in R 2 with ϑ ∈ [0, 1]. To obtain a stricter criterion describing the deviation between a deformed closed curve s and a circle, the inequality can be provided with an exponent n, which reads then N represents the set of natural numbers. Since the present phenomenological pseudo-elastic approach is intended to model cavitation in rubbers and rubber-like materials, it is decisive to describe the geometrical evolution of an initially spherical cavity with one criterion. For this purpose, the isoperimetric inequality is converted into three-dimensional space in which ϑ is evaluated in three unique geometrical planes. Since the shape of the spherical cavity changes only due to an external deformation, it is obvious to describe the isoperimetric inequality as a function of a macroscopic deformation measure. Therefore the eigenvalues are used to evaluate the isoperimetric inequality. The geometrical planes, in which ϑ is evaluated, are uniquely defined by the cross product of the eigenvectors − → n of the left stretch tensor v, which reads − → 2,3]. Duplicate planes are excluded in the following considerations. For reasons of simplicity, a geometric mean value is calculated to evaluate the isoperimetric inequality for the three geometric planes (ϑ 1−2 , ϑ 2−3 and ϑ 3−1 ), i.e.
This simplification is permissible, since for the present pseudo-elastic Helmholtz free energy function, on the one hand, the strain energy contribution through anisotropic pore growth is considered purely phenomenologically and, on the other hand, a smeared isotropic material behavior is assumed. Both assumptions are merely working hypotheses that can be improved by integrating structural tensors to describe anisotropic material behavior.
As an example, the isoperimetric inequality is determined for the − → P 1−2 geometrical plane. Within this plane, the surface area of a cut void simply reads A s,1−2 = πλ 1 λ 2 , whereas the perimeter can only be calculated using an approximation. Utilizing the approximate solution in accordance to Ramanujan (Villarino 2005), the perimeter can be calculated with To obtain the values for A s and L s for all other principal stretch planes (1-2, 2-3 and 3-1), the indices of the principal axis have to be interchanged cyclically without considering duplicates. Consequently, the isoperimetric inequality can easily be calculated with respect to Eq. (22) depending on the desired principal stretch plane (1-2, 2-3 and 3-1).
To better understand the behavior of the isoperimetric inequality Π , a three-dimensional representation of it is shown in Fig. 5 depending on the exponent n. In addition, deformation paths according to uniaxial, biaxial and hydrostatic loading conditions are added in this plot to clearly show the dependency between Π and the applied deformation. It is obvious that the inequality assumes the value one for the reference configuration regardless of the deformation applied. Following the hydrostatic loading path, the inequality is always fulfilled and remains at the value one, even if the pore grows spherically. If the isochoric loading paths are followed, the value moves from one to zero for the inequality with increasing deformation. This behavior can be enhanced rapidly by increasing the exponent n.

Equivalent void growth measure
In this section, an equivalent void growth measure Θ is introduced that distinguishes between isochoric and volumetric deformations in terms of a scalar value. To differentiate between these deformations, a so-called void growth criterion ε eqn formulated in Hencky strain space according to Drass et al. (2018b) is used and mathematically modified, so that Θ takes the value zero for volumetric deformations, whereas it becomes one under isochoric loading. The requirement that the criterion equals zero for volumetric deformations results from Eq. (15), which means that the pseudo-elastic Helmholtz free energy function degenerates to Eq. (6). That is, the Helmholtz free energy function is able to represent stress softening due to cavitation under pure hydrostatic loading.
Returning to the void growth criterion, it is defined by where I 1,ε , I 2,ε and I 3,ε describe the first, second and third invariants of the Hencky strain tensor ε with As already mentioned, the equivalent void growth measure is a kind of Heaviside step function, which is generally defined by where [0, 1, c = 0] ∈ K. Keeping this in mind, the equivalent void growth measure can be fully described by This special formulation shows that Θ = 1 applies to isochoric deformations since ε eqn < 0 or takes a value of Θ = 0 for the reference configuration and for purely volumetric deformations since ε eqn > 0 .

Constitutive modelling and algorithmic settings
The calculation of the Cauchy stress tensor for the pseudo-elastic material model has already been presented in Sect. 3.1 and is summarized in Eq. (13). Due to the additional equilibrium equation resulting from the inclusion of an internal variable in the constitutive model according to Lazopoulos and Ogden (1998), the calculation of the stress tensor could be significantly simplified.
For reasons of completeness the spatial tangent stiffness will be additionally derived in the following. As already shown in the derivation of the Cauchy stress tensor, the components of the additional strain energy coupling term depend on J and b, but only the components of the derivative with respect to J are included in the final equation of the Cauchy stress tensor due to the additional equilibrium equation of ∂Ψ couple /∂Ω = 0. Therefore, only the derivation terms according to J are described below. According to the proposal of Holzapfel (2000, p. 265), the volumetric spatial tangent of the strain energy couple term c couple can be written directly in the Eulerian configuration by where I ⊗ I represents the dyadic product of two second order identity tensors and I describes a fourthorder identity tensor. Based on the simplifications due to the additional equilibrium equation of Eq. (12), the hydrostatic pressure of the strain energy couple term reads Since p couple is now defined by Eq. (33), the parameter s couple must be also computed. It describes the second order derivative of Ψ couple with respect to J , which reads Keeping in mind that Ω was already defined in Eq. (18), the calculation of s couple can be simplified to An algorithmic box for the numerical treatment is given in Table 1.

Parameter studies of isoperimetric extension of volumetric Helmholtz free energy function
To get more insight into the function of Ω or Θ and Π , and to understand how the algorithm behaves with any deformations, reference Fig. 6. Here Ω, Θ and Π are shown under different deformation modes with respect to the principal stretch λ 1 . In this context, λ 1 is displayed on a logarithmic scale to investigate the behavior of Ω or Θ and Π at small deformations. For this purpose, a single element test under uniaxial, biaxial and hydrostatic load is numerically investigated. Since Ω, Π and Θ depend only on the principal stretches λ i and an acceleration parameter n, the results displayed depend on the deformation λ 1 and n. As mentioned at the beginning, the parameter n controls how fast one moves from the initial value of ϑ = 1 of the isoperimetric inequality to the value zero due to an isochoric load. However, to obtain the classical isoperimetric inequality, n must be set to one.
To evaluate the described parameters Ω, Π and Θ numerically, a normalized Neo-Hooke material model with μ = 1 MPa is used. The parameters of Ψ vol,ND are selected as examples to perform the numerical calculations of the single-element benchmark tests (κ 0 = 0.02181, κ 1 = 0.67900 and κ 2 = 0.09791).
As shown in Fig. 6, Θ starts at Θ = 0 and very quickly reaches Θ = 1 for isochoric deformations (UT, BT). In contrast, for triaxial loading, Θ = 0 with the result that only Ψ vol,ND is active during numerical calculation. It should also be noted that the equivalent void growth measure Θ starts at a value of Θ = 0 for zero deformations. The structural reaction of the isoperimetric inequality Π starts at zero deformations (λ 1 = λ 2 = λ 3 = 1) with a value of Π = 1. For the homogeneous hydrostatic load state, Π = 1 for the entire deformation regime, while for isochoric deformations the isoperimetric inequality decreases from one to zero depending on the parameter n. By increasing the parameter n, the proposed constitutive model converges to the classical formulation of Ψ vol,classic much faster, as the influence of Ψ vol,ND is reduced significantly.
In summary, in this section the volumetric Helmholtz free energy function proposed by Drass et al. (2018b) was successfully transferred into the context of pseudoelasticity in order to prevent anomalies in structural behaviour, especially in isochoric deformation (see "Appendix B"). Therefore, the volumetric component of the model was equipped with the nonlinear Helmholtz free energy function that accounts for isotropic void growth under hydrostatic loading. An energy coupling term was then added that numerically explicates strain energy under isochoric deformation, while also guaranteeing physical material behaviour. The energy contribution is calculated internally by analysing the geometric evolution of inherent voids. On the basis of one-element tests at any deformation, parameter studies were carried out to confirm the physical behaviour and to analyse the influence of individual parameters.

Numerical validation of constitutive model
This chapter aims at the validation of the proposed material model. First of all, the quasi-static uniaxial tensile test and the bulge test are investigated with the pseudo-elastic cavitation model presented in Sect. 3. To show the general validity of this model, pancake tests with three different diameters are numerically calculated and compared with the experimental ones. Since, according to literature (Dal et al. 2018), it has not yet been possible to validate cavitation material models Table 1 An algorithmic box of the FE procedure during numerical simulations for the pseudo-elastic cavitation model 1.
compute shape function Ω based on (2) and (3) 5. compute strain energy couple term Ψ couple (J, Ω) and the resulting pseudo-elastic Helmholtz free energy function 6. compute Cauchy stress tensor σ at t n+1 compute tangent moduli c iso and c vol at t n+1 c = c iso + c vol with

8.
compute tangent moduli c couple = p couple + J s couple I ⊗ I − 2 p couple I, with over full-scale small component simulations, this is presented here for the first time.

Transparent structural silicone adhesive
The product DOWSIL TM TSSA, a transparent structural silicone adhesive, promises however due to its material properties, its high resistance to external environmental influences and its transparent appearance a higher application of silicones in this industry (Drass et al. 2019b). An interesting behavior can be observed when this material is exposed to isochoric or hydrostatic loading conditions. Under such loading conditions, the initially transparent silicone shows an increasingly opaque or white appearance, which usually is taken as an indicator of material softening.
The micro-structure of TSSA consists of polydimethylsiloxanes (PDMS) and a high amount of nanosilica particles (ca. 20-30%). The nano-silica particles form aggregates, where one nano-silica particle has a diameter of approximately 1 nm (Drass et al. 2019b). At the boundaries of these aggregates a multitude of nanocavities can develop due to the manufacturing process, which may lead to an increased porosity or also called finite porosity. In addition, the initially constant porosity can increase during deformation, when vacuoles arise at the boundaries of the nano-silica surfaces during stretching due to weak or non-existing fillermatrix interaction. Especially in flat-bonded transparent silicone adhesive joints under axial load, the socalled cavitation effect is present, in which a very dominant effective stress softening can be detected due to strong void growth under hydrostatic tensile load. In Fig. 6 Evaluation of Ω, Π and Θ as function of n dependent on the applied deformation experimental investigations on the cavitation effect on TSSA, which are usually carried out in so-called pancake tests, the effect of the white coloring was detected and analyzed by Drass et al. (2018b). In the pancake test, flat-bonded steel cylinders are pulled axially. As a result of the hindered lateral contraction, a triaxial stress state is generated in the material, which leads to excessive void growth.

Simulation of dumbbell-shaped tensile test
Dumbbell-shaped tensile tests, also known as uniaxial tensile tests, are the most common tests when analyzing the structural behavior of any materials, since they are relatively simple to perform. To proof that the pseudoelastic cavitation model of Sect. 3 is able to approximate this behavior correctly-not leading into a nonphysical behavior due to disproportional void growth-a full-scale simulation of the dumbbell-shaped tests presented by Drass et al. (2018c) is performed.
Since the present constitutive model only considers isothermal processes at a constant strain rate of , this study investigates the stress-strain behavior under a constant temperature of T = 23 • C and a quasistatic displacement rate v UT = 5 mm/s. Based on the results of Drass et al. (2018c), TSSA exhibits a nonlinear elastic behavior until fracture occurs-a characteristic property for rubber-like materials-however, TSSA does not exhibit strain stiffening due to locking stretches of the macro-molecules in the polymer matrix. Consequently, when considering isochoric deformations, simple isochoric constitutive models, like the Neo-Hookean model, can be used for TSSA (Drass et al. 2018d). The parameters of the pseudo-elastic cavitation model were obtained using inverse numerical optimization methods. Since the focus of this work is not on determining the parameters, reference is made to the work of Drass et al. (2017) to determine them. The material parameters used to characterize TSSA are given separately in Table 2 for the isochoric and volumetric Helmholtz free energies. The parameters have been determined simultaneously based on the uniaxial, biaxial and pancake tests. Since the meaning of the shear modulus is simple, the fully metric material parameters of Ψ vol,ND are described vividly in "Appendix C" in order to better understand and illustrate their meaning.
Based on the determined material parameters for the pseudo-elastic cavitation model, the critical hydrostatic pressure p cr at which cavitation ensues is reached for a critical relative volume of J cr = 1.041. During the uniaxial tension test the deformations are mainly volumeconstant because the voids drastically change shape but not in volume. This means the numerical calculation of the test should not give way to substantial void growth.
This can be confirmed regarding Fig. 7, where the results of the numerical calculation are visually compared with the experimental results. Figure 7 a compares the structural response of both the experiment and the finite element calculations. The min. / max. values of the experiments are gray curves, the mean values are red circles and the standard deviation (1 x SD) is an error bar. All numerical calculations were performed with an increasing amount of finite elements (328-23,142 Finite Elements) to analyze a potential mesh dependency. The results show no mesh-dependency for the simulations of the uniaxial tensile tests.
The numerical model in Fig. 7b displays the boundary conditions at a constant mesh density. Since threedimensional numerical calculations were performed, three-dimensional, higher-order volume elements with purely quadratic displacement behaviors were used. A solid element consists of 20 nodes with three degrees of translational freedom. In order to shorten the calculation time, three symmetry axes were used to analyze one eighth of a sample. A displacement of u = 75.00 mm was applied incrementally to the top of the numerical model. Thus, the simulation was carried out displacement-controlled. The numerical results of the pseudo-elastic cavitation model and the experimental structural response calculations correspond well. Since the constitutive approach takes void growth and the geometric evolution of cavities into account, the cavity evolution at various positions and load-steps can be seen in Fig. 7c. The geometrical evolution of the cavities is calculated by evaluating the isoperimetric, volumetric shape function Ω (Π, Θ). Based on numerically calculated values of the equivalent void growth measure Θ and the isoperimetric inequalities t 1−2 , t 2−3 and t 3−1 , the geometric evolution of an initially spherical cavity is calculated and visualized. Figure7c illustrates how the initially spherical cavities deform with minimal change in volume: they actually morph into a "cigar-like" shape. This is confirmed by evaluating the actual void fraction f in the current configuration. The void fraction f can be calculated by The relative volume J is calculated at each Gaussian point and at each load step. TSSA has an approximate initial pore content of approximately f 0 ≈ 3.0 % (Drass et al. 2019b). Figure 7 demonstrates that when a compressible hyperelastic material is loaded isochoric, there is a slight increase of void fraction, especially at large strains. This is natural, since the bulk modulus is not infinite for the analyzed material. Stress softening due to excessive void growth is not observed in the experiment or in the simulation. Consequently, the presented constitutive modeling approach, where an isochoric Neo-Hookean material model was coupled with the pseudo-elastic cavitation model from Sect. 3, is able to validate a simple tensile test without leading to nonphysical results.

Simulation of bulge-test
So-called bulge tests on TSSA were presented by Drass et al. (2018c), where the foil-like TSSA was inflated in a balloon shape using a special test device. In order to measure local strains, a stochastic speckle pattern was applied to the material to be tested. By this experimental set-up, a strain state near equi-biaxial tension can be examined.
The numerical model of the bulge tests was build up according to the test set-up. Rotational symmetries in the test set-up were employed to reduce the calculation time (see Fig. 8b). In accordance with Sect. 4.2, higher order, three-dimensional solid elements were used to calculate the three-dimensional numerical model of the bulge test. The anchor plate and the edge of the TSSA were fixed supported. A pressure p was applied to the bottom of the TSSA film to inflate it like a balloon. The finite element calculation was accomplished forcecontrolled. As a result of the force-controlled simulation of the bulge test, an instability problem occurs at very large deformations, so that the arc length method must be applied in the simulation. With increasing deformation, the TSSA film leans against the anchor plate so that a contact formulation must be used in the numerical model. The normal Lagrangian detection method was utilized for contact modeling The contact formulation defines a surface contact that can only be open or closed, and the normal Lagrangian formulation closes all gaps and eliminates penetration between the contact and target surfaces. The investigation of the effects of friction determined that the coefficient of friction has no influence on the force-displacement behavior. Therefore, the contact algorithm is formulated with a frictionless approach. Figure 8a compares the structural response between the bulge test and the finite element analyses. Here it becomes clear that the structural responses of the numerical models lie within the standard deviation of the experiment. This proves that the novel constitutive model approach is ideally suited to take biaxial deformations into account. The increasing mesh density shows no significant deviation between the structural responses of the finite element calculations. Figure 8c references the cavity evolution during biaxial deformation. According to the results of the uniaxial tensile tests (see Fig. 7c), the initially spherical cavities deform under almost volume constant conditions. This is obvious, since the bulge test is supposed to represent a nearly isochoric deformation. Softening in the experimental and numerical structural responses is a typical phenomenon for the bulge test. However, the observed softening has nothing to do with stresssoftening due to cavitation. This can be explained by the void volume's lack of change, even at large deformations. Due to the isochoric, biaxial deformation during the bulge test, the spherical cavities deform into a penny shape. As already mentioned, stress softening due to cavitation can be excluded when analysing the bulge test. Considering Fig. 8c the actual void volume is shown at different positions and at different load steps. Here, a slight increase in void volume can be seen at very large deformations, but this is primarily due to the compressible material formulation and not to cavitation.

Simulation of pancake tests
In this section numerical simulations of pancake tests from literature (Hagl et al. 2012a, b;Drass et al. 2018b) are carried out. Starting with the simulation of the pancake tests with a diameter of 20 mm based on the investigations by Hagl et al. (2012a, b), two flat-bonded steel cylinders were loaded in axial tension. Due to the lack of lateral contraction in the experimental set-up and the almost incompressible material behavior, large volumetric deformations occur in the adhesive layer. These conditions lead to the cavitation effect, and consequently to severe stress softening in the structural behavior.
The aim of this study is to numerically reproduce the structural response of the pancake test including the pronounced cavity growth that leads to a strong softening of the material. Symmetry boundary conditions were used to numerically model an angle ratio ϕ between a full cylindrical model and a specific section (see Fig. 9b). The lower surface of the numerical model is only supported in the z-direction. Because of this, r and ϕ are free in their deformation. The radial and tangential degrees of freedom (r and ϕ) were fixed for the upper side of the numerical model, while an axial displacement corresponding to the deformation of the pancake test was applied in the vertical z-direction. As a result, the finite element calculation was controlled by a continuous increase of the displacement. In accordance with the three-dimensional models for the uniaxial tensile test and the bulge test, the pancake test was also modeled three-dimensionally with Fig. 8 Validation of biaxial tension tests and schematic visualization of cavity evolution three-dimensional volume elements. Since the structural response of the numerical model and the experimental results agree, the novel constitutive modeling approach can adequately represent the strongly nonlinear structural behavior of the pancake test (see Fig. 9a). Furthermore, no mesh-dependent results are observed in the numerical analyses. Due to the initially almost incompressible material behavior, void growth must be taken into account for the numerical analysis of the pan-cake test. Cavity evolution is easy to prove: The initially spherical cavities begin to grow in both axial and lateral directions, and a closer look at Fig. 9c reveals that the void volume inside the adhesive increases significantly because of the constrained axial loading. At the edge of the adhesive, where the disability of lateral contraction is not very pronounced, the void volume increases only slightly. This demonstrates that there is a clear softening inside of the material due to excessive void Fig. 9 Validation of pancake tests performed by Hagl et al. (2012a, b) and schematic visualization of cavity evolution growth, while the marginal areas are not influenced by the cavitation effect. Therefore, the deformation of the pancake test is not comparable with the almost volumepreserving deformation of uniaxial tensile tests and the biaxial tensile test.
Since it has already been shown that the pseudoelastic cavitation model is well suited to approximate the pancake tests performed by Hagl et al. (2012a, b), additional pancake tests with other diameters are sim-ulated in the following. All simulation models are provided with the same set of material parameters, which are summarized in Table 2. The simulation results for three different diameters are shown in Fig. 10a. It is obvious that the pseudo-elastic cavitation model is very well suited to approximate the behavior of TSSA in the pancake test for different diameters. The model is able to represent the increase in initial stiffness for larger diameters of the pancake test specimen. Furthermore, Fig. 10 Validation of pancake tests performed by Hagl et al. (2012a, b) and Drass et al. (2018b) with (a) pseudo-elastic cavitation model and (b) classical volumetric Helmholtz free energy function Ψ vol,classic the occurrence of the cavitation effect and the associated strong material softening is well approximated. With increasing deformation, the load can be increased, which can also be modeled by the new model.
In order to further show that the classical formulation of the volumetric Helmholtz free energy function Ψ vol,classic , which is generally implemented by default in commercial FE codes, is not sufficient to represent the material softening in the pancake test, these were simulated again with the classical formulation (Fig. 10b). Here, only the initial stiffness can be approximated, which is not sufficient to represent the stress softening and accordingly the ductility of the bond.

Conclusions
This paper analyses cavitation in rubber-like materials, commonly understood as sudden cavity growth under hydrostatic stress. The experimental, analytical and numerical treatment of this phenomenon is an integral part of materials and mechanics research. However, many numerical approaches lack a user-friendly environment in terms of finite element calculations and numerical robustness. This paper presented a novel constitutive modeling approach for cavitation in rubbers and rubber-like materials, in which the volumetric part is equipped with a strongly non-linear Helmholtzfree energy function (Drass et al. 2018b), which is responsible for isotropic void growth under hydrostatic loading conditions.
The classic formulation of Ψ vol,ND was extended to a generalized form (see Eq. 43) so that the set of volumet-ric material parameters can be controlled by the running index n. Furthermore, an extension of the classic volumetric-isochoric split was carried out, where the addition of a coupling term Ψ vol,couple makes it possible to consider the cavitation effect in numerical calculations using a simple approach. The coupling term creates the consideration of cavitation in the form of stress softening under volumetric loads. However, if isochoric deformations occur, such as uniaxial or biaxial tension, the coupling term prevents extreme pore growth and consequently material softening due to pore growth.
The key element of this constitutive modeling approach is the geometric analysis of the development of initially spherical cavities under any deformations. For this purpose, Ω was defined as an isoperimetric volumetric form function that takes into account whether a void grows and how the void grows geometrically. Since Ω must explain if and how a void grows, it was constructed as a function of equivalent volumetric strain and isoperimetric inequality.
To validate the novel constitutive model approach, three experiments from the literature were simulated with the novel model. Three-dimensional finite element calculations were then performed on a simple uniaxial tensile test and a biaxial tensile test (in the form of a bulge test) Drass et al. (2018c) to validate the material law for isochoric deformations. The numerical analyses presented that the material behavior, as well as the evolution of inherent voids, could be well approximated. Furthermore, it was shown that the new material model does not soften the material due to void growth in isochoric deformations, in accordance with physical behavior. In contrast to the isochoric tests, the so-called pancake test was also modeled numerically, in which two flat-bonded steel cylinders are pulled axially (Hagl et al. 2012a, b). Due to the almost incompressible material behavior in the heavily constrained tensile test, large volumetric deformation occurs due to void growth, which leads to an enormous material softening. The experimental structural response as well as the geometry evolution of inherent pores could be simulated very well with the novel material model. Furthermore, the new material model made it possible to calculate material softening, which was previously not possible with commercial FE codes.
Future work will deal with an extension of the developed model with respect to the viscoelastic and temperature-dependent behavior in that direction that the models can be applied for structure calculation of thin silicone adhesive joints in high temperature environments such as the Middle East.
Acknowledgements Open Access funding provided by Projekt DEAL. We would like to gratefully thank The Dow Chemical Company ("DOW") for their support during our studies. Additionally, we would like to thank Julia Schultz and Joseph Yost for the help with the linguistic fine-tuning of the paper.
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/.

A Volumetric hyperelastic constitutive models
In contrast to the brief description of isochoric hyperelastic models and the restriction of the detailed description to only three essential models due to the enormous number of isochoric models, almost all volumetric energy functions occurring in the literature are presented in detail in the following, since only few approaches exist here. An overview of the various volumetric Helmholtz free energy functions for hyperelastic materials is given by Hartmann and Neff (2003). The simplest form of a volumetric Helmholtz free energy is dependent on the bulk modulus K and the relative volume J only. Due to its simple form, the functional is denoted as the classical volumetric Helmholtz free energy function, which reads Since the classical volumetric Helmholtz free energy function violates the postulate of convexity for J → 0, which reads Pence and Gou (2015) suggested a volumetric strain energy function that does not exhibit this disadvantage resulting in Simo and Taylor (1982) also proposed an improvement of the classical volumetric Helmholtz free energy function by adding the natural logarithm to the power of two. The enhanced volumetric strain energy function reads Regarding variational methods for the volume-constraint in finite deformation elasto-plasticity, Simo et al. (1985) proposed a variation of Eq. (40) by considering the natural logarithm only, i.e. Miehe (1994) proposed also an alteration of the abovementioned volumetric Helmholtz free energies, which reads Following Hartmann and Neff (2003) and Li et al. (2007), most of the models were developed without a physical background. The main focus was set on enforcing the incompressibility constraint under ensuring the element stability (Bonet and Wood 2008, p. 171). In contrast, Drass et al. (2018b) proposed a micro-mechanically motivated volumetric Helmholtz free energy function accounting for cavitation in structural silicones. The novel approach reads in a generalized form A represents an integration constant ensuring a stressfree reference configuration. In Fig. 11, the hydrostatic pressure p depending on the relative volume J is illustrated with respect to all proposed volumetric Helmholtz free energy functions. In addition, Ψ vol,ND is presented as a fitting result of all conventional volumetric strain energy functions to demonstrate its flexibility. From Fig. 11 it is obvious that Ψ vol,ND includes all classical energy functions, since Ψ vol,ND adequately approximates the structural responses of all classical models. Regarding Fig. 11, all proposed models show a similar structural response. The classic Helmholtz free energy function shows a linear dependence between p and J , while all other strain energy functions show a moderate to a strong nonlinearity. Regarding hydrostatic tensile loading (J > 1), moderate softening of structural behavior for Ψ vol,ext , Ψ vol,ST and Ψ vol,M can be observed. In contrast, Ψ vol,STP shows a distinct softening of the structural behavior, which is not due to the cavity growth under hydrostatic loading and the associated macroscopic material softening, but only due to the chosen mathematical approach. In summary, the classical volumetric approaches are good in the sense of a stable finite element calculation, but processes such as stress softening due to the cavitation effect in rubbers under hydrostatic tensile loading cannot be reproduced. However, with the micro-mechanical approach of Drass et al. (2017Drass et al. ( , 2018b void growth due to cavitation can be modeled.

B Remarks on volumetric Helmholtz free energy functions
In order to study the behavior of the proposed volumetric strain energy density function not only in pure hydrostatic load cases, simple tension problems using a Flory-type hyperelastic material model are analyzed in the following. The compressible Neo-Hookean material model combined with the novel nonlinear volumetric Helmholtz free energy function reads Dealing with a compressible Flory-type hyperelastic material formulation, anomalies within the structural behavior may occur due to the volumetric-isochoric split of an arbitrary Helmholtz free energy function (Ehlers and Eipper 1998). Since the proposed nonlinear volumetric Helmholtz free energy function accounts for stress softening due to the cavitation effect, it is logical that it may lead to an abnormal structural performance under isochoric deformation modes.
To do so, simple tension problems for a nearly incompressible behavior (κ 0 → 0) and highly compressible material behavior have been analyzed, where classical formulations of the volumetric Helmholtz free energy function were applied. In contrast, in this study we analyze the novel, strongly nonlinear volumetric strain energy density function, where the initial bulk modulus is varied between 10 ≤ K ≤ 1000 MPa. Since the present study does not investigate the influence of the bifurcation point for the onset of cavitation and softening behavior, these parameters are set constant. Thus, the bifurcation point is governed by κ 1 = 0.7 and the stress softening parameters read κ 2 = 0.3 and κ 3 = 0.
In Fig. 12 the results for the simple tension problem are illustrated under different viewpoints. Figure 12a shows the first and second principal stretch of a unit cell under uniaxial tension. Keeping in mind that the second principal stretch λ 2 must be calculated as implicit function of the first principal stretch λ 1 , an anomaly in the structural behavior is obvious. With an increase of the first principal stretch λ 1 , a decrease of the lateral stretch can be regarded. However, reaching a critical deformation level, the lateral stretch increases with an increasing longitudinal deformation. This means that the specimens becomes thicker under uniaxial tension deformation, which is not physical. This can be also observed regarding Fig. 12b, where a strong increase of the relative volume can be regarded for large longitudinal stretches. The strong increase of relative volume leads to stress softening effects under a uniax-ial deformation. Considering Fig. 12c, where the axial Cauchy stress is plotted against the Hencky strain, the structural responses show a physical behavior until a critical strain is reached. After reaching the critical point, a decrease of stress can be observed for all analyzed case studies. This can be also confirmed with regard to Fig. 12d, where a distinct stress softening can be observed for the Cauchy stress plotted against the Fig. 12 Parameter studies on novel volumetric Helmholtz free energy function in accordance to Eq. (44) analyzing a unit cell under uniaxial tension loading: a λ 1 − λ 2 , b λ 1 − J , c σ − ε and d σ − (J − 1) Fig. 13 Presentation of the experimental results of the pancake tests in the engineering stress-strain diagram, annotation of the meaning of the material parameters in the diagram and approximate equations for the volumetric parameters of the new volumetric Helmholtz free energy function dilatation J − 1. Since this behavior shows evidently an anomaly in the structural performance in accordance to Ehlers and Eipper (1998), the proposed Helmholtz free energy function must be enhanced to exclude such softening effect under isochoric deformations. For pure volumetric deformations, however, the new volumetric Helmholtz free energy function is well suited to represent a stiff structural response followed by excessive stress softening. Even stress softening with subsequent stress stiffening can be easily adjusted with the newly proposed volumetric Helmholtz-free energy function.
To avoid this non-physical mechanical behavior, Drass et al. (2018b) presented a so-called void growth criterion. The void growth criterion differentiates between isochoric and volumetric deformations based on case differentiation an "if-request" during finite element calculations. Another possibility is to formulate a condition for the start of volumetric softening. This approach was followed by Drass et al. (2017Drass et al. ( , 2018a in terms of a continuum damage formulation. However, it is advantageous to avoid case differentiation in the algorithm in order to obtain an efficient and robust numerical calculation without obtaining meshdependent results (Waffenschmidt et al. 2014). This procedure is applied and presented in detail in Sect. 3, where a continuous pseudo-elastic formulation of the non-linear volumetric Helmholtz free energy function is formulated without an if-request in the algorithmic setting.

C Definition of volumetric material parameters for Ψ vol,ND
In this section, the volumetric material parameters of Ψ vol,ND are vividly explained based on the structural response of a pancake test under axial loading. Converting the force-displacement curves of the pancake tests into engineering stress strain curves results in the following diagram. As can be seen from Fig. 13, the initial stiffness or the volumetric material parameter κ 0 respectively can be approximately calculated by the reciprocal of the initial engineering stress σ 0 divided by the initial strain ε 0 . This is indicated in the diagram by the gradient at the origin. The second parameter describes the bifurcation point or bifurcation stress p cr in the stress-strain diagram. From this the second volumetric material parameter κ 1 can be read approximately by forming the reciprocal of the bifurcation stress p cr . Since the last two volumetric material parameters are related to each other, they cannot be read directly from the engineering stress-strain diagram of the pancake tests. However, for such pancake diagrams the material parameter κ 2 must be within the limits κ 2 ∈ [0, 1] and the parameter κ 3 within the limits κ 3 ∈ [−1, 0]. At this point it should be noted that the material parameter κ 2 describing the inclination after the occurrence of cavitation failure reflects the effective stress softening due to excessive void growth, from which it is positive according to the definition of Nelder (1966). In contrast, the material parameter κ 3 describes the effective stress stiffening behaviour after cavitation has occurred, so that this is given a negative sign.