Internal resonance and bending analysis of thick visco-hyper-elastic arches

In this study, a comprehensive analysis of visco-hyper-elastic thick soft arches under an external time-independent as well as time-dependent loads is presented from bending and internal resonance phenomenon perspectives. Axial, transverse and rotation motions are considered for modelling the thick and soft arch in the framework of the Mooney–Rivlin and Kelvin–Voigt visco-hyper-elastic schemes and third-order shear deformable models. The arch is assumed to be incompressible and is modelled using von Kármán geometric nonlinearity in the strain–displacement relationship. Using a virtual work method, the bending equations are derived. For the vibration analysis, three, coupled, highly nonlinear equations of motions are obtained using force-moment balance method. The Newton–Raphson method together with the dynamic equilibrium technique is used for the bending and vibration analyses. A detailed study on the influence of having visco-hyper-elasticity and arch curvature in the frequency response of the system is given in detail, and the bending deformation due to the applied static load is presented. The influence of having thick, soft arches with different slenderness ratios is shown, and the forced vibration response is discussed. Moreover, internal resonance in the system is studied showing that the curvature term in the structure can lead to three-to-one internal resonances, showing a rich nonlinear frequency response. The results of this study are a step forward in studying the visco-hyper-elastic behaviour of biological structures and soft tissues.

et al. [7], it was shown that by using arch structures, the total harvested power was 200% times more power compared to straight structures. As another example, shallow arches have also been used and examined for microelectromechanical systems (MEMS) [8].
During the last few years, researchers have studied the mechanical behaviour of both flat structures (including beams and plates) [9][10][11][12][13][14][15] and curved structures (including shells and arches) [16][17][18][19][20][21][22][23] to be able to use them efficiently for specific purposes. Focusing on arch structures, their behaviour in different engineering conditions has been examined lately. To name a few, Yang et al. [24] reinforced elastic arch structures using graphene nanoplatelets and studied the linear vibration behaviour of graphene-reinforced arches. Different fibre distribution was modelled using Halpin-Tsai micromechanics model and Hamilton's principle showing that out-of-plane and in-plane linear frequencies of the elastic arch are highly affected by the fibre distribution and percentage. For dynamic instability analysis, Zhao et al. [25] studied the influence of having porosity in modelling graphene-reinforced elastic arches; it was shown that the best dynamic instability resistance can be seen in the symmetric porosity distribution model. The nonlinear mechanics of elastic shallow arches with nonuniform cross section made of functionally graded materials (FGMs) have been investigated by Ghayesh and Farokhi [26]. Coupled axial-transverse equations of motion were obtained using Euler-Bernoulli theory and Hamilton's principle; it was shown that as the arch curvature increases, the maximum amplitude of transverse vibration decreases. For MEMS structures, Farokhi et al. [27] studied the pull-in phenomena in micro-scale arch structures. Using a high-dimensional reduced-order Galerkin model and generalised Hamilton's principle, it was shown that by increasing the curvature of the arch, the DC voltage for pull-in phenomena increases.
The previous studies were mainly focused on linear-elastic structures under small strains. However, soft materials undergo large strains, which makes the classic linear modelling of the behaviour inaccurate for many hyper-elastic structures. There has been a significant attention on hyper-elastic structures in the past few years focusing on the statics [28][29][30][31] and dynamics [32][33][34][35]. In some recent studies, the anti-clastic bending of hyper-elastic beams in finite elasticity has been investigated by Lanzoni and Tarantino [36][37][38]. Khaniki et al. [39,40] have examined the nonlinear vibration and mechanics of axially moving hyper-elastic structures and layered hyper-elastic beams. In another study, Khaniki et al. [41] investigated the effect of having porosity in hyper-elastic structures providing a porous-hyper-elastic strain energy model and modelling the vibration behaviour of porous soft beams. A detailed review on the nonlinear dynamic behaviour of different hyper-elastic structures can be found in Ref. [42].
However, since some soft structures show a combination of viscoelasticity and hyper-elasticity, researchers have developed different models for such structures [43][44][45][46][47]. In a recent study presented by Li et al. [48], a wide range of visco-hyper-elastic constitutive models were developed. The viscoelasticity was modelled via Maxwell and Kelvin-Voigt models, and the hyper-elasticity was modelled via neo-Hookean, Mooney-Rivlin, and Ogden models. It was shown that the proposed models are capable of modelling both hyper-elasticity and viscoelasticity of the biological structure. Since visco-hyper-elastic structures are widely applicable for modelling human body organs, soft robots and prosthesis, the visco-hyper-elastic constitutive models are useful for studying and manufacturing accurate mechanical and biomechanical structures in different environments.
Besides, since most of the hyper-elastic structures are significantly softer than linear-elastic structures, it is very likely to have them bent (like arches and shells) in real-life applications. Accordingly, in this study, viscohyper-elastic thick soft shallow arches are examined in the framework of nonlinear bending and vibrations (with and without internal resonances). Since the soft arch is assumed to be thick, a third-order shear deformable theory is used and the visco-hyper-elasticity is modelled following Ref. [48] by using Mooney-Rivlin and Kelvin-Voigt models. The bending behaviour of the soft arch is formulated using the virtual work method and solved using the Newton-Raphson method. For the vibration analysis, the equations of motion are obtained using the force-moment balance method and solved using a dynamic equilibrium technique. The effect of having curvature in the structure is analysed, and the internal resonance phenomena caused by having this nonlinearity in the system are discussed. A detailed discussion of the influence of curvature term, the slenderness ratio and three-to-one internal resonance is given.

