Processing math: 60%
Research article Special Issues

Free vibration of summation resonance of suspended-cable-stayed beam

  • Received: 29 April 2019 Accepted: 01 July 2019 Published: 08 August 2019
  • Free vibration of summation resonance of suspended-cable-stayed beam is investigated in the article. A 3-DOF model of the coupled structure is built, with the main cable and sling (vertical cable) considered to be geometrically nonlinear, and the beam is taken as linear Euler beam. Hamilton's principle is used to derive the dynamic equilibrium equations of the coupled structure. Then, the dynamic equilibrium equations are solved by means of multiple scales method, the second order approximation solutions of single-modal motion of the coupled structure are obtained. Numerical examples are presented to discuss time history of free vibration of the summation resonance, with and without damping. Additionally, fourth-order Runge-Kutta method is directly used for the dynamic equilibrium equations to complement and verify the analytical solutions. The results show that the coupled structure performs strongly nonlinear and coupled characteristics, which is useful for engineering design.

    Citation: Chunguang Dong, Zhuojie Zhang, Xiaoxia Zhen, Mu Chen. Free vibration of summation resonance of suspended-cable-stayed beam[J]. Mathematical Biosciences and Engineering, 2019, 16(6): 7230-7249. doi: 10.3934/mbe.2019363

    Related Papers:

    [1] Abigail Gbemisola Adeyonu, Olubunmi Lawrence Balogun, Babatunde Oluseyi Ajiboye, Isaac Busayo Oluwatayo, Abiodun Olanrewaju Otunaiya . Sweet potato production efficiency in Nigeria: Application of data envelopment analysis. AIMS Agriculture and Food, 2019, 4(3): 672-684. doi: 10.3934/agrfood.2019.3.672
    [2] Saeed Hosseini Nejad, Reza Moghaddasi, Amir Mohammadi Nejad . On the role of credit in agricultural growth: An Iranian panel data analysis. AIMS Agriculture and Food, 2018, 3(1): 1-11. doi: 10.3934/agrfood.2018.1.1
    [3] Lindie V. Maltitz, Yonas T. Bahta . Empowerment of smallholder female livestock farmers and its potential impacts to their resilience to agricultural drought. AIMS Agriculture and Food, 2021, 6(2): 603-630. doi: 10.3934/agrfood.2021036
    [4] Ogunniyi Adebayo, Omonona Bolarin, Abioye Oyewale, Olagunju Kehinde . Impact of irrigation technology use on crop yield, crop income and household food security in Nigeria: A treatment effect approach. AIMS Agriculture and Food, 2018, 3(2): 154-171. doi: 10.3934/agrfood.2018.2.154
    [5] Raulston Derrick Gillette, Norio Sakai, Godfrid Erasme Ibikoule . Role and impact of contract farming under various pricing standards: A case of Guyana's rice sector. AIMS Agriculture and Food, 2024, 9(1): 336-355. doi: 10.3934/agrfood.2024020
    [6] Hien Thi Vu, Ke-Chung Peng, Rebecca H. Chung . Efficiency measurement of edible canna production in Vietnam. AIMS Agriculture and Food, 2020, 5(3): 466-479. doi: 10.3934/agrfood.2020.3.466
    [7] Abbas Ali Chandio, Yuansheng Jiang . Factors influencing the adoption of improved wheat varieties by rural households in Sindh, Pakistan. AIMS Agriculture and Food, 2018, 3(3): 216-228. doi: 10.3934/agrfood.2018.3.216
    [8] Marcello Mastrorilli, Raffaella Zucaro . Towards sustainable use of water in rainfed and irrigated cropping systems: review of some technical and policy issues. AIMS Agriculture and Food, 2016, 1(3): 294-314. doi: 10.3934/agrfood.2016.3.294
    [9] Radhwane Derraz, Farrah Melissa Muharam, Noraini Ahmad Jaafar . Uncertainty sources affecting operational efficiency of ML algorithms in UAV-based precision agriculture: A 2013–2020 systematic review. AIMS Agriculture and Food, 2023, 8(2): 687-719. doi: 10.3934/agrfood.2023038
    [10] Boris Boincean, Amir Kassam, Gottlieb Basch, Don Reicosky, Emilio Gonzalez, Tony Reynolds, Marina Ilusca, Marin Cebotari, Grigore Rusnac, Vadim Cuzeac, Lidia Bulat, Dorian Pasat, Stanislav Stadnic, Sergiu Gavrilas, Ion Boaghii . Towards Conservation Agriculture systems in Moldova. AIMS Agriculture and Food, 2016, 1(4): 369-386. doi: 10.3934/agrfood.2016.4.369
  • Free vibration of summation resonance of suspended-cable-stayed beam is investigated in the article. A 3-DOF model of the coupled structure is built, with the main cable and sling (vertical cable) considered to be geometrically nonlinear, and the beam is taken as linear Euler beam. Hamilton's principle is used to derive the dynamic equilibrium equations of the coupled structure. Then, the dynamic equilibrium equations are solved by means of multiple scales method, the second order approximation solutions of single-modal motion of the coupled structure are obtained. Numerical examples are presented to discuss time history of free vibration of the summation resonance, with and without damping. Additionally, fourth-order Runge-Kutta method is directly used for the dynamic equilibrium equations to complement and verify the analytical solutions. The results show that the coupled structure performs strongly nonlinear and coupled characteristics, which is useful for engineering design.


    As history indicated, infectious diseases have been becoming a main enemy affecting human's health and economic development. Mathematical modelling is an useful tool to investigate the mechanisms of transmission of diseases and make optimal control measures [1]. There are many ways to suppress the disease transmission, such as media propagation, vaccination, quarantine and so on [2]. As we know, vaccination is one of effective methods to control and prevent disease prevalence. Indeed, vaccination has succeeded in slowing down transmission of diseases such as tuberculosis, hepatitis, and some children diseases [3]. However, it has been reported that vaccination immunity waning has caused some diseases reemergence such as measles, rubella and pertussis. There is no doubt that vaccine waning has great effects on understanding the evolution of diseases. Based on the epidemic compartment knowledge in Kermack and McKendrick [4,5,6], an SIVS epidemic model can be written as follows [7]:

    dS(t)dt=bβSI(μ+p)S+ϵV,dI(t)dt=βSI(μ+γ)I,dR(t)dt=γI(t)μR,dV(t)dt=pS(μ+ϵ)V(t), (1)

    where the total population is splited into four classes (Susceptible, Infected, Recovered and Vaccinated). b is the birth rate, μ is the natural death rate, β is the transmission rate, γ is the recovery rate, p is the vaccinated rate, ϵ is the vaccine waning rate. In [7], Li and Yang investigated two vaccine strategies consisting of continuous and impulsive styles, and they obtained the global stability of equilibria by constructing Lyapunov functionals. Zaleta and Hernández proposed an SIVS model with a standard incidence rate and a disease-induced death rate [8]. They showed that their model exhibits a backward bifurcation. Based on model (1), many researchers have evolved many different structures and successfully captured the key characters of diseases transmission and evaluated the risk of their prevalence (see, for examples, [9,10,11]).

    We note that all models mentioned above are based on the homogeneous assumptions for host population. However, host heterogeneity plays an important role in exploring their dynamics. Many diseases such as tuberculosis, hepatitis C, HIV/AIDS and so on infect their host for a long time and sometimes for the duration of lifespan. During the long infectious period, the variability of infectivity with age-since-infection has been studied most extensively in HIV infection [12]. The immunity waning process of pertussis satisfies Gamma distributions [13], which can be expressed by the vaccinated age. In this case, the densities of the infected and vaccinated in time t and age a are denoted by i(t,a), and v(t,a), respectively. The parameters β,γ, and ϵ in system (1) are associated with age a. Based on system (1), the model can be described as follows:

    dS(t)dt=bS(t)0β(a)i(t,a)da(μ+p)S(t)+0ϵ(a)v(t,a)da,i(t,a)t+i(t,a)a=(μ+γ(a))i(t,a),v(t,a)t+v(t,a)a=(μ+ϵ(a))v(t,a),i(t,0)=S(t)0β(a)i(t,a)da,v(t,0)=pS(t),dR(t)dt=0γ(a)i(t,a)daμR(t),S(0)=S0R+,i(0,a)=i0(a)L1+(R+),v(0,a)=v0(a)L1+(R+), (2)

    where R+=(0,+), and L1+(R+) denotes the space of all the integral functions in L1 and maintaining positivity after integral. Obviously, system (2) is a hybrid system combining an ordinary differential equation and two partial differential equations. Global dynamics of such a system has been becoming a challenging issue due to lack of well posed mathematical techniques.

    On the other hand, some diseases such as mumps, measles, gonorrhea, HIV/AIDS etc exhibit heterogeneity in host populations. Groups can be geographical such as counties, cities, communities, or epidemiological as different infectivity and multi-stain agents. Many authors introduced an irreducible matrix and summing elements of the matrix together as a kernel function to describe the inter-group and intra-group infections. In [14], Lajmanovich and Yorke proposed the earliest multi-group model for gonorrhea spread in a community and investigated the global stability. They assumed that the total size of population doesn't change and maintain a constant. Under this assumption their model can be simplified to just consider the infected classes. Since then, many multi-group epidemic models have been studied (see, for example [15,19,20]). In [21], Guo et al proposed a multi-group SIR epidemic model and used the graph-theoretic approach to investigate the global stability of endemic equilibrium. This method combining with nonnegative irreducible matrix is an effective tool in solving the global behavior of endemic equilibrium. In [20], Kuniya considered a multi-group SVIR model to explore the global behavior of equilibria by constructing Lyapunov functional and using a developed graph-theoretic method. There are few literatures incorporating age structure into multi-group epidemic models [15].

    Motivated by the discussions above, we separate the total population into n groups and four compartments: susceptible, infected, and vaccinated, recovered, denoted by Sk(t),ik(t,a) and vk(t,a) and Rk(t), respectively. ik(t,a) denotes the infected individuals at time t and infection age a in group k. vk(t,a) denotes the vaccinated individuals at time t and vaccinated age a in group k. Rk(t) represents the recovered individuals at time t. Susceptible individuals in group k can be infected by infected individuals in group j at rate βkj(a). Hence, we denote the incidence rate in group k in the form of

    λk(t)=nj=10βkj(a)ij(t,a)da.

    Susceptible individuals in group k can be vaccinated at rate pk and become vaccinated individuals with immunity. The vaccinated individuals in group k lose its immunity at rate ϵk(a) and become susceptible individuals. A multi-group SIRVS epidemic model is formulated by the following differential equations:

    dSk(t)dt=bkSk(t)λk(t)(μk+pk)Sk(t)+0ϵk(a)vk(t,a)da,ik(t,a)t+ik(t,a)a=(μk+γk(a))ik(t,a),vk(t,a)t+vk(t,a)a=(μk+ϵk(a))vk(t,a),ik(t,0)=Sk(t)λk(t),vk(t,0)=pkSk(t),dRk(t)dt=0γk(a)ik(t,a)daμkRk(t), (3)

    where bk denotes the birth rate, γk is the recovery rate with respect to infection age a, μk is the natural death rate in group k. In order to satisfy the biological meaning, all the parameters are assumed to be nonnegative and bk>0 and μk>0. Note that the total population Nk(t)=Sk(t)+0ik(t,a)da+Rk(t)+0vk(t,a)da satisfies the following equation:

    dNk(t)dt=bkμkNk(t), (4)

    which yields limtNk(t)=bkμk. Without loss of generality, we assume that the total population is Nk=bkμk. Since the first four equations in (3) do not contain the variable Rk, we can consider the following closed subsystem:

    dSk(t)dt=bkSk(t)λk(t)(μk+pk)Sk(t)+0ϵk(a)vk(t,a)da,ik(t,a)t+ik(t,a)a=(μk+γk(a))ik(t,a),vk(t,a)t+vk(t,a)a=(μk+ϵk(a))vk(t,a),ik(t,0)=Sk(t)λk(t),vk(t,0)=pkSk(t),Sk(0)=Sk0R+,ik(0,a)=ik0(a)L1+(R+),vk(0,a)=vk0(a)L1+(R+). (5)

    Once behaviors of Sk(t), ik(t,a) and vk(t,a) are known, those of Rk(t) can be derived from the fourth equation in (3). For convenience, we make the following assumption:

    Assumption 1.1. For system (5), we assume

    (ⅰ) For each j,k{1,2,,n}, βjk(a),γk(a),ϵk(a)L+(0,). That is, there exist positive constants β+jk and ϵ+k such that

    esssupa[0,+)βjk(a)=β+jk,esssupa[0,)γk(a)=γ+k,esssupa[0,)ϵk(a)=ϵ+k.

    (ⅱ) For each j,k{1,2,,n}, βjk(a) satisfies the following property:

    limh00|βjk(a+h)βjk(a)|da=0.

    (ⅲ) For each j,k{1,2,,n}, there exists an ε0>0 such that for almost all a[0,+),βjk(a)ε0.

    (ⅳ) For each j,k{1,2,,n}, the matrix (βjk)n×n is irreducible.

    The assumptions above on the parameters in system (5) are naturally satisfied for some real diseases. Evidence exists that the infectivity βjk(a) has been addressed by Gamma and Log-normal distributions for smallpox [16,17], Weibull distribution for Ebola [18]. It is easy to find that all of these probability distribution functions have peak values and they are continuous. Hence, (ⅰ) and (ⅱ) in Assumption 1.1 readily hold. As for (ⅲ) in Assumption 1.1, we can modify it in a more generalized form:

    (ⅲ)' There exists a positive constant}aβ such that βjk(a)(k,jN) is positive in a neighbourhood of aβ.

    Obviously, the distribution functions mentioned above have this property. If we assume that every group keeps up close exchanges in mutual contact, then the generated graph is strongly connected. (ⅳ) in Assumption 1.1 is satisfied automatically.

    In order to investigate the dynamic behavior of system (5), define the functional spaces

    X=(R×L1(R+,R))n,Y=X×X,Z=Rn×X×X,

    and

    X0=({0}×L1(R+,R))n,Y0=X0×X0,Z0=Rn×X0×X0

    with the norm

    ϕRn=nj=1|ϕj|,ϕ=(ϕ1,ϕ2,,ϕn)TRn,
    ϕX0=nj=10|ϕj(a)|da,ϕ=(ϕ1,ϕ2,,ϕn)TL1(R+,Rn),
    ψY0=ψ1X0+ψ2X0,ψZ=ψ1Rn+ψ2X0+ψ3X0,

    where ψi=(ψi1,ψi2,,ψin)TRn or L1(R+,Rn)(i=1,2,3). In addition, denote X+,Y+ and Z+ as the positive cones of X,Y and Z. We define

    X0+=X0X+,Y0+=Y0Y+,Z0+=Z0Z+.

    Under Assumption 1.1 we see that the set

    Ω={(S(t),0,i(t,),0,v(t,))Z0+|Sk(t)+0ik(t,a)da+0vk(t,a)dabkμk}

    is invariant, where S=(S1,S2,,Sn)T, i=(i1,i2,,in)T and v=(v1,v2,,vn)T. In the following, we just assume that all the initial values are taken from Ω.

    Next, we will show that system (5) has a globally classic solution in Ω. Let b,p,μ,γ(a), and ϵ(a) be diagonal matrixes given by

    b=diag(b1,b2,,bn),p=diag(p1,p2,,pn),μ=diag(μ1,μ2,,μn),γ(a)=diag(γ1(a),γ2(a),,γn(a)),ϵ(a)=diag(ϵ1(a),ϵ2(a),,ϵn(a)). (6)

    Then let us define a linear operator Ai:D(Ai)XX as

    Ai(0ψ)=(0Aiψ)=(0ddaψ(μ+γ(a))ψ),

    and

    Av(0ψ)=(0Avψ)=(0ddaψ(μ+ϵ(a))ψ),

    where D(Aj)={(0,ψ)X0|ψ is obsolutely continuous,ψL1(R+,Rn),ψ(0)=0}, j=i,v. If λρ(Ai) (ρ(Ai) denotes the resolvent set of Ai), for any initial value (θi,ϕi)TX, we have

    ψi(a)=ea0[μi+γi(s)+λ]dsθi+a0ϕi(s)eas[μi+γi(ξ)+λ]dξds. (7)

    Similarly, for any λρ(Av) and any initial value (θv,ϕv)TX, we have

    ψv(a)=ea0[μi+ϵi(s)+λ]dsθv+a0ϕi(s)eas[μi+ϵi(ξ)+λi]dξds. (8)

    Furthermore, define two nonlinear operators as

    Fj(0ψ)=(Bj(ψ)0),j=i,v,

    where Bi(ϕ)=Snj=10βj(a)ϕij(a)da and Bv(ϕ)=pϕS. We further define another nonlinear operator

    FS(ψ)=bpψSBi(ψi)+0ϵ(a)ψv(a)da.

    If we set u=(S,(0i),(0v))TZ, A=diag(μ,Ai,Av), and F=(Fs,Fi,Fv), (5) can be rewritten as the following abstract Cauchy problem

    du(t)dt=Au(t)+F(u(t)),u(0)=u0Ω. (9)

    Proposition 2.1. There exists a uniquely determined semiflow {U(t)}t0 on Z0+ such that, for each u=(S(t),(0i(t,)),(0v(t,)))Z0+, there exists a unique continuous map UC(R+,Z0+) which is an integral solution of the Cauchy problem (9), that is, for all t0,

    U(t)u=u0+At0U(s)uds+t0F(U(s)u)ds. (10)

    Proof. By the definition of A and equations (7) and (8), A is dense in part of X0 (see Page 1117, , [22]), and the resolvent operator is bounded. From Proposition 3.2 in, [23], we need only to verify that F satisfies Lipschitz condition. Denote

    xk=(ψSk,(0ψik),(0ψvk))TZ,ˉxk=(ˉψSk,(0ˉψik),(0ˉψvk))TZ,

    and

    LSk=max{pk+pk0ϵk(a)da+bkμknj=1β+kj,bkμknj=1β+kj}.

    By calculation, we have

    |FSk(xk)FSk(ˉxk)|=|pk(ψSkˉψSk)+ψSknj=10βkj(a)ψij(a)daˉψSknj=10βkj(a)ψij(a)da+0ϵk(a)dapk(ψSkˉψSk)|pk|ψSkˉψSk|+nj=10βkj(a)ψij(a)da|ψSkˉψSk|+ψSknj=10βkj(a)|ψij(a)ˉψij(a)|da+p0ϵk(a)da|ψSkˉψSk)|pk|ψSkˉψSk|+bkμknj=1β+kj|ψSkˉψSk|+bkμknj=1β+kjψijˉψijL1+pk0ϵk(a)da|ψSkˉψSk)|LSkxˉx. (11)

    Similarly, we have Flk(xk)Flk(ˉx)Llkxˉx,l=i,v, where Lik=bkμknj=1β+kj, and Lvk=pk. Therefore, F(x)F(ˉx)nk=1Fk(xk)F(ˉxk)Lxˉx, where L=max{LSk,Lik,Lvk}.

    In order to prove the positivity of the solution of (5), we first rewrite the operator Av and the nonlinear function Fvk as follows:

    Av(0ψ)=(0Avψ)=(0ddaψμψ)

    and

    Fv(0ψ)=(pS0ϵ(a)ψ(a)da0).

    Denote A=diag(μ,Ai,Av), and F=(FS,Fi,Fv). Then system (9) can be written as

    du(t)dt=Au(t)+F(u(t)),u(0)=u0. (12)

    It is obvious that eAtu0Ω+ if u0Ω+, where Ω is the positive cone of Ω. Then we need to show that

    t0eA(ts)Fk(Ω+)dsΩ+.

    For any ψlkR and ϕlkΩ+, it follows from the definitions of A and F that

    ψ1+ψ2+ψ3=t0eμ(ts)FS(ϕ(s))ds+t0eμ(ts)tsγ(a)daFi(ϕ(s))ds+t0eμ(ts)Fv(ϕ(s))dst0bkeμ(ts)ds=bkμk(1eμkt)bkμk. (13)

    On the other hand, choose a κkR+ and redefine the operator in (12) as follows: Aκk=diag((μk+κk),Aiκk,Avκk), where

    Aiκk(0ψ)=(0ddaψ(μk+γk(a)+κk)ψ),

    and

    Avκk(0ψ)=(0ddaψ(μk+κk)ψ).

    The nonlinear functions in (12) are defined as Flκk(ϕ)=Flk+κkϕk(l=S,i,v) for any ϕZ. Hence, we have

    t0e(μk+κk)(ts)FSκk(ϕk(s))ds=t0e(μk+κk)(ts)[bkpkϕSk(s)ϕSk(s)nj=10βkj(a)ϕik(t,a)da+t0ϵk(a)ϕvk(a)da+κkϕSk]dst0e(μk+κk)(ts)[bk+(κkpknk=1β+kj)ϕij)]ϕSk(s)ds, (14)
    t0e(μk+κk)(ts)es0γk(a)daFiκk(ϕk(s))ds=t0e(μk+κk)(ts)[ϕSk(s)nj=10βkj(a)ϕij(t,a)da+κkϕik]dst0e(μk+κk)(ts)[ϕSk(s)ε0+κk]ϕikds, (15)

    and

    t0e(μk+κk)(ts)Fvκk(x)ds=t0e(μk+κk)(ts)[pkϕSk(s)0ϵk(a)ϕvk(s,a)da+κϕvk]dst0e(μk+κk)(ts)[pkϕSk(s)+(κkϵ+k)ϕvk]ds. (16)

    If we choose κk>max{pk+nj=1β+kjϕij,ϵ+k}, it follows from (14)-(16) that

    t0eAκk(ts)Fκk(ϕ(s))ds>0,ifϕ>0.

    From what has been discussed above, we have the following result.

    Proposition 2.2. If Assumption 1.1 holds, system (5) has a unique positive solution in Ω.

    In this section, we show the computation process of the basic reproduction number R0, which is the average number of secondary cases produced by a classical infected individual during its infectious period in a fully susceptible population. Many literatures have given different methods to solve this problem. In this paper, it follows from the renewal process mentioned in Diekman et al in [24] that its value is determined by a next generation operator. Note that E0=(ˆS0,0,0,0,ˆv0(a)) is the disease-free equilibrium of system (5) where ˆS0=(S10,S02,,S0n) and ˆv0(a)=(v01(a),v02(a),,v0n(a)) with S0k=bkμk+pk(1K1k),v0k(a)=pkS0kπ1k(a), and π1k(a)=ea0(μk+ϵk(s))ds, K1k=0ϵk(a)π1k(a)da. Linearizing system (5) at the disease-free equilibrium E0, the subsystem (5) can be rewritten as

    i(t,a)t+i(t,a)a=(μ+γ(a))i(t,a),i(t,0)=S0Λ(i(t,)), (17)

    where

    Λ(i(t,))=diagjN(nk=10βjk(a)ik(t,a)da),

    and

    S0=diag(S01,S02,,S0n).

    Let B:D(B)XX be a linear operator defined by

    B(0ϕ(a))=(ϕ(0)ddaϕ(a){μ+γ(a)}ϕ(a)),D(B)={ϕX:ϕis absolutely continuous, ϕL1(R+)}, (18)

    where μ and γ are defined in (6). Based on the above definition, (7) can be rewritten as the following linear Cauchy problem in Z0:

    ddtI(t)=BI(t)+S0ΛI(t),I(0)=I0X, (19)

    Let u(t)=etB be the C0 semigroup generated by the generator B. By the variation of constants formula, we have

    I(t)=u(t)I0+t0u(ts)S0ΛI(s)ds.

    Mapping S0Λ on both sides of the above equation yields

    v(t)=f(t)+t0Ψ(s)v(ts)ds, (20)

    where v(t)=S0ΛI denotes the density of newly infectives in the linear invasion phase, f(t)=S0Λu(t)I0 and Ψ(s)=ΛS0u(s). Then the next generation operator is defined by

    K=0Ψ(s)ds=Λ(S0B)1,

    where we used the relation (zB)1=0ezsV(s)ds, for zρ(B)(ρ(B) denotes the resolvent set of B) and 0ρ(B). In fact, the next generation operator is calculated as follows:

    K=(K1,K2,,Kn)T, (21)

    where

    Kj=nk=1S0k0βjk(a)πk(a)da,j=1,2,,n. (22)

    Based on Diekmann et al [24], the basic reproduction number R0=r(K) is the spectral radius of the next generation operator K, where r(A) denotes the spectral radius of a bounded operator A. It follows from Inaba [25] that Malthusian parameter or asymptotic growth rate of infectives is positive if R0>1, otherwise it is negative.

    In this section, we focus on the existence of the endemic equilibrium E of system (5). It follows from Section 2 that equilibria of system (5) satisfy Au+F(u)=0. Actually, it satisfies the following equations

    0=bk(μk+pk)SkSkλk+0ϵk(a)vk(a)da,dik(a)da=(μk+γk(a))ik(a),ik(0)=Skλk,λk=nj=10βkj(a)ik(a)da,dvk(a)da=(μk+ϵk(a))vk(a),vk(0)=pkSk. (23)

    From the last two equations of (23), we have

    vk(a)=pkSkπ1k(a),π1k(a)=ea0[μk+ϵk(s)]ds. (24)

    Substituting (24) into the first equation of (23) yields

    Sk=bkμk+pk(1K1k)+λk. (25)

    It follows from the second and the third equations of (23) that

    ik(a)=ik(0)πk(a)=Skλkπk(a),πk(a)=ea0[μk+γk(s)]ds.

    Note that

    λk=nj=10βkj(a)ij(a)da=nj=1ij(0)Kkj=nj=1bjλjμj+pj(1K1j)+λjKkj, (26)

    where Kkj=0βkj(a)πj(a)da. Hence, we can define a nonlinear operator

    H(ϕ):=(H1(ϕ),H2(ϕ),,Hn(ϕ))TX,ϕR, (27)

    where

    Hk(ϕ)=nj=1bjϕjμj+pj(1K1j)+ϕjKkj.

    In fact, it follows from (26) that the endemic equilibrium of (23) is a positively nontrivial fixed point of the operator H. Note that the Fréchet derivative of H at ϕ=0 is given by

    Hk[0]=limh0Hk[h]Hk[0]h=limh0nj=1hbjKkjh(μj+pj(1K1j)+h)=limh0nj=1bjKkjμj+pj(1K1j)+h=Kk.

    Let

    H[0]=(ˆH1,ˆH2,,ˆHn),

    where

    ˆHk=nj=1bjμj+pj(1K1j)Kkj=nj=1S0jKkj=Kk,k=1,2,,n.

    Consequently, H[0] is equal to the next generation operator K defined by (21).

    Next, we show that R0 determines the existence of the positive fixed point of operator H. Under Assumption 1.1, the following lemma holds.

    Lemma 4.1. Let K be defined by (21). We have

    (a) K is compact.

    (b)K is nonsupporting.

    Proof. Assume that B0 is an arbitrary bounded subset of R. Then there exists a positive constant c0 such that ϕc0 for all ϕB0. Note that

    K(ϕ)=nj=1|Kj(ϕ)|nj=1nk=1S0k0βjk(a)|ϕk(a)|danj=1nk=1S0kβ+jk0|ϕk(a)|da=nj=1nk=1S0kβ+jkϕkL1. (28)

    (28) implies that the operator K is bounded. It follows from Fréchet-Kolmogorov Theorem [31] that the operator K is compact. It is obvious that K is nonsupporting under (ⅲ) or (ⅲ)' of Assumption 1.1.

    Lemma 4.1, together with the monotonicity of the operator H with respect to ϕ manifests that the following lemma holds.

    Lemma 4.2. Let H be defined (27). H is compact and H(R+) is bounded.

    Employing Theorem 4.11 in [26] (Krasnoselskii fixed theorem) and Krein-Rutman Theorem in [27], Lemma 4.2 implies that r(K) is the uqiue eigenvalue of the operator K associated with a positive eigenvector and there is no eigenvector of K associated with eigenvalue 1. As a consequence of Corollary 5.2 in [15], the following result is immediate.

    Proposition 4.1. If R0>1, then H has at least one nontrivial fixed point in Z0+.

    Corollary 4.1. If R0>1, then system (5) has at least one positive endemic equilibrium E=(S,0,i,0,v)Z0+, where S=(S1,S2,,Sn),i=(i1,i2,,in),v=(v1,v2,,vn).

    In order to determine the uniqueness of the solution of system (23), we define two operators by (ⅱ) and (ⅲ) of Assumption 1.1 as follows:

    H+k(ϕ)=nj=1bjβ+kjϕjμj+pj(1K1j)+ϕj0πj(a)da,kN,

    and for a small positive value δ

    Hk(ϕ)=nj=1bjε0ϕjμj+pj(1K1j)+ϕjaβ+δaβδπj(a)da,kN.

    For convenience, define

    H+(ϕ)=diag(H+1(ϕ),H+2(ϕ),,H+n(ϕ)),

    and

    H(ϕ)=diag(H1(ϕ),H2(ϕ),,Hn(ϕ)).

    It is obvious that 0H(ϕ)eH(ϕ)eH+(ϕ)e, where e=(1,1,,1)T.

    Theorem 4.3. If R0>1, then the operator H defined in (27) has at most one nontrivial fixed point in Z0+.

    Proof. By Proposition 4.1, we assume that there are two different nontrivial fixed points in R+, denoted by ϕ and ˆϕ. Borrowing the definitions of H and H+, we have

    ϕ=H(ϕ)H(ϕ)=H(ϕ)H+(ϕ)H+(ϕ)H(ϕ)H+(ϕ)H(ϕ)=H(ϕ)H+(ϕ)ϕ.

    Hence, there exists a positive constant q=sup{r0,ϕrˆϕ}>0. Suppose 0<q<1, for all ϕR+ and kN, then

    Hk(qϕ)=qHk(ϕ)+ξk(ϕ,q), (29)

    where ξk(ϕ,q)=q(1q)nj=1bjϕj(μj+pj(1K1j)+ϕj)(μj+pj(1K1j)+qϕj)Kkj. It follows from 0<q<1 and (ⅲ) of Assumption 1.1 that ξk is positive for all kN+ and ϕR+. Furthermore,

    H(qϕ)qH(ϕ)+ξ(q,ϕ)e,ξ=diag(ξ1,ξ2,,ξn). (30)

    Inequality (30), together with the monotonicity of H implies that

    ϕ=H(ϕ)qH(ϕ)+ξ(ˆϕ,q)e=qˆϕ+ξ(ˆϕ,q){H+(ˆϕ)}1H+(ˆϕ)eqˆϕ+ξ(ˆϕ,q){H+(ˆϕ)}1H(ˆϕ)=qˆϕ+ξ(ˆϕ,q){H+(ˆϕ)}1ˆϕ. (31)

    This contradicts the definition of q. Therefore, q1 and ϕqˆϕˆϕ. Exchanging the role of ϕ and ˆϕ, we can prove ϕˆϕ. This implies that ϕ=ˆϕ.

    In this section, we will show the relative compactness of the positive orbit {U(t,u0)}t0 defined by (5). This process is spurred by Lemma 19 in [28] and Theorem 2.46 in [29]. To apply them, we define

    ˜ϕk(t,a):={0,t>a,ik(t,a),ta,˜ψk(t,a):={0,t>a,vk(t,a),ta,

    and

    ˜ik(t,a)=ik(t,a)˜ϕk(t,a),˜vk(t,a)=vk(t,a)˜ψk(t,a).

    Then we can divide the solution of semigroup U(t) into two parts:

    V(t)u0=(0,0,˜ϕ(t,),0,˜ψ(t,)) (32)

    and

    W(t)u0=(S(t),0,˜i(t,),0,˜v(t,))T (33)

    where

    ˜ϕ=(˜ϕ1,˜ϕ2,,˜ϕn)T,˜ψ=(˜ψ1,˜ψ2,,˜ψn)T,˜i=(˜i1,˜i2,,˜in)T,˜v=(˜v1,˜v2,,˜vn)T.

    This implies that U(t)u0=V(t)u0+W(t)u0.

    Theorem 5.1. The semiflow U:R+×X0X0 is asymptotically smooth if there are maps V(t),W(t):R+×X0X0 such that U(t)u0=V(t)u0+W(t)u0, and the following statements hold for any bounded closed set Ω that is forward invariant under U:

    (i) For any u0Ω, there exists a function δ:R+×R+R+ such that for any r>0 limt+δ(t,r)=0 with u0Ωr, then V(t,u0)Ωδ(t,r);

    (ii) there exists a tΩ0 such that W(t)(Ω) has a compact closure for each ttΩ.

    Lemma 5.2. There exists a function δ:R+×R+R+ such that for any r>0,

    limt+δ(t,r)=0 (34)

    and

    V(t)u0Ωδ(t,r),u0Z,u0Ωr,t0. (35)

    Proof. Integrating ik and vk equations along the characteristic line ta= const. yields

    ˜ϕ(t,a)={0,t>a,ik0(at)πk(a)πk(at),ta,˜ψ(t,a)={0,t>a,vk0(at)π1k(a)π1k(at),ta

    for all kN. Hence for any u0Y and u0Yr,

    V(t)x0Ω=0+˜ϕ(t,)X+˜ψ(t,)X=nj=1{tij0(at)πj(a)πj(at)da+tvj0(at)π1j(a)π1j(at)da}eμ_tnj=10[ij0(a)+vj0(a)]da=eμ_t{i0X+v0X}eμ_tr

    where δ(t,r)=eμ_tr and μ_=minkN{μk}. Obviously, V(t)u0Ω approaches 0 as t goes to infinity.

    Lemma 5.3. W(t) maps any bounded subsets of Ω into sets with compact closure in Y.

    Proof. It follows from the first equation of (5) that Sk(t) remains in the compact set {ϕRn+|0ϕkbkμk,kN} for all t0. Therefore, we need to show that ˜i and ˜v still remain in pre-compact subsets of Y0 which is independent of u0. Note that

    ˜ik(t,a)={ik(ta,0)πk(a),t>a,0,ta,˜vk(t,a)={vk(ta,0)π1k(a),t>a,0,ta.

    Since for all kN+, 0˜ik(t,a)=ik(ta,0)πk(a)bkμ_2nj=1β+kjbjeμ_a:=Δkeμ_a where Δk=bkμ_2nj=1β+kjbj, and 0˜vk(t,a)bkpkμkeμ_a. It is easy to see that (ⅰ) - (ⅲ) of Theorem B.2 in [29] hold.

    In what follows, we need to show that (iv) of Theorem B.2 in [29] also holds. Assume that h(0,t) without loss of generality. Then

    nj=10|˜ij(t,a+h)˜ij(t,a)|da=nj=1tth|0˜ij(t,a)|da+nj=1th0|˜ij(t,a+h)˜ij(t,a)|danj=1Δjh+nj=1th0˜ij(t,a)|π(h)1|danj=1Δjh+nj=1Δj(th)h=nk=1Δjh(1+th). (36)

    Similarly,

    nj=10|˜vj(t,a+h)˜vj(t,a)|dank=1bjpjμjh(1+th). (37)

    Obviously, both (36) and (48) uniformly converge to 0 as h0 which is independent of u0.

    Lemmas 5.2 and 5.3, together with U(t)u0Z0nk=1bkμk, imply that U(t) has compact closure in Z for u0Ω. It follows from Proposition 3.13 in [32] that the solution orbit is relatively compact and the semiflow U(t) is asymptotically smooth.

    Proposition 5.1. The semiflow U(t) defined in (10) is asymptotically smooth.

    In this section, we establish the uniform persistence of (5) when R0>1. This property guarantees the well-definition of the Lyapunov functionals in Section 7. For some kN, define ρk:ΓR+ and u0=(S0,0,i0(),0,v0())Ω0 by

    ρk(S(t),0,i(t,),0,v(t,))=λk(t)=nj=10βkj(a)ij(t,a)daforu0Ω.

    Let

    Ω0={(S0,0,i0(),0,v0())Ω:ρk(U(t0,u0))>0for somet0R+}.

    Obviously, if u0ΩΩ0, then (S(t),0,i(t,),0,v(t,))E0 as t. Hence, if U has a global compact attractor in Ω0, then it also has a global compact attractor in Ω.

    Assumption 6.1. The support of the initial age-since-infection value ik0(a)L1(R+) lies to the right of the support the infectivity function βkj for all k,jN.

    Proposition 6.1. If R0>1 and Assumption 6.1 hold, then system (5) is uniformly weakly ρ-persistent for some kN.

    Proof. Since R0>1, there exists an ηj0>0 and ϕjR+ for k,jN such that

    nj=1˜S0j(ηj0)ˆKkj(ηj0)ϕj>ϕk, (38)

    where ˜S0k(ηk0)=bkμ+ηk0+pk(1K1k)ηk0(>0) and ^Kkj()=0eλaβkj(a)πk(a)da.

    Suppose that, for any ηk0>0, there exists an u0Ω0 such that

    lim suptρk(U(t,u0))ηk0

    and show a contradiction. Therefore, there exists a t0R+ such that ρk(U(t,u0))ηk0 for tt0 and all kN. Without loss of generality, we shift the time to t0=0. Then λk(t)ηk0 for tt0=0 and kN.

    Next, we show that SkS0k(ηk0):=bkμk+ηk0+pk(1K1k), where Sk=lim inftSk(t). By the Fluctuation Lemma [35], we can pick up a sequence {tn} such that Sk(tn)Sk, dSk(tn)dt0 as n. Then from the first equation of (5), it follows that

    dSk(tn)dtbk(μk+pk)Sk(tn)ηk0Sk(tn)+0εk(a)v(tn,a)da.

    This, combined with (40), gives

    dSk(tn)dtbk(μk+pk)Sk(tn)ηk0Sk(tn)+tn0εk(a)pkSk(tna)π1k(a)da.

    Letting n\to\infty leads to

    0\ge b_k-(\mu+\phi)S_{k{\infty}}-\eta_{k0}S_{k\infty}+p_kS_{k\infty} K_k^1.

    This implies that S_{k\infty}\ge S_k^0(\eta_{k0}).

    Finally, since S_{k\infty}\ge S_k^0(\eta_{k0}), there exists a t_1\in\mathbb{R}_+ such that S_k(t)\ge \tilde S_k^0(\eta_{k0}) for t\ge t_1. Again, without loss of generality, we can assume t_1 = 0. Solving i_k(t, a) by the characteristic line method yields

    i_k(t, a) = \left\{ \begin {array}{ll} b_k(t-a)\pi_k (a), \quad & t\ge a, \\ i_{k0}(a-t)\frac {\pi_k(a)}{\pi_k(a-t)}, \quad & t \lt a, \end {array} \right. (39)

    where b_k(t) = S_k(t)\lambda_k(t). Then

    \begin{eqnarray*} \lambda_k(t) & = &\sum\limits_{j = 1}^n\int_0^\infty\beta_{kj}(a)i_j(t, a)da \\ &\ge&\sum\limits_{j = 1}^n\int_0^t\beta_{kj}(a)S_j(t-a)\lambda_j(t-a)\pi_j (a)da \\ &\ge&\sum\limits_{j = 1}^n\tilde S_j(\eta_{j0})\int_0^t\beta_{kj}(a)\pi_j (a)\lambda_j(t-a)da. \end{eqnarray*}

    Taking Laplace transforms on both sides of the above inequality gives

    \widehat{\lambda_k(\xi)}\ge \sum\limits_{j = 1}^n\tilde S_j(\eta_{j0})\widehat{K_k(\xi)}\widehat{\lambda_j(\xi)}.

    This inequality holds for any \xi and \eta_{k0}. If we take \xi and \eta_{k0} small enough, then this is a contradiction with (38). This completes the proof.

    Clearly, for all k\in\mathbb N, \rho_k is a continuous function on \mathbb R_+. Proposition 5.1 implies that \{\mathcal U\}_{t\ge0} has a global attractor. From Theorem A.34 in [34], we need to show that for any bounded total orbit h(t+s) = \mathcal U(s, u(t)) of \mathcal U_t such that \rho_k(h(t))>0 for all t\in\mathbb R and s\in\mathbb R_+. For the total trajectory, we have

    i_k(t, a) = b_k(t-a)\pi_k(a), v_k(t, a) = pS_k(t-a)\pi_k^1(a)

    for t\in\mathbb R and a\in\mathbb R_+. In order to prove the strongly uniform persistence, the following lemma is helpful.

    Lemma 6.1. Let ({\bf S}(t), {\bf 0, i}(t, \cdot), {\bf0, v}(t, \cdot)) be a solution of (5). Then S_k^{\infty}\le S_k^0, where S_k^{\infty} = \limsup\limits_{t\rightarrow\infty}S_k(t) for k\in\mathbb N.

    Proof. By Fluctuate Lemma in [35], there exists \{t_n\} such that t_n\to\infty, S_k(t_n)\to S_k^{\infty}, and \frac{dS_k(t_n)}{dt}\to 0 as n\to\infty. Integrating v_k equation in (5) along the characteristic equation t-a = \mathrm{const.}, we have

    v_k(t, a) = \left\{ \begin{array}{ll} p_k S_k(t-a)\pi_k^1(a), \quad & t\ge a, \\ v_{k0}(a-t)\frac{\pi_k^1(a)}{\pi_k^1(a-t)}, & t \lt a, \end {array} \right. \textrm{for all} k\in\mathbb N. (40)

    Substituting v_k into S_k equation yields

    \begin {eqnarray*} \frac{dS_k(t_n)}{dt} & = & b_k-(\mu_k+p_k)S_k(t_n)-S_k(t_n)\lambda_k(t_n)+\int_0^{\infty}\varepsilon_k(a)v_k(t_n, a)da \\ & \le & b_k-(\mu_k+p_k)S_k(t_n) +\int_0^{t_n}\varepsilon_k(a) p_k S_k(t_n-a)\pi_k^1(a)da \\ &&+\int_{t_n}^{\infty}\varepsilon_k(a)v_{k0}(a-t_n)\frac{\pi_k^1(a)}{\pi_k^1(a-t_n)}da \\ & \le & b_k -(\mu_k+p_k)S_k(t_n)+\int_0^{t_n}\varepsilon_k(a) p_k S_k(t_n-a)\pi_k^1(a)da \\ &&+\int_{0}^{\infty}\varepsilon_k(a+t_n)v_{k0}(a)\frac{\pi_k^1(a+t_n)}{\pi_k^1(a)}da. \end {eqnarray*}

    Applying Fluctuate Lemma, we immediately obtain

    0\le b_k-(\mu_k+p_k)S_k^{\infty}+p_k S_k^{\infty}K_k^1

    or S_k^{\infty}\le \frac{b_k}{\mu_k+p_k(1-K_k^1)} = S_k^0 as required. This completes the proof.

    Lemma 6.2. Let us define \Theta(b(t)) = \sum\limits_{k = 1}^nb_k(t). Suppose that Assumption 6.1 holds, then \Theta(b(t)) for total trajectory h(t) is identically zero on \mathbb R, or it is strictly positive on \mathbb R.

    Proof. Suppose that there exists a t_1>0 such that b_k(t) = 0 for all t\le t_1. By the definition of \lambda_k, we have

    \begin{aligned} b_k(t) = &S_k(t)\sum\limits_{j = 1}^n\int_0^\infty\beta_{kj}(a)i_j(t, a)da \\ = &S_k(t)\sum\limits_{j = 1}^n\int_0^\infty\beta_{kj}(a)b(t-a)\pi_j(a)da\le S_k^0\bar\beta\int_0^\infty \Theta(b(t-a))da \\ = &S_k^0\bar\beta\int_0^{t-t_1} \Theta(b(t-a))da\le S_k^0\bar\beta\int_0^t\Theta(b(a))da, \end{aligned} (41)

    where \bar \beta = \max\limits_{j, k\in\mathbb N}\{\beta_{kj}^+\}. Summing k from 1 to n on both sides of (41) yields

    \Theta(b(t))\le \int_0^t\Theta(b(a))daS^0_{\bar\beta}, S_{\bar\beta}^0 = \sum\limits_{k = 1}^n\bar\beta S_k^0

    for all t> t_1. It follows from Gronwall inequality that \Theta(b(t)) = 0 for all t\ge t_1.

    Suppose there doesn't exist a t_1 such that \Theta(b(t)) = 0 for all t\le t_1. Thus, there exists a sequence \{t_m\} towards -\infty such that \Theta(b(t_m))>0 for each m. That means that b_{k}(t_m)>0 for each m. Moreover, there exists a sequence a_m such that i_k(t_m, a_m) = i_k(t_m-a_m, 0)\pi_k(a_m)>0 for each m. In view of the first equation, with the dissipative property of system (5), we obtain

    S_k'(t)\ge b_k-(\mu_k+p_k)S_k(t).

    Hence, there exists a positive constant \zeta>0 such that S_k(t)>\zeta>0 holds for all t\in\mathbb R. Let b_{km}(t) = b_k(t+t_m^*) for each n, where t_m^* = t_m-a_m. Recalling equation (41), we arrive at

    \begin{aligned} b_{km}(t)\ge&\zeta\left[\sum\limits_{j = 1}^n\int_0^t\beta_{kj}(a)b_{jm}(t-a)\pi_j(a)da+\tilde\Theta_{k m}(t)\right] \\ \ge& \zeta\left[\varepsilon_0\int_0^t \Theta(b_m(t-a))da+\tilde\Theta_{k m}(t)\right], \end{aligned} (42)

    where \varepsilon_0 is defined in (ⅲ) of Assumption 1.1 and

    \tilde \Theta_{k m}(t) = \sum\limits_{j = 1}^n\int_t^\infty \beta_{kj}(a) i_{jm0}(a-t)\frac{\pi_j(a)}{\pi_j(a-t)}da.

    Consequently,

    \Theta(b_m(t))\ge\zeta\left[\varepsilon_0n\int_{0}^t \Theta(b_m(a))da+\tilde\Theta(t)\right], \tilde\Theta_m(t) = \sum\limits_{k = 1}^n\tilde\Theta_{k m}(t). (43)

    Since \tilde\Theta_m(0) = \sum\limits_{k = 1}^n\sum\limits_{j = 1}^n\int_0^\infty \beta_{kj}(a) i_{jm0}(a)da>0 and \tilde\Theta_m(t) is continuous at t = 0, it follows from Assumption [29] and Gronwall inequality that \Theta(b_m(t)) and \tilde\Theta_m(t) are positive for sufficiently small t. Furthermore, from Corollary B.6 [29], we conclude that there exists a constant l>0 such that \Theta(b_m(t))>0 for all t>l. Since \Theta(b_m(t)) is a time shift of \Theta(b(t)) by t_m^* with t_m^*\rightarrow-\infty as m\rightarrow\infty, it follows that \Theta(b(t))>0 for all t\in\mathbb R. This completes the proof.

    Corollary 6.1. Suppose that Assumption 6.1 holds. Then for all k\in\mathbb N, b_k(t) for the total trajectory h(t) is strictly positive for every t\in\mathbb R.

    From Propositions 5.1 and 6.1, together with Corollary 6.1, we apply Theorem 5.2 in [29] to illustrate the \rho- strongly uniform persistence of system (5).

    Lemma 6.3. Suppose that \mathcal R_0>1 and the assumption of Corollary 6.1 hold. Then system (5) persists uniformly strongly, in this sense, there exists some \eta_{k0} such that

    \liminf\limits_{t\rightarrow\infty} \lambda_k(t)\ge \eta_{k0}

    for some k\in\mathbb N and \lambda_k(0)\neq0.

    Proof. From Corollary 6.1, we readily see that \inf \lambda_k(0)>0 for some k\in\mathbb N. Applying Theorem A.34 in [34], we conclude that there exists a constant \eta_{k0}>0 such that \liminf\limits_{t\rightarrow+\infty}\rho_k(\mathcal U(t, u_0))>\eta_{k0}.

    By Lemma 6.3, we present the following theorem to state the uniformly strong persistence of system (5).

    Theorem 6.4. Suppose that \mathcal R_0>1 and Assumption 6.1 hold. There exist some positive constants \eta_{k0}>0(k\in\mathbb N) such that for all t\in\mathbb R and a\in\mathbb R_+

    \liminf\limits_{t\rightarrow\infty}S_k(t)\ge\eta_{k0}, \liminf\limits_{t\rightarrow\infty}i_k(t, a)\ge\eta_{k0}\pi_k(a), \liminf\limits_{t\rightarrow\infty}v_k(t, a)\ge\eta_{k0}\pi_k^1(a).

    In this section, we will show the global behavior of the equilibria of system (5). To achieve this goal, we employ a Volterra type functional defined by g(x) = x-1-\ln x in [22], which is positive and attains minimum value 0 at x = 1. In what follows, we check this Volterra type functional is well-defined in infinite dimension and make the following assumption.

    Assumption 7.1. For all j\in\mathbb N_+, S_{j0}\in \mathbb R_+, \displaystyle\int_0^\infty|\ln h_{j0}(a)|da < +\infty, h = i, v.

    Lemma 7.1. If Assumption 7.1 holds, then \displaystyle\int_0^\infty v_j^0(a)\ln\frac{v_j(t, a)}{v_j^0(a)}da is bounded.

    Proof. For t>a,

    \begin{aligned} \left|v_j^0\ln\frac{v_j(t, a)}{v^0_j(a)}\right| = &|v_j^0(a)\ln v_j(t, a)-v_j^0(a)\ln v_j^0(a)|\\ \leq& |v_j^0(a)\ln v_j(t, a)|+|v_j^0(a)\ln v_j^0(a)|\\ = &|p_jS_j^0\pi_j^1(a)\ln p_jS_j^0(t-a)\pi_j^1(a)+|p_jS_j^0\pi_j^1(a)\ln p_jS_j^0\pi_j^1(a)|\\ \le& 2p_jS_j^0\ln \frac{p_j\Lambda_j}{\mu_j}e^{-\mu_ja}+2p_jS_j^0e^{-\mu_ja}\mu_ja. \end{aligned} (44)

    For t\le a,

    \begin{aligned} \left|v_j^0\ln\frac{v_j(t, a)}{v^0_j(a)}\right|\leq& |v_j^0(a)\ln v_j(t, a)|+|v_j^0(a)\ln v_j^0(a)|\\ = &|p_jS_j^0\pi_j^1(a)\ln v_{j0}(a-t)\pi_j^1(t)|+|p_jS_j^0\pi_j^1(a)\ln p_jS_j^0\pi_j^1(a)|\\ \le& p_jS_j^0|\ln v_{j0}(a-t)\pi_j^1(t)|+p_jS_j^0\ln \frac{p_j\Lambda_j}{\mu_j}e^{-\mu_ja} \\ &+2p_jS_j^0e^{-\mu_ja}\mu_ja. \end{aligned} (45)

    It follows from (44) - (45) that

    \begin{align} &\displaystyle\int_0^\infty\left|v_j^0(a)\ln\frac{v_j(t, a)}{v^0_j(a)}\right|da\notag \\ = &\displaystyle\int_0^\infty|v_j^0(a)\ln v_j(t, a)-v_j^0(a)\ln v_j^0(a)|da\notag\\ \leq& p_jS_j^0\ln \frac{p_j\Lambda_j}{\mu_j}\int_0^t e^{-\mu_ja}da+p_jS_j^0\displaystyle\int_t^\infty|\ln v_{j 0}(a-t)|e^{-\mu_jt}da\notag\\ &+p_jS_j^0\ln \frac{p_j\Lambda_j}{\mu_j}\int_0^\infty e^{-\mu_ja}da+2p_jS_j^0\int_0^\infty ae^{-\mu_ja}da\label{eqnE5.3}\\ = &p_jS_j^0\ln \frac{p_j\Lambda_j}{\mu_j}(1-e^{-\mu_jt})+p_jS_j^0 e^{-\mu_jt}\displaystyle\int_0^\infty|\ln v_{j 0}(a)|da\notag\\ &+\frac{p_jS_j^0}{\mu_j}\ln \frac{p_j\Lambda_j}{\mu_j}+\frac{2p_jS_j^0}{\mu_j}\notag\\ \leq&2p_jS_j^0\ln \frac{p_j\Lambda_j}{\mu_j}+p_jS_j^0 e^{-\mu_jt}\displaystyle\int_0^\infty|\ln v_{j 0}(a)|da +\frac{2p_jS_j^0}{\mu_j}.\notag \end{align} (46)

    Therefore, it follows from Assumption 7.1 that \displaystyle\int_0^\infty v_j^0(a)\ln\frac{v_j(t, a)}{v^0_j(a)}da is bounded.

    Theorem 7.2. Let Assumption 7.1 hold. If \mathcal R_0 = r(\mathcal K) < 1, the disease-free equilibrium E^0 is a global attractor in \Omega.

    Proof. For j\in\mathbb N, define

    V_j(t) = \sum\limits_{k = 1}^nK_{jk}S_k^0g\left(\frac{S_k(t)}{S_k^0}\right)+\sum\limits_{k = 1}^n\int_0^\infty\alpha_{jk}(a)i_k(t, a)da+\sum\limits_{k = 1}^nK_{j k}\int_0^\infty\delta_k(a)g\left(\frac{v_k(t, a)}{v^0_k(a)}\right)da,

    where \alpha_{kj}(a) = \int_a^\infty\beta_{kj}(s)\frac{\pi_j(s)}{\pi_j(a)}ds and \delta_k(a) = \int_a^\infty\varepsilon_k(s)v^0_k(s)da. Lemma 7.1 ensures V_j(t) is well-defined. Then

    \alpha'_{kj}(a) = -\beta_{kj}(a)+(\mu_k+\epsilon_k(a))\alpha_{kj}(a), \quad \delta'_k(a) = -\epsilon_k(a)i_k^0(a). (47)

    From Lemma 7.1, together with 6.4, V_j(t) is well-defined. Deviating it along the solution of (5) yields

    \begin{array}{rl} \frac{dV_j(t)}{dt}|_{(5)} = &\sum\limits_{k = 1}^n\left\{K_{jk}\left(1-\frac{S_k^0}{S_k(t)}\right)S_k'+\int_0^\infty\alpha_{jk}(a)\frac{\partial i_{jt}(t, a)}{\partial t}da+K_{j k}\int_0^\infty\delta_k(a)\frac{\partial v_{k}(t, a)}{\partial t}da\right\}\\ = &\sum\limits_{k = 1}^n\left[K_{j k}\left(1-\frac{S_k^0}{S_k(t)}\right)\left(b_k-(\mu_k+p_k)S_k-i_k(t, 0)+\int_0^\infty\epsilon_k(a)v_k(t, a)da\right)\right]\\ &+\sum\limits_{k = 1}^n\left[K_{j k}i_k(t, 0)-\sum\limits_{j = 1}^n\int_0^\infty\beta_{kj}(a)i_k(t, a)da\right]\\ &+\sum\limits_{k = 1}^n\left[K_{j k}\int_0^\infty\epsilon_k(a)\left(g\left(\frac{v_k(t, 0)}{v_k^0}\right) -g\left(\frac{v_k(t, a)}{v_k^0(a)}\right)\right)da\right]\\ = &-\sum\limits_{k = 1}^n\left[K_{j k}(\mu_k+p_k(1-K_k^1))\left(2-\frac{S_k(t)}{S_k^0}-\frac{S_k^0}{S_k(t)}\right)\right]\\ &-\sum\limits_{k = 1}^n\left[K_{j k}\int_0^\infty\epsilon_k(a)g\left(\frac{v_k(t, a)S_k^0}{v_k^0(a)S_k(t)}\right)da\right]\\ &+\sum\limits_{k = 1}^n(S_k^0K_{j k}-1)\int_0^\infty\beta_{jk}(a)i_k(t, a)da. \end{array} (48)

    Note that \mathcal R_0 = r(\mathcal K) < 1 implies that \sum\limits_{k = 1}^nS_k^0K_{jk} < 1. Therefore, V'_j(t)\le0 and it is easy to see that the equality holds if and only if ({\bf S}(t), {\bf 0}, {\bf i}(t, \cdot), {\bf 0}, {\bf v}(t, \cdot)) = ({\bf S^0}, {\bf 0}, {\bf 0}, {\bf 0}, {\bf v^0}(\cdot)). This implies that the largest positive invariant subset of \{u(t)\in \Omega|V'_j(t) = 0\} is the singleton \{({\bf S}^0, {\bf 0}, {\bf 0}, {\bf 0}, {\bf v^0}(\cdot))\}. This shows that the disease-free equilibrium E_0 is a global attractor.

    In the following, we give a lemma to show the boundedness of the Lyapunov functional for proving the global attractivity of the endemic equilibrium E^*.

    Lemma 7.3. If Assumption 7.1 holds, then \displaystyle\int_0^\infty h^*_j(a)\ln\frac{h_j(t, a)}{h^*_j(a)}da (h = i, v) is bounded.

    Proof. Prior to this proof, denote \beta_j^+ = \max\limits_{k\in\mathbb N}ess.\sup\limits_{a\in\mathbb R_+}\beta_{jk}(a) for all j\in\mathbb N. For t>a, it follows from b_j(t)\le\beta^+_j\frac{\Lambda_j^2}{\mu_j^2} that

    \begin{aligned} \left|i_j^*(a)\ln\frac{i_j(t, a)}{i^*_j(a)}\right| = &|i_j^*(a)\ln i_j(t, a)-i_j^*(a)\ln i_j^*(a)|\\ \leq& |i_j^*(a)\ln i_j(t, a)|+|i_j^*(a)\ln i_j^*(a)|\\ = &|i_j^*(0)\pi_j(a)\ln b_j(t-a)\pi_j(a)+|i_j^*(0)\pi_j(a)\ln i_j^*(0)\pi_j(a)|\\ \le& 4\left(\frac{\Lambda_j}{\mu_j}\right)^2 e^{-\mu_ja}\beta_j^+\ln\beta_j^+\frac{\Lambda_j}{\mu_j}+2\left(\frac{\Lambda_j}{\mu_j}\right)^2 e^{-\mu_ja}\mu_ja. \end{aligned} (49)

    For a>t,

    \begin{aligned} \left|i_j^*(a)\ln\frac{i_j(t, a)}{i^*_j(a)}\right| = &|i_j^*(a)\ln i_j(t, a)-i_j^*(a)\ln i_j^*(a)|\\ \leq& |i_j^*(a)\ln i_j(t, a)|+|i_j^*(a)\ln i_j^*(a)|\\ = &|i_j^*(0)\pi_j(a)\ln i_{j 0}(a-t)\pi_j(a)+|i_j^*(0)\pi_j(a)\ln i_j^*(a)\pi_j(a)|\\ \le&2\left(\frac{\Lambda_j}{\mu_j}\right)^2 e^{-\mu_ja}|\ln i_{j 0}(a-t)|+ 2\left(\frac{\Lambda_j}{\mu_j}\right)^2 e^{-\mu_ja}\beta_j^+\ln\beta_j^+\frac{\Lambda_j}{\mu_j}\\ &+2\left(\frac{\Lambda_j}{\mu_j}\right)^2 e^{-\mu_ja}\mu_ja. \end{aligned} (50)

    Then,

    \begin{align} &\displaystyle\int_0^\infty\left|i_j^*(a)\ln\frac{i_j(t, a)}{i^*_j(a)}\right|da\notag \\ = &\int_0^\infty|i_j^*(a)\ln i_j(t, a)-i_j^*(a)\ln i_j^*(a)|\notag\\ \leq& 2\frac{\Lambda_j^2}{\mu^3_j}\beta_j^+\ln\beta_j^+\frac{\Lambda_j}{\mu_j}\int_0^t e^{-\mu_ja}da+2\left(\frac{\Lambda_j}{\mu_j}\right)^2\displaystyle\int_t^\infty|\ln i_{j 0}(a-t)|e^{-\mu_jt}da\notag\\ &+2\frac{\Lambda_j^2}{\mu^3_j}\beta_j^+\ln\beta_j^+\frac{\Lambda_j}{\mu_j}\int_0^\infty e^{-\mu_ja}da+2\left(\frac{\Lambda_j}{\mu_j}\right)^2\mu_j\displaystyle\int_0^\infty ae^{-\mu_ja}da\notag\\ = &2\frac{\Lambda_j^2}{\mu^3_j}\beta_j^+\ln\beta_j^+\frac{\Lambda_j}{\mu_j}(1-e^{-\mu_jt})+2\left(\frac{\Lambda_j}{\mu_j}\right)^2\displaystyle\int_0^\infty|\ln i_{j 0}(a)|e^{-\mu_jt}da\label{eqnD5.2}\\ &+2\frac{\Lambda_j^2}{\mu^3_j}\beta_j^+\ln\beta_j^+\frac{\Lambda_j}{\mu_j}+2\frac{\Lambda_j^2}{\mu_j^3}\notag\\ \leq&4\frac{\Lambda_j^2}{\mu^3_j}\beta_j^+\ln\beta_j^+\frac{\Lambda_j}{\mu_j}+2\left(\frac{\Lambda_j}{\mu_j}\right)^2 e^{-\mu_jt}\displaystyle\int_0^\infty|\ln i_{j 0}(a)|da+\frac{2\Lambda_j^2}{\mu_j^3}.\notag \end{align} (51)

    By the assumption, it follows that \displaystyle\int_0^\infty i^*_j(a)\ln\frac{i_j(t, a)}{i^*_j(a)}da is bounded.

    Similarly, \displaystyle\int_0^\infty v^*_j(a)\ln\frac{v_j(t, a)}{v^*_j(a)}da is also bounded.

    Assume that f(t, a)\in\mathbb R\times L^1(\mathbb R_+) is a solution of the following system

    \begin{aligned} \displaystyle\frac{\partial f(t, a)}{\partial t}+\frac{\partial f(t, a)}{\partial a} = &-m(a)f(t, a), \\ f(t, 0) = &L(t)\displaystyle\int_0^\infty\eta(a)f(t, a)da, \\ f(0, a) = &f_0(a), \end{aligned} (52)

    where m(a), \eta(a)\in L^1(\mathbb R_+) and L(t)\in\mathbb R. Obviously, the nontrivial equilibrium E_1^* of system (52) satisfies the following equations:

    \begin{aligned} \displaystyle\frac{d f^*(a)}{d a} = &-m(a)f^*(a), \\ f^*(0) = &L^*\displaystyle\int_0^\infty\eta(a)f^*(a)da. \end{aligned} (53)

    Define a Lyapunov functional V_1(t) = \int_0^\infty M(a)g\left(\frac{f(t, a)}{f^*(a)}\right)da, where M(a) = \int_a^\infty\beta(s) f^*(s)ds.

    Lemma 7.4. Suppose that \displaystyle\int_0^\infty|\ln f_0(a)|da is bounded. There exists a positive value \eta_0>0 such that f(t, a)>\eta_0e^{-\int_0^am(s)ds}, then

    \displaystyle\frac{dV_1(t)}{dt} = \int_0^\infty M(a)f^*(a)\left[g\left(\frac{f(t, 0)}{f^*(0)}\right)-g\left(\frac{f(t, a)}{f^*(a)}\right)\right]da, (54)

    where M(a) = \displaystyle\int_a^\infty\eta(s)f^*(s)da.

    Proof. Note that

    \begin{aligned} \frac{\partial g\left(\frac{f(t, a)}{f^*(a)}\right)}{\partial a} = &\displaystyle\frac{\partial }{\partial a}\left(\frac{f(t, a)}{f^*(a)}-1-\ln \frac{f(t, a)}{f^*(a)}\right)\\ = &\displaystyle\frac{\partial }{\partial a}\frac{f(t, a)}{f^*(a)}-\frac{\partial }{\partial a}\ln \frac{f(t, a)}{f^*(a)}\\ = &\displaystyle\frac{f'_a(t, a)f^*(a)-f(t, a)f^*_a(a)}{(f^*(a))^2}-\frac{f^*(a)}{f(t, a)}\frac{f'_a(t, a)f^*(a)-f(t, a)f^{*'}_a(a)}{(f^*(a))^2}\\ = &\displaystyle\left(\frac{1}{f^*(a)}-\frac{1}{f(t, a)}\right)\left(f'_a(t, a)-f(t, a)\right)\frac{f_a^{*'}(a)}{f^*(a)}\\ = &\displaystyle\left(\frac{1}{f^*(a)}-\frac{1}{f(t, a)}\right)f'_a(t, a)+\left(\frac{1}{f^*(a)}-\frac{1}{f(t, a)}\right)m(a)f(t, a), \end{aligned} (55)

    where we denote h'_l(s, l) = \frac{\partial h(s, l)}{\partial l}.

    It follows from the proving process of Lemma 7.3, together the assumption of Lemma 7.4 that V_1(t) is well-defined. Deviating it along the solution of (52), we obtain

    \begin{aligned} \frac{dV_1(t)}{dt}|_{(52)} = &\displaystyle\int_0^\infty M(a)\frac{\partial g\left(\frac{f(t, a)}{f^*(a)}\right)}{\partial t}da\\ = &\displaystyle\int_0^\infty M(a)\left(1-\frac{f^*(a)}{f(t, a)}\right)\frac{f'_t(t, a)}{f^*(a)}da\\ = &\displaystyle\int_0^\infty M(a)\left(\frac{1}{f^*(a)}-\frac{1}{f(t, a)}\right)(-f'_a(t, a)-m(a)f(t, a))da\\ = &-\displaystyle\int_0^\infty M(a)\left(\frac{1}{f^*(a)}-\frac{1}{f(t, a)}\right)f'_a(t, a)da \\ &-\int_0^\infty M(a)m(a)f(t, a)da. \end{aligned} (56)

    Observing (55), we have

    \frac{dV_1(t)}{dt}\Big|_{(52)} = -\int_0^\infty M(a)\frac{\partial g\left(\frac{f(t, a)}{f^*(a)}\right)}{\partial a}da. (57)

    With the help of integral by parts, we obtain

    \begin{aligned} \displaystyle\frac{dV_1(t)}{dt}\Big|_{(52)} = &-M(a)g\left(\displaystyle\frac{f(t, a)}{f^*(a)}\right)\Big|_0^\infty+\int_0^\infty M'_a(a)g\left(\frac{f(t, a)}{f^*(a)}\right)da\\ = &M(0)g\left(\displaystyle\frac{f(t, 0)}{f^*(0)}\right)-\int_0^\infty\eta(a)f^*(a)g\left(\frac{f(t, a)}{f^*(a)}\right)da\\ = &\displaystyle\int_0^\infty\eta(a)f^*(a)\left[g\left(\frac{f(t, 0)}{f^*(0)}\right)-g\left(\frac{f(t, a)}{f^*(a)}\right)\right]da, \end{aligned} (58)

    here we used the fact M(0) = \int_0^\infty\eta(a)f^*(a)da, and M'(a) = -\eta(a)f^*(a).

    Next, we will give the global attractivity of the endemic equilibrium E^*.

    Theorem 7.5. Suppose \mathcal R_0>1, (iv) of Assumption 1.1 and Assumption 7.1 hold. Then the endemic equilibrium E^* is a global attractor in \Omega_0.

    Proof. Define

    V(t) = \sum\limits_{j = 1}^n\kappa_j\left\{S_j^*g\left(\frac{S_j(t)}{S_j^*}\right)+\int_0^\infty\alpha_j(a)g\left(\frac{i_j(t, a)}{i_j^*(a)}\right)da +\int_0^\infty\delta_j(a)g\left(\frac{v_j(t, a)}{v_j^*(a)}\right)da\right\},

    where \alpha_j(a) = \sum\limits_{k = 1}^n\int_a^\infty\beta_{jk}(s)i^*_k(s)da and \delta_j(a) = \int_a^\infty\varepsilon_j(s)v^*_j(s)ds. \kappa_j will be determined later and considered as a weighted coefficient. The well-definition of V(t) follows from Lemma 7.3. From the definitions of \alpha_j(a) and \delta_j(a), it follows that

    \alpha_j'(a) = -\sum\limits_{k = 1}^n\beta_{jk}(a)i^*_k(a),

    and

    \delta_j'(a) = -\varepsilon_j(a)v^*_j(a).

    Assisting with Lemma 7.4 and deviating V(t) along the solution of (5) yield

    \begin{align} \displaystyle\frac{dV(t)}{dt}\Big|_{(5)} = & \displaystyle\sum\limits_{j = 1}^n\kappa_j\left\{\left(1-\frac{S_j^*}{S_j(t)}\right)S_j'(t)\right.\notag \\ &+\sum\limits_{k = 1}^n\int_0^\infty\beta_{jk}(a)i_k^*(a)\left[g\left(\frac{i_j(t, 0)}{i^*_j(0)}\right)-g\left(\frac{i_j(t, a)}{i^*_j(a)}\right)\right]da\notag\\ &\left.+\displaystyle\int_0^\infty\epsilon_j(a)v_j^*(a)\left[g\left(\frac{v_j(t, 0)}{v^*_j(0)}\right)-g\left(\frac{v_j(t, a)}{v^*_j(a)}\right)\right]da\right\}\notag\\ = &\displaystyle\sum\limits_{j = 1}^n\kappa_j\left\{(1-\frac{S_j^*}{S_j(t)})\left[-(\mu_j+p_j)(S_j(t)-S_j^*)\right.\right.\notag\\ &-\left(S_j(t)\displaystyle\sum\limits_{k = 1}^n\int_0^\infty\beta_{jk}(a)i_k(t, a)da-S_j^*\sum\limits_{k = 1}^n\int_0^\infty\beta_{jk}(a)i_k^*(a)da\right)\notag\\ &\left.+\displaystyle\int_0^\infty\epsilon_j(a)(v_j(t, a)-v_j^*(a))da\right]\notag\\ &+\displaystyle\sum\limits_{k = 1}^n\int_0^\infty\beta_{jk}(a)i_k^*(a)\left[g\left(\frac{i_j(t, 0)}{i^*_j(0)}\right)-g\left(\frac{i_j(t, a)}{i^*_j(a)}\right)\right]da\label{eqn5.14}\\ &\left.+\displaystyle\int_0^\infty\epsilon_j(a)v_j^*(a)\left[g\left(\frac{v_j(t, 0)}{v^*_j(0)}\right)-g\left(\frac{v_j(t, a)}{v^*_j(a)}\right)\right]da\right\}\notag\\ = &\displaystyle\sum\limits_{j = 1}^n\kappa_j\left\{-(\mu_k+p_k)\left(2-\frac{S_k(t)}{S_k^*}-\frac{S_k^*}{S_k(t)}\right)+S_j^*\displaystyle\sum\limits_{k = 1}^n\int_0^\infty\beta_{jk}(a)i_k^*(a)\right.\notag\\ &\displaystyle\times\left[\frac{i_j(t, 0)}{i^*_j(0)}-\frac{i_j(t, a)}{i^*_j(a)}-\frac{i_k(t, a)S_j(t)}{i^*_k(a)S_j^*}+\frac{i_k(t, a)}{i^*_k(a)} -\frac{S_j^*}{S_j(t)}\right]da\notag\\ &\left.+\displaystyle\int_0^\infty\epsilon_j(a)v^*_j(a)\left[\frac{S_j^*}{S_j(t)}+\frac{S_j(t)}{S_J^*}-1-\frac{S_j^*v_j(t, a)}{S_j(t)v_j^*(a)}+\ln\frac{S_j^*v_j(t, a)}{S_j(t)v_j^*(a)}\right]da\right\}.\notag \end{align} (59)

    Note that

    \begin{align} S_j^*\displaystyle\sum\limits_{k = 1}^n\int_0^\infty\beta_{jk}(a)i_k^*(a)\frac{i_k(t, a)S_j(t)}{i^*_k(a)S_j^*}da = &S_j(t)\displaystyle\sum\limits_{k = 1}^n\int_0^\infty\beta_{jk}(a)i_k(t, a)da\notag\\ = &i_j(t, 0)\label{eqna5.15}\\ = &S_j^*\displaystyle\sum\limits_{k = 1}^n\int_0^\infty\beta_{jk}(a)i_k^*(a)\frac{i_j(t, 0)}{i^*_j(0)}da.\notag \end{align} (60)

    Equation (60) implies that

    S_j^*\displaystyle\sum\limits_{k = 1}^n\int_0^\infty\beta_{jk}(a)i_k^*(a)\frac{i_k(t, a)S_j(t)i_j^*(0)}{i^*_k(a)S_j^*i_j(t, 0)}da = S_j^*\displaystyle\sum\limits_{k = 1}^n\int_0^\infty\beta_{jk}(a)i_k^*(a)da. (61)

    Noting that \int_0^\infty\epsilon_k(a)v_k^*(a)da = p_kS_k^*\int_0^\infty\epsilon_k(a)\pi_k^1(a)da = p_kS_k^*K_k^1, it follows from (60) and (61) that

    \begin{array}{rl} \displaystyle\frac{dV(t)}{dt}\Big|_{(5)} = &\displaystyle\sum\limits_{j = 1}^n\kappa_j\left\{-(\mu_j+p_j(1-K_j^1))\left(2-\frac{S_j(t)}{S_j^*}-\frac{S_j^*}{S_j(t)}\right)\right.\notag\\ &-S_j^*\displaystyle\sum\limits_{k = 1}^n\int_0^\infty\beta_{jk}(a)i_k^*(a)\left[g\left(\frac{S_j^*}{S_j(t)}\right)+g\left(\frac{i_k(t, a)S_j(t)i_j^*(0)}{i^*_k(a)S_j^*i_j(t, 0)}\right)\right]da\notag\\ &+S_j^*\displaystyle\sum\limits_{k = 1}^n\int_0^\infty\beta_{jk}(a)i_k^*(a)[H_k(i(t, a))-H_j(i(t, a))]da\notag\\ &\left.-\displaystyle\int_0^\infty\epsilon_j(a)v^*_j(a)g\left(\frac{S_j^*v_j(t, a)}{S_j(t)v_j^*(a)}\right)da\right\}, \end{array} (62)

    where H_j(i) = \displaystyle\frac{i_j(t, a)}{i^*_j(a)}-\ln\frac{i_j(t, a)}{i^*_j(a)}. Define \Theta_{jk} = S_j^*\int_0^\infty\beta_{jk}(a)i_k^*(a)[H_k(i)-H_j(i)]da, (j, k = 1, 2, \cdots, n) and a Laplacian matrix

    \Theta = \left( \begin{array}{cccc} \displaystyle\sum\limits_{k\neq1}\theta_{1k} & -\theta_{21} & \cdots & -\theta_{n1} \\ -\theta_{12} & \sum\limits_{k\neq2}\theta_{2k} & \cdots & -\theta_{n2} \\ \vdots & \vdots & \ddots & \vdots \\ -\theta_{1n} & -\theta_{2n} & \cdots & \sum\limits_{k\neq n}\theta_{nk} \\ \end{array} \right).

    By Lemma 2.1 in [21], we have that the solution space of the linear system \Theta\kappa = 0 is 1 and one of its basis is given by

    \kappa = (\kappa_1, \kappa_2, \cdots, \kappa_n)^T = (c_{11}, c_{22}, \cdots, c_{nn})^T,

    where c_{jj}>0 (j = 1, 2, \cdots, n) denotes the cofactor of the j-th diagonal element of matrix \Theta. This implies that \sum\limits_{k = 1}^n\theta_{jk}\kappa_j = \sum\limits_{j = 1}^n\theta_{kj}\kappa_k and

    \begin{align}\label{eqn5.15} &\displaystyle\sum\limits_{j = 1}^n\kappa_j\sum\limits_{k = 1}^n\int_0^\infty\beta_{jk}(a)i_k^*(a)[H_k(i(t, a))-H_j(i(t, a))]da\notag\\ & = \displaystyle\sum\limits_{j = 1}^n\sum\limits_{k = 1}^n\kappa_k\int_0^\infty\beta_{kj}(a)i_j^*(a)[H_j(i(t, a))-H_k(i(t, a))]da. \end{align} (63)

    Employing the graph-theoretic approach mentioned in [21], (63) is equal to zero. Therefore, V'(t)\le0. The equality holds if and only if \frac{i_k(t, a)}{i_k^*(a)} = \frac{i_j(t, 0)}{i_j^*(0)} = \frac{i_j(t, a)}{i_j^*(a)}. This implies that \mathbf i = \mathbf ci^*. It follows from the first equation of (5) with respect to the monotonicity of \mathbf c that \mathbf c = \mathbf 1. Hence, the largest invariant set of \{u(t)\in \Omega_0| V'(t) = 0\} is the singleton E^*. Combining the relative compactness of the solution orbit (see Lemma 5.3) with the invariance principle (see [Theorem 4.2, [30]]), we see that the endemic equilibrium E^* is a global attractor in \Omega_0.

    In this section, we perform some numerical experiments to illustrate our theoretical results. For the experimental operability, we set a_{max} = 10 instead of infinity. For simplicity, we assume system (5) describes some sexually transmission diseases, such as Zika, Ebola and genital warts etc, which consist of two groups - male group and female group. As some news reported in [33], human papillomaviruses (HPV) is an effective and safe vaccine to control some sexually transmitted diseases inducing by virus. In order to illustrate system (5), we firstly fix the demographic parameters as follows:

    b_1 = 20.5, b_2 = 40.5, \mu_1 = \mu_2 = 0.01.

    Hence the total male population maintains the size N_1(t) = 2050 and the size of the total female population is N_2(t) = 4050. It has been reported that nearly half of newly infections are diagnosed as females in age range between 15 to 24. Therefore, we assume that the vaccination rates p_1 = 0.21 for male and p_2 = 0.22 for female. Other parameters are associated with infection age in the form of

    y(x, A, B) = \frac{1}{B^A\Gamma(A)}x^{A-1}e^{-\frac{x}{B}}

    where \Gamma(n) = (n-1)!, or \Gamma(z) = \int_0^\infty x^{z-1}e^{-x}dx (z\in\mathbb N). Obviously, y(x, A, B) satisfies (ⅲ)' of Assumption 1.1. We fix \gamma(a) = y(a, 2, 5), \varepsilon(a) = y(a, 1, 1) and verify the transmission rate \beta(a) = y(a, 6, B). If we choose B = 9, this implies that the mean of the transmission function is 54 and the variance is 486. It follows from Figure 1 that the disease-free equilibrium E_0 is asymptotically stable. Then we decrease B = 8 associated with the variance changed as 384. Figure 2 shows that the endemic equilibrium E^* is asymptotically stable.

    Figure 1.  The solution of (5) with initial conditions S_{j0} = 200, i_{j0}(a) = 10, V_0 = 10 and all the parameters except B = 9 are enclosed in the text.
    Figure 2.  The solution of (5) with initial conditions S_{j0} = 200, i_{j0}(a) = 10, V_0 = 10 and all the parameters except B = 8 are enclosed in the text.

    Comparing infected male and female populations in Figure 2, we find that the total infected number for female population is larger than that for male population. However, the first peak time for female population is later than the time for male population. Although the vaccination rate for female designed is higher than such rate for male, the level of infected peak and the total number for female population are still larger than those for male population. Consequently, the government should pay more attentions to female population.

    In this paper, we proposed a multi-group SIVS epidemic model with age-since-infection. We calculated the basic reproduction number \mathcal R_0 by the renewal equation, which is the spectral radius of the next generation operator \mathcal K. From Theorems 7.2 and 7.5, we see that the global attractivity of system (5) is totally determined by \mathcal R_0. This implies that the basic reproduction number is a sharp threshold determining that the disease prevails or vanishes. Implicit Euler method performed illustrates the theoretical results. Numerical experiment shows that the number of the infected females is larger than the number of the infected male although the vaccination rate for female group is higher than that for male group. The government should pay much more concerns on female group for suppressing sexual diseases prevalence.

    Irreducibility of the transmission matrix \beta_{jk}, j, k\in\mathbb N has been highly influential in analyzing the global stability of equilibria. This assumption is still a basic assumption of epidemics spreading on scale free networks. Therefore, we hope our analysis method proving theoretical results and numerical method can be generalized to investigate the dynamics of some epidemic models on complex networks [36]. Besides, oscillation is one of the important phenomena in diseases transmission. What mechanisms resulting in oscillation has been becoming an increasing trend in investigating epidemic models with age-since-infection [37]. We leave this for our future work.

    The authors are very grateful to the Editor-in-Chief, Professor Yang Kuang, and the anonymous referees for their valuable comments and suggestions, which helped us to improve the presentation of this work significantly.

    This work was supported by the National Natural Science Foundation of China (Nos. 61573016, 11871316), Shanxi Scholarship Council of China (2015-094), Shanxi Scientific Data Sharing Platform for Animal Diseases, and Startup foundation for High-level Personal of Shanxi University.

    The authors declare there is no conflict of interest.



    [1] K. Takahashi and Y. Konishi, Non-linear vibrations of cables in three dimensions, Part II: out-of-plane vibrations under in-plane sinusoidally time-varying load, J. Sound Vib., 118 (1987), 85–97.
    [2] J. L. Lilien and A. P. Da Costa, Vibration amplitudes caused by parametric excitation of cable stayed structures, J. Sound Vib., 174 (1994), 69–90.
    [3] H. N. Arafat and A. H. Nayfeh, Non-linear responses of suspended cables to primary resonance excitations, J. Sound Vib., 266 (2003), 325–354.
    [4] H. Wang, T. Y. Tao, R. Zhou, et al., Parameter sensitivity study on flutter stability of a long-span triple-tower suspension bridge, J. Wind Eng. Ind. Aerod., 128 (2014), 12–21.
    [5] T. Grigorjeva and Z. Kamaitis, Numerical analysis of the effects of the bending stiffness of the cable and the mass of structural members on free vibrations of suspension bridges, J. Civ. Eng. Manag., 21 (2015), 948–957.
    [6] C. G. Koh and Y. Rong, Dynamic analysis of large displacement cable motion with experimental verification, J. Sound Vib., 272 (2004), 183–206.
    [7] M. H. EI Ouni and N. B. Kahla, Nonlinear dynamic analysis of a cable under first and second order parametric excitations, J. Civ. Eng. Manag., 18 (2012), 557–567.
    [8] D. Xue, Y. Liu, J. He, et al., Experimental study and numerical analysis of a composite truss joint, J. Constr Steel Res., 67 (2011), 957–964.
    [9] K. Takahashi, Dynamic stability of cables subjected to an axial periodic load, J. Sound Vib., 144 (1991), 323–330.
    [10] A. P. Da Costa and J. A. C. Martins, Oscillations of bridge stay cables induced by periodic motions of deck and/or towers, J. Eng. Mech., 122 (1996), 613–622.
    [11] H. Chen, D. Zuo, Z. Zhang, et al., Bifurcations and chaotic dynamics in suspended cables under simultaneous parametric and external excitations, Nonlinear Dyn., 62 (2010), 623–646.
    [12] D. Zulli and A. Luongo, Nonlinear energy sink to control vibrations of an internally nonresonant elastic string, Meccanica, 50 (2014), 781–794.
    [13] B. N. Sun and Z. G. Wang, Cable oscillation induced by parametric excitation in cable-stayed bridges, J. Zhejiang Univ-sc A, 4 (2003), 13–20.
    [14] Y. Fujino and P. Warnitchai, An experimental and analytical study of auto parametric resonance in a 3-DOF model of cable-stayed-beam, Nonlinear Dyn., 4 (1993), 111–138.
    [15] Y. Xia and Y. Fujino, Auto-parametric vibration of a cable-stayed-beam structure under random excitation, J. Eng. Mech., 132 (2006), 279–286.
    [16] Y. Xia, Q. X. Wu, Y. L. Xu, et al., Verification of a cable element for cable parametric vibration of one-cable-beam system subject to harmonic excitation and random excitation, Adv. Struct. Eng., 14 (2011), 589–595.
    [17] V. Gattulli and M. Lepidi, Nonlinear interactions in the planar dynamics of cable-stayed beam, Int J. Solids Struct., 40 (2003), 4729–4748.
    [18] V. Gattulli and M. Morandini, A parametric analytical model for non‐linear dynamics in cable-stayed beam, Earthq. Eng. Struct. D., 31 (2002), 1281–1300.
    [19] V. Gattulli and M. Lepidi, One-to-two global-local interaction in a cable-stayed beam observed through analytical, finite element and experimental models, Int. J. Nonlin Mech., 40 ( 2005), 571–588.
    [20] Z. Wang, C. Sun, Y. Zhao, et. al., Modeling and nonlinear modal characteristics of the cable-stayed beam, Eur. J. Mech- A-Solid, 47 (2014), 58–69.
    [21] M. H. Wei, K. Lin, L. Jin, et. al., Nonlinear dynamics of a cable-stayed beam driven by sub-harmonic and principal parametric resonance, Int. J. Mech. Sci., 110 (2016), 78–93.
    [22] Z. W. Wang and T. J. Li, Nonlinear dynamic analysis of parametrically excited space cable-beam structures due to thermal loads, Eng. Struct., 83 (2015), 50–61.
    [23] H. J. Kang, H. P. Zhu, Y. Y. Zhao, et al., In-plane non-linear dynamics of the stay cables, Nonlinear Dyn., 73 (20136), 54–59.
    [24] Y. Y. Zhao and H. J. Kang, In-plane free vibration analysis of cable-arch structure, J. Sound Vib., 312 (2008), 363–379.
    [25] H. J. Kang, Y. Y. Zhao and H. P. Zhu, Out-of-plane free vibration analysis of a cable-arch structure, J. Sound Vib., 332 (2013), 907–921.
    [26] A. H. Nayfeh and D. T. Mook, Nonlinear Oscillations, John Wiley & Sons, (1979), 71–79.
    [27] J. R. Zhao, Y. G. Tang and L. Q. Liu, Study on the sum type combination resonance responses of a classic Spar platform, Ocean Eng., 27 (2009), 23–30.
    [28] L. Q. Chen, J. W. Zu and J. Wu, Steady-state response of the parametrically excited axially moving string constituted by the Boltzmann superposition principle, Acta Mech., 162 (2003), 143–155.
    [29] Q. X. Wu, L. Liu and B. C. Chen, Theoretical equations of in-plane natural vibration for cables considering bending stiffness, Eng. Mech., 27 (2010), 9–27.
    [30] M. Irvine, Cable Structures, The MIT Series in Structural Mechanics, (1981).
    [31] Japan Society of Civil Engineers, Basic and application of cable and space structures, Maruzen, Japan: JSCE, (1999).
    [32] G. H. Li, H. F. Xiang and Z. Y. Shen, Stability and vibration of bridge structures, China Railway Publishing House, (2010), 132–134.
    [33] N. J. Gimsing and C. T. Georgakis, Cable supported bridges: concept and design. 3th ed, John Wiley & Sons, (2011), 220–225.
    [34] A. H. Nayfeh and D. T. Mook, Nonlinear Oscillation, John Wiley & Sons, (1979).
    [35] M. Al-Qassab and S. Nair, Wavelet-Galerkin method for the free vibrations of an elastic cable carrying an attached mass, J. Sound Vib., 270 (2004), 191–206.
    [36] S. H. Zhou, G. Q. Song, Z. H. Ren, et al., Nonlinear analysis of a parametrically excited beam with intermediate support by using Multi-dimensional incremental harmonic balance method, Chaos, Solitons & Fractals, 93 (2016), 207–222.
    [37] P. J. Li, R. H. Wang and N. J. Ma, Analytic algorithm of transverse vibration frequency of cables with intermediate flexible supports, J. South China Univ T., 40 (2012), 7–18.
  • This article has been cited by:

    1. Wei Lv, Hanfei He, Kezan Li, Robust optimal control of a network-based SIVS epidemic model with time delay, 2022, 161, 09600779, 112378, 10.1016/j.chaos.2022.112378
    2. Jing Yang, Shaojuan Ma, Juan Ma, Jinhua Ran, Xinyu Bai, Stochastic Analysis for the Dual Virus Parallel Transmission Model with Immunity Delay, 2024, 1557-8666, 10.1089/cmb.2024.0662
  • Reader Comments
  • © 2019 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)
通讯作者: 陈斌, bchen63@163.com
  • 1. 

    沈阳化工大学材料科学与工程学院 沈阳 110142

  1. 本站搜索
  2. 百度学术搜索
  3. 万方数据库搜索
  4. CNKI搜索

Metrics

Article views(3953) PDF downloads(361) Cited by(1)

Figures and Tables

Figures(10)  /  Tables(3)

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog