Stability analysis of a SAIR epidemic model on scale-free community networks

1.   Introduction

Mechanical metamaterials are engineered systems designed to create exotic materials, often prepared in laboratories using two or more components or phases. These materials exhibit unconventional mechanical properties not found in conventional materials ^[1]. Some notable examples include bandgap metamaterials composed of phononic crystals ^[2,3], which can be designed to exhibit adjustable frequency intervals where the transmission of mechanical waves is not permitted ^[4,5]. Tensegrity bandgap metamaterials are formed by alternating lumped masses with tensegrity units of different types, allowing for the creation of mass–spring systems that exhibit tunable frequency bandgaps. These bandgaps can be adjusted by internal and external prestress, as well as by contrasts in mass and spring properties between their elements (refer, e.g., to ^[6] and references tharein).

In the case of a forced mass–spring chain, when the length of such a system is finite and/or when the presence of a boundary breaks the periodicity condition, it is well known that the Bloch–Floquet wave form solutions of the dynamic problem do not longer satisfy the equations of motion. Two alternative strategies are available to overcome such an issue. The first one consists of solving the ordinary differential equations of motion of the discrete system with aid of the competent initial conditions, and accounting for any kinematic restrictions ^[6]. Such a strategy is generally accurate, but has a computational disadvantage if the number of unit cells is very large. The second strategy consists of deriving a continuum model of the mass-spring system equipped with competent Partial Differential Equations (PDEs) and Boundary Conditions (BCs). Such an approach can be accurate, in the presence of a correct mathematical formulation, only when the wave-length is larger than the size $\varepsilon$ of the unit cell. However it offers the advantage that the associated computational effort is either negligible in presence of analytical solutions to the equations of motion or dependent on the adopted mesh size in the case of finite element approximations.

This work develops a continuum model of the tensegrity mass–spring system recently proposed in ^[6] (see Section 2 for a review of this model). The approach that is followed moves on from the introduction of kinematic descriptors at the continuum level for the chain domain (i.e., 'hosting' and 'resonant' displacement fields) and next introduces appropriate Piola's ansatzes to link such descriptors to those corresponding to the discrete model ^[7,8,9,10]. The adopted ansatz guarantees that the continuum descriptors correspond to the limits for $\varepsilon \rightarrow 0$ of the discrete counterparts. Using the continuum displacement fields and a variational principle for the action functional system, the work derives the PDEs of the Maxwell chain as well as the dispersion relation at the continuum level (Section 3). Analytic solutions of the wave dynamics of the continuum model are presented in Section 4, which illustrate its main features and permit the validation of the dispersion relation presented in the previous section. The paper ends with concluding remarks and directions for future work in Section 5.

2.   A tensegrity Maxwell chain model

We hereafter recall the main features of the tensegrity Maxwell chain model diffusely presented in ^[6], which is graphically illustrated in Figure 1. It is composed of a number of $N$ tensegrity $\theta = 1$ prisms ('$P_H$ prisms') ^[14] arranged in parallel with $2 N$ minimal regular tensegrity prisms ('$P_R$ prisms' or T3 prisms) ^[16]. The $P_H$ prisms are interposed between $N+1$ hosting masses $M_H$, while two $P_R$ prisms are interposed between two consecutive $M_H$ masses, being separated by a resonant $M_R$ mass. The $M_H$ and $M_R$ masses exhibit characteristic sizes ($r_H$ and $r_R$) and rotational moments of inertia $J_H = n_H \, M_H \, r_H^2$ and $J_R = n_R \, M_R \, r_R^2$ with respect to the axis of the chain, respectively, where $n_H$ and $n_R$ are two scalar factors depending on the shape of such elements ($n_H = n_R = 1/2$ in the case of masses consisting of circular discs with radii $r_H$ and $r_R$).

The model analyzed hereafter describes the response of the chain in the small displacement regime from the reference configuration, by reducing the prism units to linear springs endowed with stiffness coefficients $K_H$ ($P_H$ prisms) and $K_R$ ($P_H$ prisms). The latter coincides with the tangent values of the axial stiffness coefficients of the $P_H$ and $P_R$ prisms in correspondence to the reference configuration and are assumed to be positive. Let $u^H_{i}$ denote the axial displacement of the $i-$th hosting mass and let $\vartheta^H_{i}$ denote the twisting rotation that accompanies the axial deformation of the $P_H$ prism interposed between the hosting masses $i$ and $i+1$. Similarly, let us use the symbols $u^R_{i}$ and $\vartheta^R_{i}$ to indicate the axial displacement and the twisting rotation of the $i-$th resonant mass ^[11], respectively (Figure 1). We assume the following compatibility equations between the axial and twisting motions of the hosting and resonant prisms ^[6]

Here, $p$ is a coupling parameter to be determined by studying the kinematics of the $P_H$ and $P_R$ prisms in proximity to the reference configuration. Equation (2.1) permits to reduce the kinematics of the Maxwell chain to the axial displacements $u^H_{i}$ and $u^R_{i}$ of the hosting and resonant units (see ^[6] for more details).

It is an easy task to obtain the expressions of the kinetic energy $\mathcal{K}$ and the internal elastic energy $\mathcal{U}_{int}$ of the discrete Maxwell chain model as follows:

where we have employed the dot notation for time derivatives, and we have introduced the equivalent masses

The external energy of the Maxwell chain is given by

where the following equivalent axial forces have been introduced:

${F}^{ext}_{H, i}$ and ${ F}^{ext}_{R, i}$ denote the external forces and the external twisting moments ${\cal M}^{ext}_{H, i}$ and ${\cal M}^{ext}_{R, i}$ applied to the chain.

It is worth to note that we decompose the external energy (2.5) into three addends. The first and third addends are the contributions from the external world, respectively, to the extreme left and to the extreme right-hand sides of the chain. The second addend contains the rest of the external energies. As a matter of fact, this last second addend is decomposed into two parts. The first is the contribution of the external energies due to the external forces applied to the hosting masses with $i = 1, ..., N-1$. The second is the contribution of the external energies due to the external forces applied to the resonance masses with $i = 0, ..., N-1$.

3.   A continuum model for the Maxwell chain

We now introduce a one-dimensional continuum domain of length $L$, characterized by two displacement fields $u_H(X, t)$ and $u_R(X, t)$ (see Figure 1) as kinematic descriptors, where $X \in \left[0, L\right]$ is a spatial variable (cf. Figure 1) and $t \in \left[t_1, t_2 \right]$ is the time variable. This domain represents the support of a continuum model for the Maxwell chain, which consists of a mixture of two solids: the first solid describes the hosting masses $M_H$ of the discrete model, while the second one describes the resonance masses $M_R$. Both solids are defined not only for every instant of time $t \in \left[t_1, t_2 \right]$ but also for every point $X \in \left[0, L\right]$ of the reference configuration. In addition, each of the two solids has its own dynamics, one being defined by the displacement field $u_H(X, t)$ and the other by the displacement field $u_R(X, t)$.

The following Piola's ansatz

are assumed to link the continuum kinematic descriptors to the discrete ones introduced in the previous section, where $f$ is either $u_H$ or $u_R$, and $f_i$ is either $u_i^H$ or $u_i^R$. Such ansatzes let us identify the realizations of the continuum descriptors $u_H(X, t)$ and $u_R(X, t)$ at the points of application of the $i$th discrete masses $M_H$ and $M_R$ with the discrete descriptors $u_i^H(t)$ and $u_i^R(t)$. Making use of Eq (3.1), we are led to rewrite the kinetic energy (2.2) of the discrete Maxwell chain as follows:

where the linear mass densities $\varrho_H$ and $\varrho_R$ are given by

Let us now look at the continuum model as the limit for $\varepsilon \rightarrow 0$ of a sequence of discrete chain models. In correspondence to such a limit, one easily realizes that Eq (3.3) do not degenerate only if the values $M_H$ and $M_R$ (see also Eq (2.4)) rescale with the size $\varepsilon$ of the unit cell of the discrete chain. We therefore assume that the values of the masses $M_H$ and $M_R$ tend to zero for $\varepsilon \rightarrow 0$, in such a way that the linear mass densities $\varrho_H$ and $\varrho_R$ remain finite.

Before we move on to identify the internal energy from Eq (2.3) of the continuum model, we need to write the displacement $u^H_{i+1}$ of the cell $i+1$ in terms of that $u^H_i$ of the cell $i$, with the use of both Piola's ansatz (3.1), adapted for the field $u^H(X, t)$, and the following Taylor's series expansion

Here, we have employed the prime notation for the derivatives with respect to the spatial variable $X$. Eq (3.4) allows us to write from Eq (2.3) the internal energy of the continuum model into the following form:

which, making use of Eq (3.1), leads us finally to

The continuum stiffness constants $\kappa_H$ and $\kappa_R$, which appear in Eq (3.8), are related to the discrete counterparts $K_H$ and $K_R$ through

We now go back to the identification of the continuum model as the limit for $\varepsilon \rightarrow 0$ of a sequence of discrete models. In such a limit, Eq (3.9) do not degenerate only if the values $\kappa_H$ and $\kappa_R$ rescale with $\varepsilon$ according to (3.9). This implies that it must result in $K_H \rightarrow \infty$ for $\varepsilon \rightarrow 0$, in such a way that $\kappa_H$ remains finite. Similarly, it must result in $K_R \rightarrow 0$ for $\varepsilon \rightarrow 0$, so that $\kappa_R$ remains finite. It is worth noting that inclusion of the last term in (3.5) would include higher-order gradient terms in the final PDEs that are here neglected. Besides, quadratic terms in (2.3) have implied linear PDEs. However, for large deformation, nonlinear PDEs would be obtained by considering exponents in (2.3) higher than $2$. In this general case, the difference between the formulations in the reference or in the present configurations should be considered via the so-called Piola's transformation ^[15]

We finally identify the external energy (2.5), making use of Piola's ansatz (3.1), for both the fields $u^H(X, t)$ and $u^R(X, t)$. We obtain the following expression of $\mathcal{U}_{ext}$ for the continuum model:

from which, making use of Eq (3.1), we get to

The distributed external forces $b^H(X, t)$ and $b^R(X, t)$, which appear in Eq (3.11), are related to the discrete counterparts ${\cal F}^{ext, eq}_{H, i}(t)$ and ${\cal F}^{ext, eq}_{R, i}(t)$ through

Considering again the limit for $\varepsilon \rightarrow 0$ of a sequence of discrete models, we realize that Eq (3.12) do not degenerate in correspondence to such a limit case only if the values $b^H$ and $b^R$ rescale with $\varepsilon$ according to (3.12). This implies that it must result in both ${\cal F}^{ext, eq}_{H, i} \rightarrow 0$ and ${\cal F}^{ext, eq}_{R, i} \rightarrow 0$ for $\varepsilon \rightarrow 0$, in such a way that $b^H$ and $b^R$ remain finite for $\varepsilon \rightarrow 0$.

The action functional of the continuum model can be formulated into the following standard form (see, e.g., ^[17])

where $t_1$ and $t_2$ are two times at which we prescribe that the displacement fields $u^H$ and $u^R$ are equal to known values $u_H^\alpha$ and $u_R^\alpha$ ($\alpha = 1, 2$), through

We now introduce the following action principle ^[17],

for any admissible variations $\mathcal{C}_H$ of $u_H$ and $\mathcal{C}_R$ of $u_R$. It is an easy task to obtain from Eqs (3.2), (3.8), (3.11), and (3.13) the following result:

Invoking the arbitrariness of the variations in (3.15), we can finally derive the following PDEs ruling the dynamic problem of the continuum model:

which are accompanied by the BCs

Let us now introduce plane wave solutions for the displacement fields $u_H(X, t)$ and $u_R(X, t)$, through

where $\mathrm{Re}$ indicates the real part operator; $I$ is the imaginary unit; $u_0^H$, $u_0^R$, $\omega$, and $k_w$, respectively, denote the complex amplitudes, the angular frequency, and the wave number of the traveling waves. Making use of Eq (3.21), we easily obtain (neglecting all the external forces and the external energy terms) the following dispersion relation of the continuum model

where the function $f_M(k_w)$ and the discriminant $\Delta_M\left(k_w\right)$ have the expressions given below

In the limit $\varepsilon \rightarrow 0$, Eq (3.22) reduce to the form

Making use of the positions (3.3) and (3.9), it is easily verified that both the Eqs (3.25) and (3.26) reduce to Eq (25), in the reference ^[6], of the discrete model, which has been presented to provide the value of the circular frequency of the optic branch for $k_w = 0$.