Bending of hyper-elastic thick arch formulation
For a thick, soft, shallow arch, using a higher-order shear deformation theory, by considering plane motion, von Kármán's geometric nonlinearity, and curvature in the thickness direction, the nonlinear static strain terms are obtained as ε 33 = 0 (should satisfy the incompressibility condition) where u 1 and u 3 are the displacements of the arch in the x 1 and x 3 directions (shown in Fig. 1a), u 30 is the curvature of the shallow arch, ε 11 and ε 33 are the axial strains in the x 1 and x 3 directions and ε 13 is the shear strain. By considering plane motion and incompressible material, the right Cauchy-Green strain tensor is written as [49] which contains nonlinear higher-order strain terms since hyper-elastic structures undergo large strains and cannot be modelled by considering the linear part alone, as is undertaken for elastic structures. The principle of virtual work is not only for linear elastic structures and can be applied for hyper-elastic structures [50]. To this end, by using the Mooney-Rivlin hyper-elastic strain energy density [51,52] and the principle of virtual work, the variation of the energy term, due to the static applied forces, can be written in three variation terms for axial (PE1), transverse (PE2) and rotation (PE3) variations as where PE1, PE2, PE3 are given in Appendix A for the sake of brevity, with the moment of area coefficients defined as Using Eqs. (1-7), three equilibrium equations for the bending analysis are obtained by having homogeneity in the arch as It can be seen that the equilibrium equations are highly nonlinear in hyper-elasticity terms, with a high coupling between the axial, rotation and transverse motions. The equilibrium equations can be validated by comparing them to those presented in Ref. [40] by neglecting the curvature terms from these equations and unsymmetrical and time-dependant terms in the equations of motion of Ref. [40]. By employing the solution procedure given in Sect. 4, the bending behaviour of soft thick shallow arches will be obtained.

Dynamics of visco-hyper-elastic thick arch formulation via force-moment balance method
For the time-dependent analysis and modelling of soft thick shallow arches, using the same assumptions given in Sect. 2, the nonlinear strain terms are written as ε 33 = 0 (should satisfy the incompressibility condition) .
For visco-hyper-elastic structures, Li et al. [48] have shown that the visco-hyper Kelvin-Voigt model has a good accuracy in modelling visco-hyper-elastic structures while using Mooney-Rivlin, neo-Hookean and Ogden strain energy density models. The visco-hyper-elastic model (Fig. 2) is presented as [48] σ visco-hyper = σ hyperelastic + ξ kε , (Talk about the feasible applicability of the approach) For an element of the thick visco-hyper-elastic arch, shown in Fig. 2b, Newton's law is applied for the rotational motion, as well as translational ones by having [53] ∂ N ( where the stress resultants are defined as [54] and the applied external periodic load with the magnitude of F is given as For isotropic visco-hyper-elastic shallow arches, by using force-moment balance method, one can reach the coupled equations of motions which are shown in Appendix B for the sake of brevity. By neglecting the viscoelasticity and the arch-related terms, the equations of motion can be validated through Ref. [40]. It can be seen that all three equations of motion are highly nonlinear in visco-elasticity and hyper-elasticity terms, with a high coupling between the axial, rotation and transverse motions. By using the solution procedure given in Sect. 5, the vibration behaviour of soft thick shallow arches will be obtained.

Solution procedure for hyper-elastic bending analysis
By assuming only transverse external load on the system (F s1 = F S2 = 0 and F S3 = F) and defining nondimensional terms as the nondimensional equilibrium equations are written as where * is neglected from the parameters for the sake of brevity. By employing a series expansion, the degrees of freedom are written as and the equations of motion are discretised as The linear stiffness coefficients (K L i j ) of the bending equilibrium equations are defined as and the nonlinear stiffness coefficients (K N L i j ) are defined as which by solving the nonlinear polynomial equations of motion using the Newton-Raphson method, the static bending due to the external static force can be obtained.

Solution procedure for visco-hyper-elastic vibrations
Using the given nondimensional terms in Sect. 4 and additional nondimensional terms as the nondimensional equations of motion are written as For time-dependant analysis, by employing the Galerkin's procedure, the degrees of freedom are written as and the equations of motion are discretised as with the mass coefficients (M i j ) as and the linear visco-elastic coefficients (C L i j ) as The nonlinear damping coefficients (C N L i j ) are shown in Appendix C for the sake of brevity, and the stiffness terms are defined in the previous sections. By solving the dynamic equilibrium equations [39,41], the nonlinear vibration response of the system is obtained.

Results and discussions for bending and vibration
In this section, the nonlinear bending and vibration responses of thick visco-hyper-elastic shallow arches are modelled for different cases of curvature, slenderness ratio and internal resonances. The hyper-elastic properties are taken from Ref [55] for vulcanised rubbers with the Mooney-Rivlin coefficients as C 1 = 0.28e6 Pa and C 2 = 0.15e6 Pa and mass density ρ = 950 kg/m 3 .

Study 1: Force-bending analysis
As the first part of the analysis, the nonlinear static deformation of thick hyper-elastic shallow arches is analysed following the solution procedure given in Sect. 4 and compared to those obtained using finite element software. The shallow arch is assumed to have the following geometrical properties: L = 0.3 m, h = 0.03 m, b = 0.03 m, with a simply supported boundary condition and assumed to undergo a static pressure (which is varied from − 1 to 1 kPa). The maximum transverse deformation is obtained and shown in Table 1 for different curvatures and external loads. By comparing the current results with those obtained from ANSYS [56], it is shown that the results are in very good agreement. Besides, it can be seen that by increasing the curvature term, the maximum deformation of the shallow arch decreases significantly especially for higher external loadings. Moreover, to understand the static behaviour of the shallow arch under different types of loadings, the curvature is assumed as following the first mode as W 0 (x 1 ) = R 0 W 1 (x 1 ) where R 0 is the nondimensional curvature term assumed as 0.6. The nondimensional external applied static force (F bending ) is modelled using a combination of sin functions as where n is varied from 1 to 70. For the given function, the static force through the length of the arch is shown in Fig. 3a for different values of n. By solving the nonlinear bending equilibrium equations, the nondimensional deformations of the hyper-elastic thick arch are shown in Fig. 3b-d. It can be seen that there is a high coupling between the degrees of freedom; also, the force type and magnitude can change the deformations in the structure, significantly. Study 2: Influence of the curvature on the bending behaviour For the same shallow arch with the given properties in the previous subsection, the influence of the curvature in resisting deformation in the structure is analysed. To do this, an external static force following Eq. (93) is given to the structure with f 1 and f 2 as 1000 and − 500, respectively. Figure 4a-c shows the axial, transverse and rotation deformation of the shallow arch by varying the curvature term as R 0 = 0.1, 0.3, 0.5, 0.7, 0.9, 1.0, and 1.5. It can be seen that increasing the shallow arch curvature leads to higher resistance against the transverse loading. Besides, as shown in Fig. 4a, the coupling between the axial and transverse deformations increases significantly by increasing the curvature term.