It is useful to apply the above formulas to the micro-scale physical models of the Maxwell chain studied in ^[6]. Such models employ $P_R$ prisms with a 5.5 mm height in the reference configuration, equilateral triangular bases with an 8.7 mm edge, 0.28 mm Spectra strings (5.48 GPa Young modulus), and 0.8 mm Ti6Al4V bars (120 GPa Young modulus). In addition, they use $P_H$ prisms showing 11 mm height in the reference configuration, equilateral triangular bases with a 6.11 mm edge, Spectra strings and Ti6Al4V bars identical to those of the $P_R$ prisms (we refer the reader to ^[6] for more detailed information about such units). The masses forming the above models are circular discs, such that $M_H = 16.03$ grams and $M_R = 4$ grams ($n_H = n_R = 1/2$). The study presented in ^[6] estimates $K_H = 61.04$ N/mm, $K_R = 13.70$ N/mm, and $p \approx 2.33$ mm/rad.

Graphical illustrations of the dispersion relations exhibited by discrete and continuum models of the chain under examination are provided in Figure 2, by distinguishing two different cases: discrete and continuum systems in which the longitudinal and twisting motions of the $P_H$ and $P_R$ prisms are coupled, according to Eq (2.1), as well as discrete and continuum systems in which the same prisms are free to tangentially slide against the lumped masses, so as that the twisting rotations of the tensegrity units are not transferred to the masses (uncoupled system). In the first case, we have $M^{eq}_H = 349.51$ grams and $M^{eq}_R = 24.85$ grams, while in the second case, it results in $M_H^{eq} = M_H$ and $M_R^{eq} = M_R$. Figure 2 shows that the acoustic and optic branches forming the dispersion relations of the continuum systems correctly reduce to the analogous branches of the discrete systems in the long wavelength regime ($k_w \rightarrow 0$). In the coupled system we predict $\omega_{Moc}(0) = 173$ Hz, while in the uncoupled system, we predict $\omega_{Moc}(0) = 465$ Hz.

4.   Analytic solutions of the wave equation

Let us derive the wave equation of the continuum model of the Maxwell chain, by making use of the results derived in the previous section. We start by rewriting Eq (3.21) into the following form:

where now $\bar{u}^H_0(X)$ and $\bar{u}^R_0(X)$ are unknown functions of the spatial coordinate $X$, and $u_S^H(X)$ and $u_S^R(X)$ are those displacement functions produced by a static load at the boundary and assuming zero external distributed forces ($b_H = b_R = 0$). Thus, Eqs (3.17) and (3.18) yield,

and $u_S^H(X)$ and $u_S^R(X)$ are solved once a particular set of boundary conditions are considered. Keeping this in mind, one can write Eqs (3.17) and (3.18) making use of Eq (4.1), to obtain

The insertion of Eq (4.5) into Eq (4.4) leads us to the second-order linear homogeneous differential equation with respect to the spatial coordinate $X$

where we have introduced the equivalent mass density

Eq (4.6) admits the following general solution:

where we introduced the integration constants $c_1$ and $c_2$, to be determined through the competent boundary conditions, and we set

Moving on to determine analytic solutions for $\bar{u}_0^R(X)$, we now insert Eq (4.8) into Eq (4.5), obtaining

Let us now determine the constants $c_1$ and $c_2$ in the particular case of a chain lying vertically in space, which is subject to a mass $M$ applied to the top end ($X = L$) and a sinusoidal displacement excitation at the bottom end ($X = 0$). The base excitation has an amplitude of $u_0$ and can be written as

It is useful to separate the dynamic and static parts of Eq (4.1). For what concerns the dynamic part, we analyze the boundary conditions

which, once inserted into Eq (4.8), leads us to compute the integration constants $c_1$ and $c_2$, and to obtain the dynamic amplitude of the hosting displacement field

and from Eq (4.10) that of the resonant one,

For what concerns the static terms $u_S^H$ and $u_S^R$, we observe that Eq (4.3) implies the following result in static equilibrium conditions:

which is accompanied by the static boundary condition

On the other hand, in the example under consideration, it is a simple task to obtain the following differential equation for $u_S^H$ ^[20]