Study 3: Influence of considering the higher-order nonlinear terms
Since incompressible hyper-elastic structures show nonlinearity in their constitutive equation as well as large deformations, the importance of considering the higher-order nonlinear terms is discussed in this subsection. Accordingly, for the same shallow arch in previous sections with the curvature of R 0 = 1.5, bending behaviour is obtained using both linear and nonlinear models when subjected to an external static force following Eq. (93) with f 1 and f 2 as 10000 and − 2500. Figure 5a-c shows the axial, transverse and rotation deformation of the shallow arch by only considering the linear part of the equilibrium equations as well as considering the full equations. It can be seen that the linear simplification of the structure leads to underestimated bending deformations for the axial deformation and inaccurate transverse and rotation motions which indicates the importance of accurately modelling the nonlinear behaviour of the structure.

Study 1: Vibration in the hyper-elasticity framework with internal resonances
Since the curvature of the thick, soft shallow arch can change the linear frequencies of the structure, this section concentrates on detecting the internal resonance behaviour of the arch due to the curvature. For the same shallow arch model given in Sect. 6.1, the ratio of the second natural frequency to the fundamental natural frequency (ω 2 /ω 1 ) is shown in Fig. 6 by having different curvature terms. It can be seen that around R 0 = 0.3, the ratio is equal to 3 and the three-to-one internal resonance in the system occurs. Accordingly, by having R 0 = 0.3 (which leads to ω 2 /ω 1 ∼ =3), the frequency-amplitude responses of the system are shown in Figs. 7, 8 and 9 for the dominant vibration coordinates in the axial, transverse and rotational directions. It can be seen that a high coupling between modes it obtained, and the system shows a rich nonlinear behaviour, with a combination of hardening, softening and modal interactions.

Study 2: Vibration in the visco-hyper-elasticity framework with internal resonances
Many soft structures show visco-elastic behaviour together with hyper-elasticity (especially biological tissues). To see the influence of this combination on the nonlinear frequency response of the structure with an internal resonance, this section analyses the visco-hyper-elasticity behaviour of thick, soft, shallow arches. Since visco-elasticity is a function of many environmental parameters, such as temperature [57], this property of hyper-elastic structures should be measured for specific working environments. In this study, visco-elasticity is modelled using the Kelvin-Voigt visco-hyper-elastic model which has been shown that has a good accuracy in modelling visco-hyper-elastic materials [48]. The influence of the visco-elasticity on the three-to-one internal resonance of the structure is shown in Figs. 10, 11 and 12 for the axial, transverse and rotational directions by varying the viscoelastic parameter as ξ k = [2e−3, 2e−4, 2e−5]. It can be seen that the visco-hyper-elasticity has a significant effect in changing the nonlinear frequency response, and increasing this term from 2e−5 to 2e−4 decreases the maximum amplitude. By increasing the viscous term to 2e−3, the complicated nonlinear internal resonance behaviour will be fully damped.