where $g$ is the gravity acceleration. Making use of Eqs (4.16) and (4.17) in Eq (4.15), we finally obtain

We are now able to cast Eq (4.1) into the following form:

Since we are assuming $\kappa_H > 0$ (cf. Section 2), Eqs (4.19) and (4.20), in association with the positions (4.9), let us conclude that one obtains periodic harmonic solutions for $u_H$ and $u_R$ only when $\alpha_{1}$ and $\alpha_{2}$ are imaginary ($\bar{ \varrho}_H > 0$). Oppositely, we obtain decaying solutions when $\alpha_{1}$ and $\alpha_{2}$ are real ($\bar{ \varrho}_H < 0$). The frequency bandgap region obviously corresponds to $\bar{ \varrho}_H < 0$, since in this case the system supports wave solutions that exponentially decay when $X$ approaches $L$ ^[18,19]. It is easily verified that the values of $\omega > 0$ such that $\bar{ \varrho}_H < 0$ are contained in the interval comprised between zero and the frequency $\omega_{Moc}$ provided by Eq (3.26). Such a result provides validation of the dispersion relation derived in the previous section, and highlights that the analyzed chain can be effectively employed as an isolation device for vertical vibrations of the top mass $M$. It is not difficult to generalize the results obtained in the present section to the case in which the examined system is under the action of a not-negligible self-weight. For such an example, indeed, the dynamic parts of $u_H$ and $u_R$ remains unchanged as compared to the present case, while the static parts of these functions will assume a quadratic expression with respect to $X$ ^[20].

5.   Conclusions

We have derived a continuum model for a Maxwell-type mass–spring chain with tensegrity architecture ^[12,13], which exhibits $\theta = 1$ tensegrity prisms ^[14] arranged in parallel with T3 prisms ^[16] and lumped masses. When compared to the discrete model presented in ^[19], the theory formulated in the present work is useful to obtain analytic solutions for the wave propagation problem of the analyzed system in the high wave-length regime, where structural applications are relevant ^[21,22]. Such a feature has been demonstrated through the analysis of a physical example, which refers to the vibration isolation problem of a mass, e.g., a device to be protected against mechanical vibrations in a hospital or another essential building. It is worth noting that the mass–spring chain studied in this work exhibits internal resonance properties that can be employed to widen the frequency bandgap width of the system, as compared to the case of a diatomic tensegrity mass–spring chain, where such a width can only be tuned by playing with mass and stiffness contrasts between the elements of the unit cell ^[23].

The research presented in this work paves the way for the design of novel tensegrity metamaterials serving as next-generation vibration isolation devices. Such systems will exhibit properties that can be easily adjusted to the structure to be protected due to the tunability of the frequency bandgap region ^[19]. The presence of geometrical and mechanical nonlinearities in such devices will be studied through future research, accounting, e.g., for the presence of superelastic cables made of shape memory alloys ^[24] and higher order terms in the internal an kinetic energies ^[25,26]. Their use in forming anisotropic and/or porous metamaterials ^[27,28], energy harvesters ^[29], and sensors and actuators for structural health monitoring ^[30] will also be investigated in future studies.

Author contributions

L.P. formulated the initial idea for the continuum model, while F.F., J.d.C.M., and R.Z.C. conceived the tensegrity modeling. L.P. and F.F. led the conceptualization, idealization, and supervision phases of the project. F.F. provided funding. J.d.C.M. led the development of the analytical results, while R.Z.C. led the numerical analysis of the physical models. All authors contributed to the first draft and the revised version of the manuscript.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Acknowledgements

This research has been funded under the National Recovery and Resilience Plan (NRRP), by the European Union–NextGenerationEU, within the project with grant number P2022CR8AJ (FF PI). This research has also been funded by the NextGenerationEU PRIN2022 research projects with grant numbers 20224LBXMZ (FF PI) and 2022P5R22A ("Ricerca finanziata dall'Unione Europea–Next Generation EU"). FF and JdCM also acknowledge the support by the Italian Ministry of Foreign Affairs and International Cooperation within the Italy-USA Science and Technology Cooperation Program 2023–2025, Project "Next-generation green structures for natural disaster-proof buildings" (grant No. US23GR15).

Conflict of interest

The authors declare that there is no conflict of interest.