Study 3: Curvature sensitivity
One of the main concentrations of this paper is on considering the curvature modelling in soft thick shallow arches. To show the importance and influence of this term on the nonlinear vibration behaviour of the structure, the shallow arch is assumed to have the following geometrical properties given in Sect. 6.1 with a simply supported boundary condition and curvature term varied as R 0 = [0.0, 0.1, 0.3, 0.5, 0.6, 0.8, 1.0]. The results for the amplitude-frequency response of the soft arch subjected to an external nondimensional load of F = 0.5 are shown in Figs. 13, 14 and 15 for the dominant vibration coordinates in the axial, transverse and rotation directions, respectively. It can be seen that increasing the curvature term from 0 to 1.0 moves the amplitude responses to higher frequencies. For higher curvature terms, transverse dynamic equilibriums show a softening effect, while for lower curvature term a combination of hardening and softening behaviour is obtained. Furthermore, by increasing the curvature term, the maximum amplitude of the first transverse dynamic equilibrium coordinate decreases significantly.

Study 4: Slenderness ratio sensitivity
Since a third-order shear deformable model has been used, different length-to-thickness models can be analysed. In this subsection, by assuming the same properties given in Sect. 6.1, the length of the arch is varied to L = S L h where S L is the slenderness ratio, which is varied as S L = [10,20,30,40,50]. The results for varying the length-to-thickness ratio are shown in Figs. 16, 17 and 18 for the axial, transverse and rotation dynamic equilibrium directions. It can be seen that increasing the slenderness ratio increases the stiffness hardening behaviour in the transverse coordinates. The first axial and rotation coordinates' amplitude decreases by increasing the slenderness ratio rate, which means the coupling between the axial and the transverse motion decreases and the rotation motion loses its importance.

Summary and conclusions
The nonlinear bending and vibration behaviours of simply supported, thick, soft, visco-hyper-elastic shallow arches are examined in this study. The visco-hyper-elasticity was modelled using an incompressible Mooney-Rivlin's hyper-elastic strain energy density model, together with Kelvin-Voigt visco-hyper-elasticity. A higherorder shear deformable model was used to model the axial, transverse and rotation motions. The curvature term was added to the model using von Kármán geometric nonlinearity, and the incompressibility condition was satisfied in the strain-displacement definition. The coupled nonlinear equations of motion were obtained and solved using force-moment balance method and the virtual work method for the vibration and bending analyses, respectively, showing that: (i) Increasing the shallow arch curvature leads to higher resistance against the transverse loading in all axial, transverse and rotation motions of the structure and changed the bending behaviour significantly. (ii) The linear simplification of the structure leads to underestimated bending deformations for all the axial, transverse and rotation motions which indicates the importance of accurately modelling the nonlinear behaviour of the structure. (iii) The coupling between the axial and transverse bending increases significantly by increasing the curvature term (iv) The curvature term can cause internal resonance in the system, leading to high coupling between modes, for which the system shows a rich nonlinear behaviour with a combination of hardening, softening and modal interactions. (v) The visco-hyper-elasticity has a significant effect in changing the nonlinear frequency response with a three-to-one internal resonance. By increasing the visco-hyper-elastic term to higher numbers, the complicated nonlinear internal resonance behaviour will be fully damped. (vi) For higher curvature terms, transverse dynamic equilibriums show a softening behaviour, while for lower curvature terms, a combination of hardening and softening behaviour is obtained. (vii) By increasing the curvature term, the maximum amplitude of the first transverse dynamic equilibrium coordinate decreases significantly. (viii) Increasing the slenderness ratio increases the stiffness hardening behaviour in the transverse coordinates.
(ix) The dominant axial and rotation coordinates' amplitude decreases by increasing the slenderness ratio rate, which means the coupling between the axial and the transverse motion decreases and the rotation motion loses its importance.
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 enabled and organized by CAUL and its Member Institutions Open Access funding enabled and organized by CAUL and its Member Institutions. No funding was received for this project.

Conflict of interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Appendix A
Variations of potential energy terms

Appendix B
For homogeneous visco-hyper-elastic shallow arches, due to the homogeneity of the visco-hyper-elastic arch, some of the inertia terms are equal to zero, simplifying the equations of motion to C N L 18 : +ξ k I 00