Loading [MathJax]/jax/element/mml/optable/MathOperators.js
Research article Topical Sections

Study of protective hard coatings of SiO2-TiO2 on aluminum substrates

  • Aluminum alloys are frequently employed in the aeronautics industry due to the remarkable mechanical properties and lightweight nature of these materials. Moreover, thin film coatings are commonly applied in order to improve the corrosion resistance under harsh environments. In this work, Al 7075-T6 substrates were coated with nanostructured SiO2-TiO2 films using a sol-gel method. The experimental approach initially consisted in the preparation of a precursor agent using tetraethyl orthosilicate (TEOS) and triethoxy(octyl)silane (ETOS). Subsequently, nanoparticles of SiO2-TiO2 were mixed in order to develop thin films using a one-step dip coating method. The roughness, nanoindentation and corrosion properties were evaluated for the coated substrates. A finite element model was created for the nanoindentation test, which determined the mechanical response between the film-contact interface during loading conditions. The average hardness, elastic modulus and critical loads leading to fracture were verified. The nanoindentation test presented a significant increase in hardness for the coated Al 7075-T6 alloy, reaching a value of 4.6 GPa. The SiO2-TiO2 thin films presented uniform and compact surface coatings with high mechanical properties. Furthermore, the performed corrosion tests indicated moderate protection by the SiO2-TiO2 thin films. The SiO2-TiO2 thin films displayed a generalized corrosion throughout the surface, presenting oxides and fractured crystals in localized regions.

    Citation: Johana Gamez, Luis Reyes-Osorio, Oscar Zapata, Roberto Cabriales, Luis Lopez, Miguel Delgado-Pamanes. Study of protective hard coatings of SiO2-TiO2 on aluminum substrates[J]. AIMS Materials Science, 2024, 11(2): 200-215. doi: 10.3934/matersci.2024011

    Related Papers:

    [1] Wahidullah Hussainzada, Han Soo Lee . Impact of land surface model schemes in snow-dominated arid and semiarid watersheds using the WRF-hydro modeling systems. AIMS Geosciences, 2024, 10(2): 312-332. doi: 10.3934/geosci.2024018
    [2] Firoz Ahmad, Laxmi Goparaju . Soil and Water Conservation Prioritization Using Geospatial Technology – a Case Study of Part of Subarnarekha Basin, Jharkhand, India. AIMS Geosciences, 2017, 3(3): 375-395. doi: 10.3934/geosci.2017.3.375
    [3] Kyle Whalley, Wei Luo . The Pressure of Society on Water Quality: A Land Use Impact Study of Lake Ripley in Oakland, Wisconsin. AIMS Geosciences, 2017, 3(1): 14-36. doi: 10.3934/geosci.2017.1.14
    [4] Maryam Khal, Abdellah Algouti, Ahmed Algouti, Nadia Akdim, Sergey A. Stankevich, Massimo Menenti . Evaluation of open Digital Elevation Models: estimation of topographic indices relevant to erosion risk in the Wadi M’Goun watershed, Morocco. AIMS Geosciences, 2020, 6(2): 231-257. doi: 10.3934/geosci.2020014
    [5] Quentin Fiacre Togbévi, Luc Ollivier Sintondji . Hydrological response to land use and land cover changes in a tropical West African catchment (Couffo, Benin). AIMS Geosciences, 2021, 7(3): 338-354. doi: 10.3934/geosci.2021021
    [6] Shailesh Kumar Singh, Nelly Marcy . Comparison of Simple and Complex Hydrological Models for Predicting Catchment Discharge Under Climate Change. AIMS Geosciences, 2017, 3(3): 467-497. doi: 10.3934/geosci.2017.3.467
    [7] Chibuike Chiedozie Ibebuchi, Itohan-Osa Abu . Relationship between synoptic circulations and the spatial distributions of rainfall in Zimbabwe. AIMS Geosciences, 2023, 9(1): 1-15. doi: 10.3934/geosci.2023001
    [8] Margherita Bufalini, Farabollini Piero, Fuffa Emy, Materazzi Marco, Pambianchi Gilberto, Tromboni Michele . The significance of recent and short pluviometric time series for the assessment of flood hazard in the context of climate change: examples from some sample basins of the Adriatic Central Italy. AIMS Geosciences, 2019, 5(3): 568-590. doi: 10.3934/geosci.2019.3.568
    [9] David W. Schwartzman . Life’s Critical Role in the Long-term Carbon Cycle: the Biotic Enhancement of Weathering. AIMS Geosciences, 2017, 3(2): 216-238. doi: 10.3934/geosci.2017.2.216
    [10] Ana Casado, Natalia C López . Comparison of synthetic unit hydrograph methods for flood assessment in a dryland, poorly gauged basin (Napostá Grande, Argentina). AIMS Geosciences, 2025, 11(1): 27-46. doi: 10.3934/geosci.2025003
  • Aluminum alloys are frequently employed in the aeronautics industry due to the remarkable mechanical properties and lightweight nature of these materials. Moreover, thin film coatings are commonly applied in order to improve the corrosion resistance under harsh environments. In this work, Al 7075-T6 substrates were coated with nanostructured SiO2-TiO2 films using a sol-gel method. The experimental approach initially consisted in the preparation of a precursor agent using tetraethyl orthosilicate (TEOS) and triethoxy(octyl)silane (ETOS). Subsequently, nanoparticles of SiO2-TiO2 were mixed in order to develop thin films using a one-step dip coating method. The roughness, nanoindentation and corrosion properties were evaluated for the coated substrates. A finite element model was created for the nanoindentation test, which determined the mechanical response between the film-contact interface during loading conditions. The average hardness, elastic modulus and critical loads leading to fracture were verified. The nanoindentation test presented a significant increase in hardness for the coated Al 7075-T6 alloy, reaching a value of 4.6 GPa. The SiO2-TiO2 thin films presented uniform and compact surface coatings with high mechanical properties. Furthermore, the performed corrosion tests indicated moderate protection by the SiO2-TiO2 thin films. The SiO2-TiO2 thin films displayed a generalized corrosion throughout the surface, presenting oxides and fractured crystals in localized regions.



    Interactions between the predators and the preys are diverse and complex in ecology. The predators increase the mortality rate of preys by direct predation. In existing literatures, many predator-prey models only involve direct predation for predator-prey interactions [1,2,3,4]. However, besides the population loss caused by direct predation, the prey will modify their behavior, psychology and physiology in response to the predation risk. This is defined as the fear effect by Cannon [5]. Zanette et al. [6] investigated the variation of song sparrows offspring reproduction when the sounds and calls of predators were broadcasted to simulate predation risk. They discovered that the fear effect alone can cause a significant reduction of song sparrows offspring reproduction. Inspired by this experimental result, Wang et al. [7] incorporated the fear effect into the predator-prey model. In their theoretical analysis, linear and Holling type Ⅱ functional response are chosen respectively. According to their results, the structure of the equilibria will not be affected by the fear effects, but the stability of equilibria and Hopf bifurcation are slightly different from models with no fear effects.

    Soon afterward, Wang and Zou [8] considered a model with the stage structure of prey (juvenile prey and adult prey) and a maturation time delay. Additionally, Wang and Zou [9] pointed out anti-predation behaviors will not only decrease offspring reproduction of prey but also increase the difficulty of the prey being caught. Based on this assumption, they derived an anti-predation strategic predator-prey model.

    Recently, the aforementioned ordinary or delayed differential equation models were extended to reaction diffusion equation [10,11,12]. Wang et al. [11] used several functional responses to study the effect of the degree of prey sensitivity to predation risk on pattern formation. Following this work, Wang et al. [10] introduced spatial memory delay and pregnancy delay into the model. Their numerical simulation presented the effect of some biological important variables, including the level of fear effect, memory-related diffusion, time delay induced by spatial memory and pregnancy on pattern formation. Moreover, Dai and Sun [13] incorporated chemotaxis and fear effect into predator-prey model, and investigated the Turing-Hopf bifurcation by selecting delay and chemotaxis coefficient as two analysis parameters.

    Denote by u1(x,t) and u2(x,t), the population of the prey and adult predator at location x and time t, respectively. We suppose juvenile predators are unable to prey on. By choosing simplest linear functional response in the model of Wang and Zou [7], the equation for prey population is given by

    tu1d1Δu1=u1(rg(K,u2)μ1mu1)pu1u2,

    where d1 is a diffusion coefficient for prey, μ1 is the mortality rate for prey, m is the intraspecies competition coefficient, p is predation rate, r is reproduction rate for prey and g(K,u2) represents the cost of anti-predator defense induced by fear with K reflecting response level. We assume g(K,u2) satisfies the following conditions.

    (H1) g(0,u2)=g(K,0)=1, limu2g(K,u2)=limKg(K,u2)=0, g(K,u2)K<0 and g(K,u2)u2<0.

    Considering the maturation period of the predator, we set b(x,t,a1) be the density of the predator at age a1, location x and time t. Establish the following population model with spatial diffusion and age-structure

    (a1+t)b(x,t,a1)=d(a1)Δb(x,t,a1)μ(a1)b(x,t,a1),b(x,t,0)=cpu1(x,t)u2(x,t), (1.1)

    where xΩ, a bounded spatial habitat with the smooth boundary Ω, t,a1>0, c is the conversion rate of the prey to predators, τ>0 be the maturation period of predator and age-specific functions

    d(a1)={d0, a1τ,d2, a1>τ,  μ(a1)={γ, a1τ,μ2, a1>τ,

    represent the diffusion rate and mortality rate at age a1, respectively. We introduce the total population of the matured predator as u2=τb(x,t,a1)da1. Thus, (1.1) together with b(x,t,)=0 yields

    tu2=d2Δu2+b(x,t,τ)μ2u2. (1.2)

    Let s1=ta1 and w(x,t,s1)=b(x,t,ts1). Along the characteristic line, solving (1.1) yields

    tw(x,t,s1)={d0Δw(x,t,s1)γw(x,t,s1),  xΩ, 0ts1τ,d2Δw(x,t,s1)μ2w(x,t,s1), xΩ, ts1>τ,w(x,s1,s1)=b(x,s1,0)=cpu1(x,s1)u2(x,s1), s10;  w(x,0,s1)=b(x,0,s1), s1<0.

    Assume linear operator d0Δγ with Neumann boundary conditions yields the C0 semigroups T1(t). Therefore,

    b(x,t,τ)=w(x,t,tτ)={T1(τ)b(,tτ,0), t>τ,T1(t)b(,0,τt), tτ.

    In particular, G(x,y,t) denotes the kernel function corresponding to T1(t). Thus

    b(x,t,τ)=T1(τ)b(,tτ,0)=cpΩG(x,y,τ)u1(y,tτ)u2(y,tτ)dy, t>τ.

    The above equation together with (1.2) yields a nonlocal diffusive predator-prey model with fear effect and maturation period of predators

    u1t=d1Δu1+u1(x,t)[rg(K,u2(x,t))μ1mu1(x,t)]pu1(x,t)u2(x,t),u2t=d2Δu2+cpΩG(x,y,τ)u1(y,tτ)u2(y,tτ)dyμ2u2(x,t). (1.3)

    Since the spatial movement of mature predators is much bigger than that of juvenile predators, we assume that diffusion rate of juvenile predators d0 approaches zero. Hence, the kernel function becomes G=eγτf(xy) with a Dirac-delta function f. Thus, the equation of u2(x,t) in (1.3) becomes

    u2t=d2Δu2+cpeγτΩf(xy)u1(y,tτ)u2(y,tτ)dyμ2u2(x,t).

    It follows from the properties of Dirac-delta function that

    Ωf(xy)u1(y,tτ)u2(y,tτ)dy=limϵ0Bϵ(x)f(xy)u1(y,tτ)u2(y,tτ)dy=u1(x,tτ)u2(x,tτ)limϵ0Bϵ(x)f(xy)dy=u1(x,tτ)u2(x,tτ)

    where Bϵ(x) is the open ball of radius ϵ centered at x. Therefore, model (1.3) equipped with nonnegative initial conditions and Neumann boundary conditions is

    u1(x,t)t=d1Δu1(x,t)+u1(x,t)[rg(K,u2(x,t))μ1mu1(x,t)]pu1(x,t)u2(x,t), xΩ,t>0,u2(x,t)t=d2Δu2(x,t)+cpeγτu1(x,tτ)u2(x,tτ)μ2u2(x,t), xΩ,t>0,u1(x,t)ν=u2(x,t)ν=0,xΩ,t>0,u1(x,ϑ)=u10(x,ϑ)0,u2(x,ϑ)=u20(x,ϑ)0,xΩ, ϑ[τ,0]. (1.4)

    If rμ1, then as t, we have (u1,u2)(0,0) for xΩ, namely, two species will extinct. Throughout this paper, suppose r>μ1, which ensures that the prey and predator will persist.

    This paper is organized as follows. We present results on well-posedness and uniform persistence of solutions and prove the global asymptotic stability of predator free equilibrium in Section 2. The nonexistence of nonhomogeneous steady state and steady state bifurcation are proven in Section 3. Hopf bifurcation analysis is carried out in Section 4. In Section 5, we conduct numerical exploration to illustrate some theoretical conclusions and further explore the dynamics of the nonlocal model numerically. We sum up our paper in Section 6.

    Denote by C:=C([τ,0],X2) the Banach space of continuous maps from [τ,0] to X2 equipped with supremum norm, where X=L2(Ω) is the Hilbert space of integrable function with the usual inner product. C+ is the nonnegative cone of C. Let u1(x,t) and u2(x,t) be a pair of continuous function on Ω×[τ,) and (u1t,u2t)C as (u1t(ϑ),u2t(ϑ))=(u1(,t+ϑ),u2(,t+ϑ)) for ϑ[τ,0]. By using [14,Corollary 4], we can prove model (1.4) exists a unique solution. Note (1.4) is mixed quasi-monotone [15], together with comparison principle implies that the solution of (1.4) is nonnegative.

    Lemma 2.1. For any initial condition (u10(x,ϑ),u20(x,ϑ))C+, model (1.4) possesses a unique solution (u1(x,t),u2(x,t)) on the maximal interval of existence [0,tmax). If tmax<, then lim supttmax(u1(,t)+u2(,t))=. Moreover, u1 and u2 are nonnegative for all (x,t)¯Ω×[τ,tmax).

    We next prove u1 and u2 are bounded which implies that tmax=.

    Theorem 2.2. For any initial condition φ=(u10(x,ϑ),u20(x,ϑ))C+, model (1.4) possesses a global solution (u1(x,t),u2(x,t)) which is unique and nonnegative for (x,t)¯Ω×[0,). If u10(x,ϑ)0(, then this solution remains positive for all (x, t)\in \overline\Omega\times(0, \infty) . Moreover, there exists a positive constant \xi independent of \varphi such that \limsup\limits_{t\rightarrow \infty}u_1\le \xi, \ \limsup\limits_{t\rightarrow \infty}u_2\le \xi for all x\in\Omega.

    Proof. Consider

    \begin{equation} \begin{aligned} {\partial w_1\over \partial t}& = d_1 \Delta w_1+rw_1-\mu_1 w_1-m w_1^2, \ x\in\Omega, t > 0,\\ {\partial w_1\over \partial \nu} & = 0, \ x\in\partial\Omega, t > 0, \ \ w_1(x,0) = \sup\limits_{\vartheta\in[-\tau,0]}u_{10}(x,\vartheta), \ x\in\Omega. \end{aligned} \end{equation} (2.1)

    Clearly, w_1(x, t) of (2.1) is a upper solution to (1.4) due to {\partial u_1/\partial t}\le d_1 \Delta u_1+ru_1-\mu_1 u_1-m u_1^2 . By using Lemma 2.2 in [16], (r-\mu_1)/m of (2.1) is globally asymptotically stable in C(\overline\Omega, \mathbb{R}^+) . This together with comparison theorem indicates

    \begin{equation} \limsup\limits_{t\rightarrow \infty}u_1\le \lim\limits_{t\rightarrow \infty} w_1 = {r-\mu_1\over m} \ \ \text{uniformly for} \ \ x\in\overline\Omega. \end{equation} (2.2)

    Thus, there exists \widetilde{\xi} > 0 which is not dependent on initial condition, such that \|u_1\|\le \widetilde{\xi} for all t > 0 . T_2(t) denotes the C_0 semigroups yielded by d_2\Delta-\mu_2 with the Neumann boundary condition. Then from (1.4),

    \begin{aligned} u_2& = T_2(t)u_{20}(\cdot,0)+cpe^{-\gamma\tau}\int_{-\tau}^{t-\tau}T_2(t-\tau-a)u_1(\cdot,a)u_2(\cdot,a)da. \end{aligned}

    Let -\delta < 0 be the principle eigenvalue of d_2\Delta-\mu_2 with the Neumann boundary condition. Then, \|T_2(t)\|\le e^{-\delta t} . The above formula yields

    \begin{aligned} \|u_2(\cdot,t)\|&\le e^{-\delta t}\|u_{20}(\cdot,0)\|+cpe^{-\gamma\tau}\widetilde{\xi}\int_{-\tau}^{t-\tau}e^{-\delta(t-\tau-a)}\|u_2(\cdot,a)\|da,\\ &\le \widetilde{B}_1+\int_{0}^{t}\widetilde{B}_2\|u_2(\cdot,a)\|da, \end{aligned}

    by choosing constants \widetilde{B}_1\ge\|u_2(\cdot, 0)\|+ cp\tau \widetilde{\xi}\sup\limits_{a\in[-\tau, 0]}\|u_2(\cdot, a)\| and \widetilde{B}_2\ge cp\widetilde{\xi} . Using Gronwall's inequality yields \|u_2(\cdot, t)\|\le \widetilde{B}_1 e^{\widetilde{B}_2t} for all 0\le t < t_{max} . Lemma 2.1 implies t_{max} = \infty . So (u_1, u_2) is a global solution. Moreover, if u_{10}(x, \vartheta)\ge 0(\not\equiv0), u_{20}(x, \vartheta)\ge 0(\not\equiv0) , then by [17,Theorem 4], this solution is positive for all t > 0 and x\in\overline\Omega .

    We next prove (u_1, u_2) is ultimately bounded by a constant which is not dependent on the initial condition. Due to (2.2), there exist t_0 > 0 and \xi_0 > 0 such that u_1(x, t)\le \xi_0 for any t > t_0 and x\in\overline\Omega . Let z(x, t) = c u_1(x, t-\tau)+u_2(x, t) , \mu = \min\{\mu_1, \mu_2\} and I_1 = \int_\Omega z(x, t)dx . We integrate both sides of (1.4) and add up to obtain

    I'_1(t)\le \int_\Omega \left(c(r-\mu_1)u_1(x,t-\tau)-\mu_2 u_2(x,t)\right)dx\le cr\xi_0|\Omega|-\mu I_1(t),\ t\ge t_0+\tau.

    Comparison principle implies

    \limsup\limits_{t\to\infty}\|u_2\|_1\le\limsup\limits_{t\to\infty}I_1\le cr\xi_0|\Omega|/\mu.

    Especially, there exist t_1 > t_0 and \xi_1 > 0 such that \|u_2\|_{1}\le \xi_1 for all t\ge t_1 .

    Now, we define V_{l}(t) = \int_{\Omega}(u_2(x, t))^ldx with l\ge1 , and estimate the upper bound of V_{2}(t) . For t > t_1 , the second equation of (1.4) and Young's inequality yield

    \begin{equation*} \label{V2} \begin{aligned} {1\over2}V'_2(t)&\le -d_2\int_{\Omega} |\nabla u_2|^2dx+cp\xi_0\int_{\Omega}u_2(x,t-\tau)u_2(x,t)dx-\mu_2 V_2(t)\\ &\le -d_2\|\nabla u_2\|_2 ^2+{cp\xi_0\over 2}V_2(t-\tau)+{cp\xi_0\over 2}V_2(t). \end{aligned} \end{equation*}

    The Gagliardo-Nirenberg inequality states:

    \begin{equation*} \label{GN} \forall\ \epsilon > 0, \exists\ \hat c > 0,\ \text{s.t.}\ \|\mathcal{P}\|_2^2\le \epsilon\|\nabla \mathcal{P}\|_2^2+\hat c\epsilon^{-n/2}\|\mathcal{P}\|_1^2,\ \forall\ \mathcal{P}\in W^{1,2}(\Omega){.} \end{equation*}

    We obtain

    V'_2(t)\le C_1+C_2 V_2(t-\tau)-(C_2+C_3)V_2(t),

    where C_1 = \hat c\epsilon^{-n/2-1}2d_2\xi_1^2 > 0 , C_2 = cp\xi_0 > 0 and C_3 = 2(d_2/\epsilon-C_2) > 0 with small \epsilon\in(0, d_2/C_2) . Using comparison principle again yields \limsup\limits_{t\to\infty}V_2(t)\le C_1/C_3 , which implies there exist t_2 > t_1 and \xi_2 > 0 such that V_2\le \xi_2 for t\ge t_2 .

    Let L_l = \limsup\limits_{t\to\infty} V_l(t) with l\ge1 , we want to estimate L_{2l} with the similar method of estimation for L_2 . Multiply the second equation in (1.4) by 2lv^{2l-1} and integrate on \Omega . Young's inequality implies

    V'_{2l}(t)\le -2d_2\int_{\Omega}|\nabla u_2^l|^2dx+cp\xi_0 V_{2l}(t-\tau)+(2l-1)cp\xi_0 V_{2l}(t).

    Then

    V'_{2l}(t)\le 2 d_2 \hat c\epsilon^{-n/2-1} V^2_{l}(t)-{2 d_2\over \epsilon}V_{2l}(t)+cp\xi_0 (V_{2l}(t-\tau)+(2l-1)V_{2l}(t)),

    via Gagliardo-Nirenberg inequality. Since L_l = \limsup\limits_{t\to\infty}V_{l}(t) , there exists t_l > 0 such that V_{l}\le 1+L_l when t > t_l . Hence,

    V'_{2l}(t)\le 2 d_2 \hat c\epsilon^{-(n/2+1)} (1+L_l)^2-{2 d_2\over \epsilon}V_{2l}(t)+lC_4 (V_{2l}(t)+V_{2l}(t-\tau))

    with C_4 = 2cp\xi_0 . We choose \epsilon^{-1} = (2C_4+1)l/(2d_2) and C_5 = 2d_2\hat c[(2 C_4+1)/(2d_2)]^{n/2+1} . Then for t > t_l , we obtain

    V'_{2l}(t)\le C_5 l^{n/2+1}(1+L_l)^2+ l C_4 V_{2l}(t-\tau)-(l C_4+l)V_{2l}(t).

    By comparison principle, the above inequality yields L_{2l}\le C_5 l^{n/2}(1+L_l)^2 , with a constant C_5 which is not dependent on l and initial conditions. Finally, prove L_{2^s} < \infty for all s\in\mathbb{N}_0 . Let B = 1+C_5 and \{b_s\}_{s = 0}^\infty be an infinite sequence denoted by b_{s+1} = B^{(1/2)^{(s+1)}}2^{sn((1/2)^{(s+2)})}b_s with the first term b_0 = L_1+1 . Clearly, L_{2^s}\le(b_s)^{2^s} and

    \lim\limits_{s\to\infty}\ln b_s = \ln b_0+\ln B+\frac{n}{2}\ln 2.

    Therefore,

    \limsup\limits_{s\to\infty}(L_{2^s})^{(1/2)^{s}}\le \lim\limits_{s\to\infty}b_s = B(1+L_1)2^{n/2}\le B(1+\xi_1)2^{n/2}\le \xi: = \max\{B(\xi_1+1)2^{n/2}, (r-\mu_1)/m\}.

    The above inequality leads to \limsup\limits_{t\rightarrow \infty}u_1\le \xi and \limsup\limits_{t\rightarrow \infty}u_2\le \xi for all x\in\Omega .

    In Theorem 2.2, we proved that the solution of model (1.4) is uniformly bounded for any nonnegative initial condition, this implies the boundedness of the population of two species. Clearly model (1.4) exists two constant steady states E_0 = (0, 0) and E_1 = ((r-\mu_1)/m, 0) , where E_0 is a saddle. Define the basic reproduction ratio [18] by

    \begin{equation} {R}_0 = {cpe^{-\gamma\tau}(r-\mu_1)\over m\mu_2}. \end{equation} (2.3)

    Thus model (1.4) possesses exactly one positive constant steady state E_2 = (u_1^*, u_2^*) if and only if {R}_0 > 1 , which is equivalent to cp(r-\mu_1) > m\mu_2 and 0\le\tau < \tau_{max}: = {1\over \gamma}\ln {cp(r-\mu_1)\over m\mu_2} . Here,

    u_1^* = {\mu_2 e^{\gamma\tau}\over pc}, \ \ u_2^* \ \ \text{satisfies} \ \ rg(K,u_2)-p u_2 = \mu_1+mu_1^*.

    The linearization of (1.4) at the positive constant steady state (\tilde u_1, \tilde u_2) gives

    \begin{equation} {\partial W/\partial t} = \mathcal{D}\Delta W+\mathcal{L}(W_t), \end{equation} (2.4)

    where domain Y: = \{(u_1, u_2)^T:u_1, u_2\in C^2(\Omega)\cap C^1(\bar\Omega), (u_1)_\nu = (u_2)_\nu = 0\ \text{on}\ \partial\Omega\} , W = (u_1(x, t), u_2(x, t))^T , \mathcal{D} = \text{diag}(d_1, d_2) and a bounded linear operator \mathcal{L}:\mathbf{C}\to X^2 is

    \mathcal{L}(\varphi) = M\varphi(0)+M_\tau\varphi(-\tau),\ \text{for}\ \varphi\in\mathbf{C},

    with

    \begin{equation*} M = \begin{pmatrix} rg(K,\tilde u_2)-\mu_1-2m\tilde u_1-p\tilde u_2 & \tilde u_1(rg'_{u_2}(K,\tilde u_2)-p)\\ 0 & -\mu_2 \end{pmatrix}, \ M_\tau = \begin{pmatrix} 0 & 0\\ cpe^{-\gamma\tau}\tilde u_2 & cpe^{-\gamma\tau}\tilde u_1 \end{pmatrix}. \end{equation*}

    The characteristic equation of (2.4) gives

    \rho \eta-\mathcal{D}\Delta \eta-\mathcal{L}(e^{\rho\cdot}\eta) = 0,\ \text{for some}\ \eta\in Y\backslash\{0\},

    or equivalently

    \begin{equation} \det(\rho I+\sigma_n \mathcal{D}-M-e^{-\rho\tau}M_\tau) = 0,\ \text{for}\ n\in\mathbb{N}_0. \end{equation} (2.5)

    Here, \sigma_n is the eigenvalue of -\Delta in \Omega with Neumann boundary condition with respect to eigenfunction \psi_n for all n\in\mathbb{N}_0 , and

    \begin{equation} 0 = \sigma_0 < \sigma_1\le\sigma_2\le\cdots\le\sigma_n\le\sigma_{n+1}\le\cdots \ \ \text{and} \ \ \lim\limits_{n\to\infty}\sigma_n = \infty. \end{equation} (2.6)

    Theorem 2.3. (i) The trivial constant steady state E_0 = (0, 0) is always unstable.

    (ii) If R_0 > 1 , then E_1 = ((r-\mu_1)/m, 0) is unstable, and model (1.4) possesses a unique positive constant steady state E_2 = (u_1^*, u_2^*) .

    (iii) If R_0\le1 , then E_1 is globally asymptotically stable in \mathbf{C}_+ .

    Proof. (ⅰ) Note, (2.5) at E_0 takes form as (\rho+\sigma_n d_1-r+\mu_1)(\rho+\sigma_n d_2+\mu_2) = 0 for all n\in\mathbb{N}_0 . Then r-\mu_1 > 0 is a positive real eigenvalue, namely, E_0 is always unstable.

    (ⅱ) The characteristic equation at E_1 gives

    \begin{equation} (\rho+\sigma_n d_1+r-\mu_1) (\rho+\sigma_n d_2+\mu_2-cpe^{-\gamma\tau}{r-\mu_1 \over m}e^{-\rho\tau}) = 0 \ \ \text{for}\ n\in\mathbb{N}_0. \end{equation} (2.7)

    Note that one eigenvalue \rho_1 = -\sigma_n d_1-r+\mu_1 remains negative. Hence, we only need to consider the root distribution of the following equation

    \begin{equation} \rho+\sigma_n d_2+\mu_2-cpe^{-\gamma\tau}{r-\mu_1 \over m}e^{-\rho\tau} = 0\ \ \text{for}\ n\in\mathbb{N}_0. \end{equation} (2.8)

    According to [19,Lemma 2.1], we obtain that (2.8) exists an eigenvalue \rho > 0 when R_0 > 1 , namely, E_1 is unstable when R_0 > 1 .

    (ⅲ) By using [19,Lemma 2.1] again, any eigenvalue \rho of (2.8) satisfies Re(\rho) < 0 when R_0 < 1 , namely, E_1 is locally asymptotically stable when R_0 < 1 . Now, consider R_0 = 1 . 0 is an eigenvalue of (2.7) for n = 0 and all other eigenvalues satisfy Re(\rho) < 0 . To prove the stability of E_1 , we shall calculate the normal forms of (1.4) by the algorithm introduced in [20]. Set

    \Upsilon = \{\rho\in\mathbb{C},\rho\ \text{is the eigenvalue of equation (2.7) and}\ Re\rho = 0\}.

    Obviously, \Upsilon = \{0\} when R_0 = 1 . System (1.4) satisfies the non-resonance condition relative to \Upsilon . Denote \overline{u}_1 = (r-\mu_1)/m , and let \mathbf{w} = ({w}_1, {w}_2)^T = (\overline{u}_1-u_1, u_2)^T and (1.4) can be written as

    \dot{\mathbf{w}}_t = \mathcal{A}_0 \mathbf{w}_t+\mathcal{F}_0(\mathbf{w}_t) \ \ \text{on} \ \ \mathbf{C}.

    Here linear operator \mathcal{A}_0 is given by (\mathcal{A}_0 \varphi)(\vartheta) = (\varphi(\vartheta))' when \vartheta\in[-\tau, 0) and

    \begin{align*} (\mathcal{A}_0\varphi)(0) = \begin{pmatrix} d_1\Delta&0\\0&d_2\Delta \end{pmatrix}\varphi(0) +\begin{pmatrix} -r+\mu_1 & (p-rg'_{u_2}(K,0))\overline{u}_1\\0&-\mu_2 \end{pmatrix}\varphi(0) +\begin{pmatrix} 0&0\\0&\mu_2 \end{pmatrix}\varphi(-\tau), \end{align*}

    and the nonlinear operator \mathcal{F}_0 satisfies [\mathcal{F}_0(\varphi)](\vartheta) = 0 for -\tau\le\vartheta < 0 . By Taylor expansion, [\mathcal{F}_0(\varphi)](0) can be written as

    \begin{align} [\mathcal{F}_0(\varphi)](0) = \begin{pmatrix} m\varphi_1^2(0)-r\overline{u}_1g''_{u_2}(K,0)\varphi_2^2(0)/2+(rg'_{u_2}(K,0)-p)\varphi_1(0)\varphi_2(0)\\ -cpe^{-\gamma\tau}\varphi_1(-\tau)\varphi_2(-\tau) \end{pmatrix} +h.o.t. \end{align} (2.9)

    Define a bilinear form

    \begin{equation*} \langle\beta,\alpha\rangle = \int_{\Omega}\left[\alpha_1(0)\beta_1(0)+\alpha_2(0)\beta_2(0)+\mu_2\int_{-\tau}^0\beta_2(\vartheta+\tau) \alpha_2(\vartheta)d \vartheta \right]dx,\ \beta\in C([0,\tau],X^2),\ \alpha\in\mathbf{C}. \end{equation*}

    Select \alpha = (1, m/(p-r g'_{u_2}(K, 0)) and \beta = (0, 1)^T to be the right and left eigenfunction of \mathcal{A}_0 relative to eigenvalue 0 , respectively. Decompose \mathbf{w}_t as \mathbf{w}_t = h\alpha+\delta and \langle\beta, \delta\rangle = 0 . Notice \mathcal{A}_0\alpha = 0 and \langle\beta, \mathcal{A}_0\delta\rangle = 0 . Thus,

    \langle \beta, \dot{\mathbf{w}}_t\rangle = \langle \beta, \mathcal{A}_0 \mathbf{w}_t\rangle+\langle \beta, \mathcal{F}_0(\mathbf{w}_t)\rangle = \langle \beta, \mathcal{F}_0(\mathbf{w}_t)\rangle.

    Moreover,

    \langle \beta, \dot{\mathbf{w}}_t\rangle = \dot{h}\langle \beta, \alpha\rangle+\langle \beta, \dot{\delta}\rangle = \dot{h}\langle \beta, \alpha\rangle.

    It follows from the above two equations that

    \dot{h}\frac{m(1+\mu_2\tau)|\Omega|}{p-rg'_{u_2}(K,0)} = \langle\beta,\mathcal{F}_0(h\alpha+\delta)\rangle = \int_{\Omega}\beta^T[\mathcal{F}_0(h\alpha+\delta)](0)dx = \int_{\Omega}[\mathcal{F}_0(h\alpha+\delta)]_2(0)dx.

    When the initial value is a small perturbation of E_1 , then \delta = O(h^2) , together with Taylor expansion yields

    \begin{equation*} [\mathcal{F}_0(h\alpha+\delta)]_{2}(0) = -cpe^{-\gamma\tau}(h\alpha_1(-\tau)+\delta_1(-\tau))(h\alpha_2(-\tau)+\delta_2(-\tau)) = -\frac{cpe^{-\gamma\tau}m}{p-rg'_{u_2}(K,0)}h^2+O(h^3). \end{equation*}

    Therefore, we obtain the norm form of (1.4) as follows

    \begin{equation} \dot{h} = \frac{-cpe^{-\gamma\tau}}{1+\mu_2\tau}h^2+O(h^3). \end{equation} (2.10)

    Then for any positive initial value, the stability of zero solution of (2.10) implies E_1 is locally asymptotically stable when R_0 = 1 .

    Next, it suffices to show the global attractivity of E_1 in \mathbf{C}_+ when R_0\le1 . Establish a Lyapunov functional V:\mathbf{C}_+\to\mathbb{R} as

    V(\phi_1,\phi_2) = \int_{\Omega}\phi_2(0)^2 dx+cpe^{-\gamma\tau}{r-\mu_1 \over m}\int_{\Omega}\int_{-\tau}^0 \phi_2(\theta)^2d\theta dx \ \ \text{for} \ \ (\phi_1,\phi_2)\in\mathbf{C}_+.

    Along solutions of (1.4), taking derivative of V(\phi_1, \phi_2) with respect to t yields

    \begin{equation*} \begin{aligned} {d V \over dt}&\le -2d_2\int_{\Omega}|\nabla u_2|^2dx+ \int_{\Omega}2cpe^{-\gamma\tau}\frac{r-\mu_1}{m}u_2(x,t-\tau)u_2(x,t)\\ &-2\mu_2 u_2^2(x,t)+cpe^{-\gamma\tau}\frac{r-\mu_1}{m}[u_2^2(x,t)-u_2^2(x,t-\tau)]dx\\ &\le\int_{\Omega}2\mu_2(R_0-1)u_2^2(x,t)dx\le0 \ \ \text{if} \ \ R_0\le1. \end{aligned} \end{equation*}

    Note \{E_1\} is the maximal invariant subset of dV/dt = 0 , together with LaSalle-Lyapunov invariance principle [21,22] implies E_1 is globally asymptotically stable if R_0\le1 .

    In Theorem 2.3, E_0 is always unstable which suggests that at leat one species will persist eventually. Moreover, if R_0\le1 , then E_1 is globally asymptotically stable in \mathbf{C}_+ , which implies that when the basic reproduction ratio is no more than one, the predator species will extinct and only the prey species can persist eventually. Next, we will prove the solution is uniformly persistent. \Theta_t denotes the solution semiflow of (1.4) mapping \mathbf{C}_+ to \mathbf{C}_+ ; namely, \Theta_t\varphi: = (u_1(\cdot, t+\cdot), u_2(\cdot, t+\cdot))\in \mathbf{C}_+ . Set \zeta^+(\varphi) = \cup_{t\ge0}\{\Theta_t\varphi\} be the positive orbit and \varpi(\varphi) be the omega limit set of \zeta^+(\varphi) . Denote

    Z: = \{(\varphi_1,\varphi_2)\in\mathbf{C}_+:\varphi_1\not\equiv0\ \text{and}\ \varphi_2\not\equiv0\},\ \partial Z: = \mathbf{C}_+\backslash Z = \{(\varphi_1,\varphi_2)\in\mathbf{C}_+:\varphi_1\equiv0\ \text{or}\ \varphi_2\equiv0\},

    \Gamma_{\partial} as the largest positively invariant set in \partial Z , and W^s((\tilde u_1, \tilde u_2)) as the stable manifold associated with (\tilde u_1, \tilde u_2) . We next present persistence result of model (1.4).

    Theorem 2.4. Suppose R_0 > 1 . Then there exists \kappa > 0 such that \liminf\limits_{t\to\infty} u_1(x, t)\ge\kappa and \liminf\limits_{t\to\infty} u_2(x, t)\ge\kappa for any initial condition \varphi\in Z and x\in\overline\Omega .

    Proof. Note \Theta_t is compact, and Theorem 2.2 implies \Theta_t is point dissipative. Then \Theta_t possesses a nonempty global attractor in \mathbf{C}_+ [23]. Clearly, \Gamma_{\partial} = \{(\varphi_1, \varphi_2)\in\mathbf{C}_+: \varphi_2\equiv0\} , and \varpi(\varphi) = \{E_0, E_1\} for all \varphi\in \Gamma_{\partial} . Define a generalized distance function \psi mapping \mathbf{C}_+ to \mathbb{R}^+ by

    \psi(\varphi) = \min\limits_{x\in\bar\Omega}\{\varphi_1(x,0),\varphi_2(x,0)\},\ \forall \varphi = (\varphi_1,\varphi_2)\in\mathbf{C}_+.

    Following from strong maximum principle [24], \psi(\Theta_t\varphi) > 0 for all \varphi\in Z . Due to \psi^{-1}(0, \infty)\subset Z , assumption (P) in [25,Section 3] holds. Then verify rest conditions in [25,Theorem 3].

    First, prove W^s(E_0)\cap \psi^{-1}(0, \infty) = \emptyset . Otherwise, there exists an initial condition \varphi\in \mathbf{C}_+ with \psi(\varphi) > 0 , such that (u_1, u_2)\to E_0 as t\to\infty . Thus, for any sufficiently small \varepsilon_1 > 0 satisfying r g(K, \varepsilon_1)-\mu_1 > p\varepsilon_1 , there exists t_1 > 0 such that 0 < u_1, u_2 < \varepsilon_1 for all x\in \Omega and t > t_1 . Note that rg(K, 0)-\mu_1 > 0 and \partial g(K, u_2)/\partial u_2 < 0 ensure the existence of small \varepsilon_1 > 0 . Then the first equation in (1.4) and (\mathbf{H}_1) lead to

    \partial_t u_1 > d_1\Delta u_1+u_1 \left[(r g(K,\epsilon_1)-\mu_1-p\epsilon_1)-mu_1\right],t > t_1.

    Notice

    \begin{equation*} \begin{aligned} &\partial_t \hat u_1-d_1\Delta \hat u_1 = \hat u_1\left[(rg(K,\epsilon_1)-p\epsilon_1-\mu_1)-m\hat u_1\right], \ x\in \Omega,t > t_1,\\ &\partial_\nu \hat u_1 = 0, \ x\in\partial\Omega,t > t_1, \end{aligned} \end{equation*}

    has a globally asymptotically stable positive steady state (r g(K, \epsilon_1)-\mu_1-p\epsilon_1)/m due to [16,Lemma 2.2], together with comparison principle yields \lim\limits_{t\to\infty}u_1\ge \lim\limits_{t\to\infty}\hat u_1 > 0 . A contradiction is derived, so W^s(E_0)\cap \psi^{-1}(0, \infty) = \emptyset .

    Next check W^s(E_1)\cap\psi^{-1}(0, \infty) = \emptyset . If not, there exists \varphi\in\mathbf{C}_+ with \psi(\varphi) > 0 such that (u_1, u_2) converges to E_1 as t\to\infty . According to (2.2), for any small \varepsilon_2 > 0 satisfying cpe^{-\gamma\tau}((r-\mu_1)/m-\epsilon_2) > \mu_2 , there exists t_2 > 0 such that u_1 > (r-\mu_1)/m-\varepsilon_2 for all x\in\overline\Omega and t > t_2-\tau . Note that R_0 > 1 ensures the existence of \varepsilon_2 > 0 . Thus, the second equation of (1.4) yields

    \partial_t u_2 > d_2\Delta u_2+cpe^{-\gamma\tau}({r-\mu_1 \over m}-\varepsilon_2)u_2(x,t-\tau)-\mu_2u_2(x,t), \ t > t_2.

    In a similar manner, we derive \lim\limits_{t\to\infty}u_2(x, t) > 0 by above inequality, cpe^{-\gamma\tau}((r-\mu_1)/m-\epsilon_2) > \mu_2 and comparison principle. A contradiction yields again. Hence, it follows from Theorem 3 in [25] that, for any \varphi\in \mathbf{C}_+ , there exists \kappa > 0 such that \liminf\limits_{t\to\infty}\psi(\Theta_t \varphi)\ge\kappa uniformly for any x\in\overline\Omega .

    Theorem 2.4 implies that when the basic reproduction ratio is bigger than one, both the predator species and prey species will persist eventually. We next investigate the stability of E_2 . The corresponding characteristic equation at E_2 gives

    \begin{equation} \rho^2+a_{1,n}\rho+a_{0,n}+(b_{1,n}\rho+b_{0,n})e^{-\rho\tau} = 0, \ n\in\mathbb{N}_0, \end{equation} (2.11)

    with

    \begin{align*} &a_{1,n} = \sigma_n(d_1+d_2)+\mu_2+mu_1^* > 0,\ a_{0,n} = (\sigma_nd_1+au_1^*)(\sigma_nd_2+\mu_2) > 0,\\ &b_{1,n} = -\mu_2 < 0, \ \ b_{0,n} = -\mu_2(\sigma_nd_1+mu_1^*+u_2^*(rg'_{u_2}(K,u_2^*)-p)). \end{align*}

    Characteristic equation (2.11) with \tau = 0 is

    \begin{equation} \rho^2+(a_{1,n}+b_{1,n})\rho+a_{0,n}+b_{0,n} = 0, \ n\in\mathbb{N}_0. \end{equation} (2.12)

    We observe that a_{0, n}+b_{0, n} = \sigma_nd_2(\sigma_nd_1+mu_1^*)-\mu_2 u_2^*(rg'_{u_2}(K, u_2^*)-p) > 0 , and a_{1, n}+b_{1, n} = \sigma_n(d_1+d_2)+mu_1^* > 0 for all integer n\ge0 , which yields that any eigenvalue \rho of (2.12) satisfies Re(\rho) < 0 . Then, local asymptotic stability of E_2 is derived when \tau = 0 which implies Turing instability can not happen for the non-delay system of (1.4). In addition, a_{0, n}+b_{0, n} > 0 for any n\in\mathbb{N}_0 leads to that (2.11) can not have an eigenvalue 0 for any \tau\ge0 . This suggests we look for the existence of simple \rho = \pm i\delta ( \delta > 0 ) for some \tau > 0 . Substitute \rho = i\delta into (2.11) and then

    \begin{equation} G_n(\delta,\tau) = \delta^4+(a_{1,n}^2-2a_{0,n}-b_{1,n}^2)\delta^2+a_{0,n}^2-b_{0,n}^2 = 0, \ \ n\in\mathbb{N}_0, \end{equation} (2.13)

    with

    \begin{align*} &a_{1,n}^2-2a_{0,n}-b_{1,n}^2 = (\sigma_nd_1+mu_1^*)^2+(\sigma_nd_2)^2+2\mu_2\sigma_nd_2 > 0,\\ &a_{0,n}+b_{0,n} = \sigma_nd_2(\sigma_nd_1+mu_1^*)-\mu_2 u_2^*(rg'_{u_2}(K,u_2^*)-p) > 0,\\ &a_{0,n}-b_{0,n} = \sigma_n^2d_1d_2+\sigma_n(2\mu_2d_1+mu_1^*d_2)+\mu_2(2mu_1^*+u_2^*(rg'_{u_2}(K,u_2^*)-p)). \end{align*}

    Thus, a_{0, n}-b_{0, n}\ge 0 for all n\in\mathbb{N}_0 is equivalent to

    (\mathbf{A}_0): \ 2mu_1^*\ge u_2^*(p-rg'_{u_2}(K,u_2^*)).

    If (\mathbf{A}_0) holds, then (2.13) admits no positive roots, together with for \tau = 0 , any eigenvalues \rho of (2.11) satisfies Re(\rho) < 0 , yields the next conclusion.

    Theorem 2.5. Suppose R_0 > 1 . Then, E_2 is locally asymptotically stable provided that (\mathbf{A}_0) holds.

    Now, we consider positive nonhomogeneous steady states. The steady state (u_1(x), u_2(x)) of (1.4) satisfies the elliptic equation

    \begin{equation} \begin{aligned} -d_1\Delta u_1& = r g(K,u_2)u_1-\mu_1 u_1-mu_1^2-pu_1u_2, \ x\in\Omega,\\ -d_2\Delta u_2& = cpe^{-\gamma\tau}u_1u_2-\mu_2u_2, \ x\in\Omega,\\ \partial_{\nu} u_1& = \partial_{\nu} u_2 = 0, \ x\in\partial\Omega. \end{aligned} \end{equation} (3.1)

    From Theorem 2.3, all the solutions converge to E_1 when R_0\le1 and the positive nonhomogeneous steady state may exist only if R_0 > 1 . Throughout this section, we assume that R_0 > 1 . In what follows, the positive lower and upper bounds independent of steady states for all positive solutions to (3.1) are derived.

    Theorem 3.1. Assume that R_0 > 1 . Then any nonnegative steady state of (3.1) other than (0, 0) , and ((r-\mu_1)/m, 0) should be positive. Moreover, there exist constants \overline{\mathcal{B}}, \underline{\mathcal{B}} > 0 which depend on all parameters of (3.1) and \Omega , such that \underline{\mathcal{B}}\le u_1(x), u_2(x)\le\overline{\mathcal{B}} for any positive solution of (3.1) and x\in\overline\Omega .

    Proof. We first show any nonnegative solution (u_1, u_2) other than E_0 and E_1 , should be u_1 > 0 and u_2 > 0 for all x\in\overline\Omega . To see this, suppose u_2(x_0) = 0 for some x_0\in\overline\Omega , then u_2(x)\equiv0 via strong maximum principle and

    0\le d_1\int_{\Omega}|\nabla(u_1-{{r-\mu_1}\over m})|^2dx = \int_{\Omega}-mu_1(x)(u_1(x)-{{r-\mu_1}\over m})^2dx\le0.

    Thus the above inequality implies u_1(x)\equiv0 or u_1(x)\equiv (r-\mu_1)/m . Now, we assume u_2 > 0 for all x\in\overline \Omega . Strong maximum principle yields u_1 > 0 for all x\in\overline \Omega . Hence, u_1 > 0 and u_2 > 0 for all x\in\overline\Omega .

    We now prove u_1 and u_2 have a upper bound which is a positive constant. Since -d_1\Delta u_1(x)\le (r-\mu_1-mu_1(x))u_1(x) , we then obtain from Lemma 2.3 in [26] that u_1(x)\le (r-\mu_1)/m for any x\in\overline\Omega .

    By two equations in (3.1), we obtain

    -\Delta(d_1cu_1+d_2u_2)\le rc (r-\mu_1)/m -\min\{{\mu_1\over d_1}, {\mu_2\over d_2}\}(d_1cu_1+d_2u_2).

    By using [26,Lemma 2.3] again, we conclude

    u_1(x),u_2(x)\le \overline{\mathcal{B}} = {rc (r-\mu_1)\over{m\min\{c\mu_1,\mu_2,{\mu_1d_2/d_1}, {\mu_2d_1c/d_2}\}}}.

    Next, we only need to prove \|u_1(x)\| and \|u_2(x)\| have a positive lower bound which is not dependent on the solution. Otherwise, there exists a positive steady states sequence (u_{1, n}(x), u_{2, n}(x)) such that either \lim\limits_{n\to\infty}\|u_{1, n}\|_{\infty} = 0 or \lim\limits_{n\to\infty}\|u_{2, n}\|_{\infty} = 0 . Integrating second equation of (3.1) gives

    \begin{equation} 0 = \int_{\Omega}u_{2,n}(cpe^{-\gamma\tau}u_{1,n}-\mu_2)dx. \end{equation} (3.2)

    If \|u_{1, n}(x)\|_\infty\to0 as n\to\infty , then cpe^{-\gamma\tau}u_{1, n}(x)-\mu_2 < -\mu_2/2 for sufficiently large n , which yields u_{2, n}(cpe^{-\gamma\tau}u_{1, n}-\mu_2) < 0 , a contradiction derived. Thus, \|u_{2, n}(x)\|_{\infty}\to0 as n\to\infty holds. We then assume that (u_{1, n}, u_{2, n})\to(u_{1, \infty}, 0) , as n\to\infty where u_{1, \infty}\ge0 . Similarly, we obtain that either u_{1, \infty}\equiv0 or u_{1, \infty}\equiv (r-\mu_1)/m . Obviously, u_{1, \infty}\not\equiv0 based on the above argument, thus u_{1, \infty}\equiv (r-\mu_1)/m and \lim\limits_{n\to \infty}cpe^{-\gamma\tau}u_{1, n}(x)-\mu_2 = \mu_2(R_0-1) > 0 . This again contradicts (3.2). Hence, we have shown \|u_1(x)\|_\infty and \|u_2(x)\|_\infty have a positive constant lower bound independent on the solution. Therefore, u_1(x) and u_2(x) have a uniform positive constant lower bound independent on the solution of (3.1) via Harnack's inequality [26,Lemma 2.2]. This ends the proof.

    Theorem 3.2. Suppose R_0 > 1 . There exists a constant \chi > 0 depending on r, \mu_1, p, c, \gamma, \tau, \mu_2, g and \sigma_1 , such that if \min\{d_1, d_2\} > \chi then model (1.4) admits no positive spatially nonhomogeneous steady states, where \sigma_1 is defined in (2.6).

    Proof. Denote the averages of the positive solution (u_1, u_2) of system (3.1) on \Omega by

    \widetilde{u_1} = \frac{\int_{\Omega}u_1(x)dx}{|\Omega|} \ \ \text{and} \ \ \widetilde{u_2} = \frac{\int_{\Omega}u_2(x)dx}{|\Omega|}.

    By (2.2), we have u_1(x)\le \overline u_1 , where \overline u_1 = (r-\mu_1)/m , which implies \widetilde{u_1}\le \overline u_1 . Multiplying the first equation by ce^{-\gamma\tau} and adding two equations of (3.1) lead to

    -(d_1ce^{-\gamma\tau}\Delta u_1+d_2\Delta u_2) = ce^{-\gamma\tau}\left(r g(K,u_2)u_1-\mu_1 u_1-mu_1^2\right)-\mu_2 u_2.

    The integration of both sides for the above equation yields

    \widetilde{u_2} = {ce^{-\gamma\tau} \over \mu_2|\Omega|}\int_{\Omega}\left(rg(K,u_2)u_1-\mu_1 u_1-mu_1^2\right) dx\le (r-\mu_1)\overline u_1{ce^{-\gamma\tau}\over \mu_2}: = M_v.

    It is readily seen that \int_{\Omega}(u_1-\widetilde{u_1})dx = \int_{\Omega}(u_2-\widetilde{u_2})dx = 0 . Note u_1 and u_2 are bounded by two constants \overline u_1 > 0 and \overline u_2: = rc \overline u_1/\min\{\mu_1d_2/d_1, \mu_2\} > 0 by Theorem 3.1. Denote M_f = \max\limits_{u_2\in[0, \overline u_2]}|g'_{u_2}(K, u_2)| and we then obtain

    \begin{equation*} \label{ineq1} \begin{aligned} d_1\int_{\Omega}|\nabla(u_1-\widetilde{u_1})|^2dx = &\int_{\Omega}(u_1-\widetilde{u_1})(rg(K,u_2)u_1-\mu_1 u_1-mu_1^2)dx-\int_{\Omega}pu_1u_2(u_1-\widetilde{u_1})dx\\ = &\int_{\Omega}(u_1-\widetilde{u_1})\left(r g(K,u_2)u_1-\mu_1 u_1-mu_1^2-(rg(K,\widetilde{u_2})\widetilde{u_1}-\mu_1 \widetilde{u_1}-m\widetilde{u_1}^2)\right)dx\\ &+\int_{\Omega}p(\widetilde{u_1}\widetilde{u_2}-u_1u_2)(u_1-\widetilde{u_1}) dx\\ \le&\left(r-\mu_1+(rM_f+p){\overline u_1\over 2}\right)\int_{\Omega}(u_1-\widetilde{u_1})^2dx +(rM_f+p){\overline u_1\over 2}\int_{\Omega}(u_2-\widetilde{u_2})^2dx,\\ \end{aligned} \end{equation*}
    \begin{equation*} \begin{aligned} d_2\int_{\Omega}|\nabla(u_2-\widetilde{u_2})|^2dx = &\int_{\Omega}(cpe^{-\gamma\tau}u_1u_2-\mu_2 u_2)(u_2-\widetilde{u_2})dx\\ = &cpe^{-\gamma\tau}\int_{\Omega}(u_1u_2-\widetilde{u_1}\widetilde{u_2})(u_2-\widetilde{u_2})dx-\int_{\Omega}\mu_2(u_2-\widetilde{u_2})^2 dx\\ \le& cpe^{-\gamma\tau}\int_{\Omega}{M_v \over 2} (u_1-\widetilde{u_1})^2dx+cpe^{-\gamma\tau} (\overline u_1+ {M_v \over 2}) \int_{\Omega} (u_2-\widetilde{u_2})^2dx. \end{aligned} \end{equation*}

    Set A_1 = (r-\mu_1)+\left((r M_f+p)\overline u_1+cpe^{-\gamma\tau} M_v\right)/2, and A_2 = \left((r M_f+p)\overline u_1+cpe^{-\gamma\tau}(2\overline u_1+M_v)\right)/2. Then, the above inequalities and Poincar \acute{e} inequality yield that

    \begin{equation*} \label{ineq3} d_1\int_{\Omega}|\nabla(u_1-\widetilde{u_1})|^2dx+d_2\int_{\Omega}|\nabla(u_2-\widetilde{u_2})|^2dx\le \chi \int_{\Omega}\left(|\nabla(u_1-\widetilde{u_1})|^2+|\nabla(u_2-\widetilde{u_2})|^2\right)dx, \end{equation*}

    with a positive constant \chi = \max\{A_1/\sigma_1, A_2/\sigma_1\} depending on r, f, \mu_1, p, c, \gamma, \tau, \mu_2 and \sigma_1 . Hence, if \chi < \min\{d_1, d_2\} , then \nabla(u_1-\widetilde{u_1}) = \nabla(u_2-\widetilde{u_2}) = 0 , which implies (u_1, u_2) is a constant solution.

    Select u_1^*: = \nu as the bifurcation parameter and study nonhomogeneous steady state bifurcating from E^* . Let u_{2, \nu} satisfy rg(K, u_2)-pu_2 = \mu_1+m\nu , E_2 = (\nu, u_{2, \nu}) , and (\widehat u_1, \widehat u_2) = (u_1-\nu, u_2-u_{2, \nu}) . Drop \widehat\cdot . System (3.1) becomes

    \begin{equation*} \mathcal{H}(\nu,u_1,u_2) = \begin{pmatrix} d_1\Delta u_1+(u_1+\nu)(rg(K,u_2+u_{2,\nu})-\mu_1-m(u_1+\nu)-p(u_2+u_{2,\nu}))\\ d_2\Delta u_2+cpe^{-\gamma\tau}(u_1+\nu)(u_2+u_{2,\nu})-\mu_2(u_2+u_{2,\nu}) \end{pmatrix} = 0, \end{equation*}

    for (\nu, u_1, u_2)\in\mathbb{R}^+\times \mathcal{Y} with \mathcal{Y} = \{(u_1, u_2): u_1, u_2\in H^2(\Omega), (u_1)_{\nu} = (u_2)_{\nu} = 0, \ \text{on}\ \partial\Omega\} . Calculating Fr \acute{e} chet derivative of \mathcal{H} gives

    \begin{equation*} D_{(u_1,u_2)}\mathcal{H}(\nu,0,0) = \begin{pmatrix} d_1\Delta-m\nu&\nu(rg'_{u_2}(K,u_{2,\nu})-p)\\ cpe^{-\gamma\tau}u_{2,\nu}&d_2\Delta \end{pmatrix}. \end{equation*}

    Then the characteristic equation follows

    \begin{equation} \rho^2+P_i(\nu)\rho+Q_i(\nu) = 0 \ \ \text{for} \ \ i\in\mathbb{N}_0, \end{equation} (3.3)

    where

    P_i(\nu) = m\nu+(d_1+d_2)\sigma_i, \ \ Q_i(\nu) = d_1d_2\sigma_i^2+d_2 m\nu\sigma_i-\mu_2u_{2,\nu}(rg'_{u_2}(K,u_{2,\nu})-p).

    Obviously, Q_i > 0 and P_i > 0 for all \nu\in\mathbb{R}^+ and i\in\mathbb{N}_0 . Therefore, (3.3) does not have a simple zero eigenvalue. According to [4], we obtain the nonexistence of steady state bifurcation bifurcating at E_2 .

    Theorem 3.3. Model (1.4) admits no positive nonhomogeneous steady states bifurcating from E_2 .

    Next, the stability switches at E_2 and existence of periodic solutions of (1.4) bifurcating from E_2 are studied. Suppose R_0 > 1 , namely, cp(r-\mu_1) > m\mu_2 and 0\le\tau < \tau_{max}: = {1\over \gamma}\ln {cp(r-\mu_1)\over m\mu_2} to guarantee the existence of E_2 .

    Recall the stability of E_2 for \tau = 0 is proved and 0 is not the root of (2.11) for \tau\ge0 . So, we only consider eigenvalues cross the imaginary axis to the right which corresponds to the stability changes of E_2 . Now, we shall consider the positive root of G_n(\delta, \tau) . Clearly, there exists exactly one positive root of G_n(\delta, \tau) = 0 if and only if a_{0, n} < b_{0, n} for n\in\mathbb{N}_0 . More specifically,

    (\mathbf{A}_1): \ 2mu_1^* < u_2^*(p-rg'_{u_2}(K,u_2^*))

    is the sufficient and necessary condition to ensure G_0(\delta, \tau) has exactly one positive zero. For some integer n\ge1 , the assumption (\mathbf{A}_1) is a necessary condition to guarantee G_n(\delta, \tau) exists positive zeros. Set

    \begin{equation} J_n = \{\tau:\tau\in[0,\tau_{max}) \ \text{satisfies}\ a_{0,n}(\tau) < b_{0,n}(\tau)\}, \ \ n\in\mathbb{N}_0. \end{equation} (4.1)

    Implicit function theorem implies G_n(\delta, \tau) has a unique zero

    \delta_n(\tau) = \sqrt{\left(\left[b_{1,n}^2+2a_{0,n}-a_{1,n}^2+\sqrt{(b_{1,n}^2+2a_{0,n}-a_{1,n}^2)^2-4(a_{0,n}^2-b_{0,n}^2)}\right]/2\right)}

    which is a C^1 function for \tau\in J_n . Hence, i\delta_n(\tau) is an eigenvalue of (2.11), and \delta_n(\tau) satisfies

    \begin{equation} \begin{aligned} \sin(\delta_n(\tau)\tau) = &\frac{\delta_n(-\mu_2\delta_n^2+\mu_2a_{0,n}+b_{0,n}a_{1,n})} {\mu_2^2\delta_n^2+b^2_{0,n}}: = h_{1,n}(\tau),\\ \cos(\delta_n(\tau)\tau) = &\frac{b_{0,n}(\delta_n^2-a_{0,n})+a_{1,n}\mu_2\delta_n^2} {\mu_2^2\delta_n^2+b_{0,n}^2}: = h_{2,n}(\tau), \end{aligned} \end{equation} (4.2)

    for n\in\mathbb{N}_0 . Let

    \begin{eqnarray*} \label{xyz} \vartheta_n(\tau) = \left\{\begin{aligned} & \arccos h_{2,n}(\tau),\ \text{if}\ \delta_n^2 < a_{0,n}+b_{0,n}a_{1,n}/\mu_2,\\ &2\pi-\arccos h_{2,n}(\tau),\ \text{if}\ \delta_n^2\ge a_{0,n}+b_{0,n}a_{1,n}/\mu_2{,} \end{aligned}\right. \end{eqnarray*}

    which is the unique solution of \sin\vartheta_n = h_{1, n} and \cos\vartheta_n = h_{2, n} and satisfies \vartheta_n(\tau)\in(0, 2\pi] for \tau\in I_n .

    According to [3,27], we arrive at the next properties.

    Lemma 4.1. Suppose that R_0 > 1 and (\mathbf{A}_1) holds.

    (i) There exists a nonnegative integer M_1 such that J_{n}\ne\emptyset for 0\le n\le M_1 , with J_{M_1}\subset J_{M_1-1}\subset\cdots\subset J_1\subset J_0 , and J_n = \emptyset for n\ge 1+M_1 , where J_n is defined in (4.1).

    (ii) Define

    \begin{equation} \mathcal{S}_n^k(\tau) = \delta_n(\tau)\tau-\vartheta_n(\tau)-2k\pi \ \ \mathit{\text{for integer}} \ 0\le n\le M_1, \ k\in\mathbb{N}_0, \ \mathit{\text{and}} \ \tau\in J_n. \end{equation} (4.3)

    Then, \mathcal{S}_0^0(0) < 0 ; for 0\le n\le M_1 and k\in\mathbb{N}_0 , we have \lim\limits_{\tau\rightarrow \widehat\tau_n^-}\mathcal{S}_n^k(\tau) = -(2k+1)\pi , where \widehat\tau_n = \sup J_n ; \mathcal{S}_n^{k+1}(\tau) < \mathcal{S}_n^{k}(\tau) and \mathcal{S}_n^k(\tau) > \mathcal{S}_{n+1}^{k}(\tau) .

    (iii) For each integer n\in[0, N_1] and some k\in\mathbb{N}_0 , \mathcal{S}_n^k(\tau) has one positive zero \overline\tau_n\in J_n if and only if (2.11) has a pair of eigenvalues \pm i\delta_n(\overline\tau_n) . Moreover,

    \begin{equation} \mathit{\text{Sign}}({Re\lambda'(\overline\tau_n)}) = \mathit{\text{Sign}}({(\mathcal{S}_n^k)'(\overline\tau_n)}). \end{equation} (4.4)

    When (\mathcal{S}_n^k)'(\overline\tau_n) < 0 , \pm i\delta_n(\overline\tau_n) cross the imaginary axis from right to left at \tau = \overline\tau_n ; when (\mathcal{S}_n^k)'(\overline\tau_n) > 0 , \pm i\delta_n(\overline\tau_n) cross the imaginary axis from left to right.

    If \sup\limits_{\tau\in I_0}\mathcal{S}_0^0\le0 , then \mathcal{S}_n^k < 0 in J_n holds for any 0\le n\le M_1 and k\in\mathbb{N}_0 ; or only \mathcal{S}_0^0 has a zero with even multiplicity in J_0 and \mathcal{S}_n^k < 0 for any positive integers n and k . Therefore, E_2 is locally asymptotically stable for \tau\in[0, \tau_{max}) . The following assumption ensures Hopf bifurcation may occur at E_2 .

    (\mathbf{A}_2) \sup\limits_{\tau\in J_0}\mathcal{S}_0^0(\tau) > 0 and \mathcal{S}_n^k(\tau) has at most two zeros (counting multiplicity) for integer 0\le n\le M_1 and k\in \mathbb{N}_0 .

    Note \sup\limits_{\tau\in J_n} \mathcal{S}_n^0 is strictly decreasing in n due to Lemma 4.1. It then follows from (\mathbf{A}_2) , and \mathcal{S}_n^k(0) < \mathcal{S}_0^0(0) < 0 , \lim\limits_{\tau\rightarrow \widehat\tau_n^-}\mathcal{S}_n^k(\tau) < 0 for any integer 0\le n\le M_1 and k\in\mathbb{N}_0 , that we can find two positive integers

    \begin{equation} M = \{n\in[0,M_1]: \ \sup \mathcal{S}_n^0 > 0 \ \text{and} \ \sup \mathcal{S}_{n+1}^0\le0\}\ge0, \end{equation} (4.5)

    and

    \begin{equation} K_n = \{j\ge1: \ \sup \mathcal{S}_n^{j-1} > 0 \ \text{and} \ \sup \mathcal{S}_n^{j}\le0\}\ge1,\ \text{for any integer}\ 0\le n\le M. \end{equation} (4.6)

    Then \mathcal{S}_n^k(\tau) admits two simple zeros \tau_n^k and \tau_n^{2K_n-k-1} for k\in[0, K_n-1] and no zeros for k\ge K_n . The above analysis, together with Lemma 4.1(ⅲ), yields the next result.

    Lemma 4.2. Suppose R_0 > 1 and (\mathbf{A}_1) and (\mathbf{A}_2) hold. Let \mathcal{S}_n^k(\tau) , M and K_n be defined in (4.3), (4.5) and (4.6).

    (i) For integer n\in[0, M] , there are 2K_n simple zeros \tau_n^j \ (0\le j\le 2K_n-1) of \mathcal{S}_n^i(\tau) \ (0\le i\le K_n-1) , 0 < \tau_n^0 < \tau_n^1 < \tau_n^2 < \cdots < \tau_n^{2K_n-1} < \widehat\tau_n , and d\mathcal{S}_n^i(\tau_n^i)/d\tau > 0 and d\mathcal{S}_n^{i}(\tau_n^{2K_n-i-1})/d\tau < 0 for each 0\le i\le K_n-1 .

    (ii) If there exist exactly two bifurcation values \tau_{n_1}^j = \tau_{n_2}^i: = \tau_* with n_1\ne n_2 and (n_1, j), (n_2, i)\in[0, M]\times[0, 2K_n-1] , then the double Hopf bifurcation occurs at E_2 when \tau = \tau_* .

    Collect all values \tau_n^j with 0\le n\le M and 0\le j\le 2K_n-1 in a set. To ensure Hopf bifurcation occurs, remove values which appear more than once. The new set becomes

    \begin{equation} \Sigma = \{\tau_0,\tau_1,\cdots,\tau_{2L-1}\},\ \text{with}\ \tau_i < \tau_j\ \text{if}\ i < j\ \text{and}\ 1\le L\le \sum\limits_{n = 0}^{M}K_n. \end{equation} (4.7)

    Lemma 4.1(ⅱ) implies \mathcal{S}_0^0(\tau) exists two simple zeros \tau_0 < \tau_{2L-1} . When \tau = \tau_i with 0\le i\le 2L-1 , the Hopf bifurcation occurs at E_2 . Moreover, E_2 is locally asymptotically stable for \tau\in[0, \tau_0)\cup(\tau_{2L-1}, \tau_{max}) and unstable for \tau\in(\tau_0, \tau_{2L-1}) . Define

    \begin{equation} \Sigma_0 = \{\tau\in\Sigma:\mathcal{S}_0^j(\tau) = 0\ \text{for integer}\ 0\le j\le K_0\}. \end{equation} (4.8)

    Theorem 4.3. Suppose R_0 > 1 . Let J_n , \mathcal{S}_n^k(\tau) , \Sigma and \Sigma_0 be defined in (4.1), (4.3), (4.7) and (4.8), respectively.

    (i) E_2 is locally asymptotically stable for all \tau\in[0, \tau_{max}) provided that either J_0 = \emptyset or \sup\limits_{\tau\in J_0} \mathcal{S}_0^0(\tau)\le0 .

    (ii) If (\mathbf{A}_1) and (\mathbf{A}_2) hold, then a Hopf bifurcation occurs at E_2 when \tau\in \Sigma . E_2 is locally asymptotically stable for \tau\in[0, \tau_0)\cup(\tau_{2L-1}, \tau_{max}) , and unstable for \tau\in(\tau_0, \tau_{2L-1}) . Further, for \tau\in\Sigma\backslash \Sigma_0 , the bifurcating periodic solution is spatially nonhomogeneous; for \tau\in\Sigma_0 , the bifurcating periodic solution is spatially homogeneous.

    Next investigate the properties of bifurcating periodic solutions by global Hopf bifurcation theorem [28]. Set y(t) = (y_1(t), y_2(t))^T = (u_1(\cdot, \tau t)-u_1^*, u_2(\cdot, \tau t)-u_2^*)^T and write (1.4) as

    \begin{equation} y'(t) = A y(t)+ Z(y_t,\tau,Q), \ (t,\tau,Q)\in\mathbb{R}^+\times[0,\tau_{max})\times \mathbb{R}^+,\ y_t\in C([-1,0],X^2), \end{equation} (4.9)

    where y_t(\nu) = y(t+\nu) for \nu\in[-1, 0] , A = \text{diag}(\tau d_1\Delta-\tau \mu_1, \tau d_2\Delta-\tau \mu_2) and

    \begin{align*} \label{F} Z(y_t) = \tau\begin{pmatrix} (y_{1t}(0)+u_1^*)\left(rg(K,y_{2t}(0)+u_2^*)-m(y_{1t}(0)+u_1^*)-p(y_{2t}(0)+u_2^*)\right)-\mu_1 u_1^* \\ cpe^{-\gamma\tau}(y_{1t}(-1)+u_1^*)(y_{2t}(-1)+u_2^*)-\mu_2 u_2^* \end{pmatrix}. \end{align*}

    \{\Psi(t)\}_{t\ge0} denotes the semigroup yielded by A in \Omega , with Neumann boundary condition. Clearly, \lim\limits_{t\to\infty}\Psi(t) = 0 . The solution of (4.9) can be denoted by

    \begin{equation} y(t) = \Psi(t) y(0)+\int_0^t \Psi(t-\sigma) Z(y_\sigma) d\sigma. \end{equation} (4.10)

    If y(t) is a a- periodic solution of (4.9), then (4.10) yields

    \begin{equation} y(t) = \int_{-\infty}^t \Psi(t-s)Z(y_\sigma) d\sigma, \end{equation} (4.11)

    since \Psi(t+na)y(0)\to0 as n\to\infty . Hence we only need to consider (4.11). The integral operator of (4.11) is differentiable, completely continuous, and G-equivariant by [24,Chapter 6.5]. The condition \min\{d_1, d_2\} > \chi ensures (1.4) admits exactly one positive steady state E_2 . Using a similar argument as in [29,Section 4.2], (H1) (H4) in [24,Chapter 6.5] hold and we shall study the periodic solution.

    Lemma 4.4. Suppose R_0 > 1 , then all nonnegative periodic solutions of (4.9) satisfies \kappa\le u_1(x, t), u_2(x, t) \le \xi for all (x, t)\in\overline\Omega\times\mathbb{R}^+ , where \xi and \kappa are defined in Theorems 2.2 and 2.4, respectively.

    Lemma 4.4 can be obtained by Theorems 2.2 and 2.4. We further assume

    (\mathbf{A}_3): \ {g(K,u_2)-g(K,u_2^*) \over u_1-u_1^*}-{m\over r}\le0 \ \ \text{for} \ \ u_1\in[\kappa, {r-\mu_1 \over m}] \ \ \text{and} \ \ u_2\in[\kappa,\xi].

    This technical condition is used to exclude the \tau- periodic solutions. Note assumption (\mathbf{A}_3) holds when K = 0 . This, together with that (g(K, u_2)-g(K, u_2^*))/(u_1-u_1^*) is continuous in K , implies there exists \varepsilon > 0 such that (\mathbf{A}_3) holds for 0\le K < \varepsilon , that is, (\mathbf{A}_3) holds for the model (1.4) with weak fear effect.

    Lemma 4.5. Assume that R_0 > 1 and (\mathbf{A}_3) holds, then model (1.4) admits no nontrivial \tau- periodic solution.

    Proof. Otherwise, let (u_1, u_2) be the nontrivial \tau- periodic solution, that is, (u_1(x, t-\tau), u_2(x, t-\tau)) = (u_1(x, t), u_2(x, t)) . Thus, we have

    \begin{equation} \begin{aligned} &\partial_t u_1 = d_1\Delta u_1+u_1\left(rg(K,u_2)-\mu_1-mu_1-pu_2\right), \ x\in\Omega, t > 0,\\ &\partial_t u_2 = d_2\Delta u_2+cpe^{-\gamma\tau}u_1u_2-\mu_2u_2, \ x\in\Omega, t > 0{,}\\ &\partial_\nu u_1 = \partial_\nu u_2 = 0,x\in\partial\Omega,t > 0{,}\\ &u_1(x,\vartheta) = u_{10}(x,\vartheta)\ge0,u_2(x,\vartheta) = u_{20}(x,\vartheta)\ge0,x\in\Omega,\vartheta\in[-\tau,0]. \end{aligned} \end{equation} (4.12)

    Claim

    (u_1,u_2)\to E_2 \ \ \text{as} \ \ t\to\infty.

    To see this, establish the Lyapunov functional \mathbb{L}_1: C(\bar\Omega, \mathbb{R}^+\times\mathbb{R}^+)\to\mathbb{R} ,

    \mathbb{L}_1(\phi_1,\phi_2) = \int_{\Omega}\left(ce^{-\gamma\tau}(\phi_1-u_1^*\ln\phi_1)+(\phi_2-u_2^*\ln\phi_2)\right)dx \ \text{for} \ (\phi_1,\phi_2)\in C(\bar\Omega,\mathbb{R}^+\times\mathbb{R}^+).

    Along the solution of system (4.12), the time derivative of \mathbb{L}_1(\phi_1, \phi_2) is

    {d\mathbb{L}_1 \over dt} = \int_{\Omega}\left[-{d_1 \mu_2|\nabla u_1|^2 \over pu_1^2}-d_2 u_2^*{|\nabla u_2|^2 \over u_2^2} +r(u_1-u_1^*)^2\left({g(K,u_2)-g(K,u_2^*) \over u_1-u_1^*}-{m \over r}\right)\right]dx.

    The assumption ( \mathbf{A}_3 ) ensures d\mathbb{L}_1/dt\le0 for all (u_1, u_2)\in C(\bar\Omega, \mathbb{R}^+\times\mathbb{R}^+) . The maximal invariant subset of d\mathbb{L}_1/dt = 0 is \{E_2\} . Therefore, E_2 attracts all positive solution of (4.12) by LaSalle-Lyapunov invariance principle [21,22] which excludes the nonnegative nontrivial \tau- periodic solution.

    To obtain the nonexistence of \tau- periodic solution for model (1.4), we must use the condition (\mathbf{A}_3) which is very restrictive. However, in numerical simulations, Lemma 4.5 remains true even (\mathbf{A}_3) is violated. Thus, we conjecture the nonexistence of \tau -periodic solution for (1.4).

    In the beginning of this section, when \tau = \tau_i with 0\le i\le 2L-1 , \pm i \delta_n(\tau_i) are a pair of eigenvalues of (2.11). Give the next standard notations:

    (i) For 0\le i\le 2L-1 , (E_2, \tau_i, 2\pi/(\delta_n(\tau_i)\tau_i)) is an isolated singular point.

    (ii) \Gamma = \mathcal{C}l\{(y, \tau, Q)\in X^2\times\mathbb{R}^+\times\mathbb{R}^+:\ y\ \text{is the nontrivial Q-periodic solution of (4.9)}\} is a closed subset of X^2\times\mathbb{R}^+\times\mathbb{R}^+ .

    (iii) For 0\le i\le 2L-1 , \mathcal{C}_i(E_2, \tau_i, Q_i) is the connected component of (E_2, \tau_i, Q_i) in \Gamma .

    (iv) For integer 0\le k\le \max\limits_{n\in[0, M]} K_n-1 , let \Sigma_H^k = \{\tau\in\Sigma: \mathcal{S}_n^k(\tau) = 0\ \text{for integer}\ 0\le n\le M\}.

    We are ready to present a conclusion on the global Hopf branches by a similar manner in [30,Theorem 4.12].

    Theorem 4.6. Assume R_0 > 1 , \min\{d_1, d_2\} > \chi and (\mathbf{A}_1) (\mathbf{A}_3) hold. Then we have the following results.

    (i) The global Hopf branch \mathcal{C}_i(E_2, \tau_i, Q_i) is bounded for \tau_i\in\Sigma_H^k with k\ge1 and i\in[0, 2L-1] .

    (ii) For any \tau\in(\min\limits_{k}\Sigma_H^k, \max\limits_{k}\Sigma_H^k) , model (1.4) possesses at least one periodic solution.

    (iii) For \tau_{i_1}\in\Sigma_H^{k_1} , \tau_{i_2}\in\Sigma_H^{k_2} with i_1, i_2\in[0, 2L-1] and k_1, k_2\in[0, \max\limits_{n\in[0, M]} K_n-1] , we have \mathcal{C}_{i_1}(E_2, \tau_{i_1}, Q_{i_1})\cap\mathcal{C}_{i_2}(E_2, \tau_{i_2}, Q_{i_2}) = \emptyset if k_1\ne k_2 .

    To verify obtained theoretical results, the numerical simulation is presented in this part. We choose the fear effect function as g(K, u_2) = e^{-K u_2} and let

    \begin{equation*} \label{par1} \Omega = (0,2\pi),d_1 = d_2 = 1, r = 8, \mu_1 = 0.1, a = 0.2, p = 1, \gamma = 0.3, \mu_2 = 0.2, c = 0.1. \end{equation*}

    Figure 1 shows the existence and stability of E_0 , E_1 and E_2 , and Hopf bifurcation curve of model (1.4) in K-\tau plane. Above the line \tau = \tau_{max} , E_2 does not exist, E_1 is globally asymptotically stable and E_0 is unstable; Below the line \tau = \tau_{max} , there exist three constant steady states E_0, E_1 and E_2 . In the region which is bounded by \tau = \tau_{max} and 2mu_1^*+u_2^*(rg'_{u_2}(K, u_2^*)-p) = 0 , no Hopf bifurcation occurs and E_2 is stable. In the region which is bounded by \tau = \tau_0 and \tau = \tau_1 , there exist periodic solutions through Hopf bifurcation bifurcating at E_2

    Figure 1.  Basic dynamics of model (1.4) in different regions with K-\tau plane.

    Fix K = 1.07 , then by simple calculation, we have \tau_{max}\approx9.944 , \tau_0\approx0.85 , \tau_1\approx3.55 , J_0 = [0, 5.25] , J_1 = [0, 3.05] , J_n = \emptyset for n\ge2 , \sup\limits_{\tau\in J_0} \mathcal{S}_0^0 > 0 and \sup\limits_{\tau\in J_1}\mathcal{S}_1^0 < 0 . Thus, all Hopf bifurcation values \tau_0 and \tau_1 are the all zeros of S_0^k(\tau) for integer k\ge0 . We summarize the dynamics of model (1.4) as follows.

    (i) For \tau\in[\tau_{max}, \infty) , we obtain E_1 is globally asymptotically stable and E_0 is unstable, see Figure 2(a).

    Figure 2.  (a) \tau = 10 > \tau_{max} , E_1 is globally asymptotically stable. (b) \tau = 0.5\in(0, \tau_0) , E_2 is locally asymptotically stable. (c) \tau = 2\in(\tau_0, \tau_1) , a homogeneous periodic solution emerges.

    (ii) For \tau\in(0, \tau_0)\cup(\tau_1, \tau_{max}) , we obtain E_2 is locally asymptotically stable, and two constant steady state E_0 and E_1 are unstable, as shown in Figure 2(b).

    (iii) For \tau\in(\tau_0, \tau_1) , we obtain E_0, E_1 and E_2 are unstable, a periodic solution bifurcates from E_2 , as shown in Figure 2(c). Further, a Hopf bifurcation occurs at E_2 when \tau = \tau_0, \text{and}\ \tau_1 .

    Next, we explore the global Hopf branches by choosing \Omega = (0, 4\pi) and

    \begin{equation} d_1 = 1,d_2 = 1,r = 10,\mu_1 = 5,a = 0.4,p = 1,\gamma = 0.05,\mu_2 = 4.75,c = 2.5,K = 0.4. \end{equation} (5.1)

    As shown in Figure 3, we collect all zeros of \mathcal{S}_n^k(\tau) for nonnegative n, k in set B_i with i = 0, 1, 2, 3 , namely,

    \begin{aligned} B_0& = \{0.06,1.88,3.86,6.14,9.54,11.8,13.36,13.8,13.98,14.04\},\\ B_1& = \{0.08,1.96,4.04,6.56,12.36,13.01,13.25,13.35\},\\ B_2& = \{0.15,2.31,4.96,10.16,10.85,11.07\}, \ B_3 = \{0.39,6.45\}. \end{aligned}
    Figure 3.  The graphs of \mathcal{S}_n^k(\tau) with integers 0\le n\le 3 and 0\le k\le 4 .

    From Theorem 4.3, E_2 is locally asymptotically stable when \tau\in(0, 0.06)\cup(14.04, \tau_{max}) and unstable when \tau\in(0.06, 14.04) , at least one periodic solution emerges for \tau\in(0.06, 14.04) . Moreover, a spatially homogeneous periodic solution bifurcates from \tau\in B_0 , see Figure 4(a); a spatially nonhomogeneous periodic solution bifurcates from \tau\in B_1\cup B_2\cup B_3 , see Figure 4(b).

    Figure 4.  (a) \tau = 1.88\in B_0 , the bifurcating periodic solution is spatially homogeneous. (b) \tau = 2.68\in(\tau_4, \tau_5) , the bifurcating periodic solution is spatially nonhomogeneous.

    In model (1.3), the kernel function takes form as \mathcal{G} = e^{-\gamma\tau}f(x-y) by reasonable assumptions and our theoretical results are derived by choosing f(x-y) as Dirac-delta function. Next, we choose

    \begin{equation} f = \frac{e^{-2|x-y|^2}}{\int_\Omega e^{-2|x-y|^2}dy}. \end{equation} (5.2)

    Here, f is the truncated normal distribution. Clearly, \int_\Omega f(x-y) dy = 1 . Let \Omega = (0, 4\pi) and the parameter values chosen according to (5.1). As shown in Figure 5, when \tau = 1.03 , a stable nonhomogeneous periodic solution emerges; when \tau = 1.5 , a homogeneous periodic solution emerges; when \tau = 16 , the positive constant steady state is stable. Numerical simulation suggests nonlocal interaction can produce more complex dynamics.

    Figure 5.  For the nonlocal model (1.3) with \mathcal{G} defined in (5.2), delay \tau induces different dynamics.

    We formulate an age-structured predator-prey model with fear effect. For R_0\le1 , the global asymptotic stability for predator-free constant steady state is proved via Lyapunov-LaSalle invariance principle. For R_0 > 1 , we prove the nonexistence of spatially nonhomogeneous steady states and exclude steady state bifurcation. Finally, we carry out Hopf bifurcation analyses and prove global Hopf branches are bounded.

    Our theoretical results are obtained by choosing a special kernel function in model (1.3). However, in numerical results, we explore rich dynamics when the nonlocal interaction is incorporated into the delayed term. The theoretical results concerning the nonlocal model are left as an open problem.

    H. Shu was partially supported by the National Natural Science Foundation of China (No. 11971285), the Fundamental Research Funds for the Central Universities (No. GK202201002), the Natural Science Basic Research Program of Shaanxi (No. 2023-JC-JQ-03), and the Youth Innovation Team of Shaanxi Universities. W. Xu was partially supported by a scholarship from the China Scholarship Council while visiting the University of New Brunswick. P. Jiang was partially supported by the National Natural Science Foundation of China (No. 72274119).

    The authors declare no conflicts of interest in this paper.



    [1] Liberini M, De Falco G, Scherillo F, et al. (2016) Nano-TiO2 coatings on aluminum surfaces by aerosol flame synthesis. Thin Solid Films 609: 53–61. https://doi.org/10.1016/j.tsf.2016.04.025 doi: 10.1016/j.tsf.2016.04.025
    [2] Immarigeon JP, Holt RT, Koul AK, et al. (1995) Lightweight materials for aircraft applications. Mater Charact 35: 41–67. https://doi.org/10.1016/1044-5803(95)00066-6 doi: 10.1016/1044-5803(95)00066-6
    [3] Liu T, Zhang F, Xue C, et al. (2010) Structure stability and corrosion resistance of nano-TiO2 coatings on aluminum in seawater by a vacuum dip-coating method. Surf Coat Technol 205: 2335–2339. https://doi.org/10.1016/j.surfcoat.2010.09.028 doi: 10.1016/j.surfcoat.2010.09.028
    [4] Soklic A, Tasbihi M, Kete M, et al. (2015) Deposition and possible influence of a self-cleaning thin TiO2/SiO2 film on a photovoltaic module efficiency. Catal Today 252: 54–60. https://doi.org/10.1016/j.cattod.2014.10.021 doi: 10.1016/j.cattod.2014.10.021
    [5] Çomakli O, Yazici M, Yetim T, et al. (2017) The effects of aging time on the structural and electrochemical properties of composite coatings on CP-Ti substrate. J Bionic Eng 14: 532–539. https://doi.org/10.1016/S1672-6529(16)60419-5 doi: 10.1016/S1672-6529(16)60419-5
    [6] Gobara M (2015) Effects of TiO2/SiO2 reinforced nanoparticles on the mechanical properties of green hybrid coating. Int Lett Chem Phys Astron 47: 56–66. https://doi.org/10.56431/p-z14m86 doi: 10.56431/p-z14m86
    [7] Krishna V, Padmapreetha R, Chandrasekhar SB, et al. (2019) Oxidation resistant TiO2-SiO2 coatings on mild steel by sol-gel. Surf Coat Technol 378: 125041. https://doi.org/10.1016/j.surfcoat.2019.125041 doi: 10.1016/j.surfcoat.2019.125041
    [8] Khosravi SH, Veerapandiyan VK, Vallant R, et al. (2020) Effect of processing conditions on the structural properties and corrosion behavior of TiO2-SiO2 multilayer coatings derived via the sol-gel method. Ceram Int 46: 17741–17751. https://doi.org/10.1016/j.ceramint.2020.04.079 doi: 10.1016/j.ceramint.2020.04.079
    [9] Jacobs M, De Vos Y, Middelkoop V (2021) Thickness controlled SiO2/TiO2 sol-gel coating by spraying. Open Ceram 6: 100121. https://doi.org/10.1016/j.oceram.2021.100121 doi: 10.1016/j.oceram.2021.100121
    [10] Widati AA, Nuryono N, Kartini I (2019) Water-repellent glass coated with SiO2-TiO2-methyltrimethoxysilane through sol-gel coating. AIMS Mater Sci 6: 10–24. https://doi.org/10.3934/matersci.2019.1.10 doi: 10.3934/matersci.2019.1.10
    [11] Zhang L, Zheng Q, Yin L, et al. (2021) Surface passivation of applying an organic-inorganic hybrid coatings toward significantly chemically stable iron powder. Colloid Surface A 610: 125910. https://doi.org/10.1016/j.colsurfa.2020.125910 doi: 10.1016/j.colsurfa.2020.125910
    [12] Zhang L, Wan W, Jiang X, et al. (2022) Enhancement of oxidation and corrosion resistance of flaky carbonyl‑iron powder via SiO2/KH560/PDMS coating applied with sol-gel. Surf Coat Technol 437: 128346. https://doi.org/10.1016/j.surfcoat.2022.128346 doi: 10.1016/j.surfcoat.2022.128346
    [13] Zhang L, Wang B, Jiang X, et al. (2022) Silicone-encapsulated carbonyl iron filler for corrosion-resistant electromagnetic shielding. Mater Chem Phys 282: 125918. https://doi.org/10.1016/j.matchemphys.2022.125918 doi: 10.1016/j.matchemphys.2022.125918
    [14] Oliver WC, Pharr GM (2004) Measurement of hardness and elastic modulus by instrumented indentation: Advances in understanding and refinements to methodology. J Mater Res 19: 3–20. https://doi.org/10.1557/jmr.2004.19.1.3 doi: 10.1557/jmr.2004.19.1.3
    [15] Kalidindi SR, Pathak S (2008) Determination of the effective zero-point and the extraction of spherical nanoindentation stress–strain curves. Acta Mater 56: 3523–3532. https://doi.org/10.1016/j.actamat.2008.03.036 doi: 10.1016/j.actamat.2008.03.036
    [16] Alaboodi AS, Hussain Z (2019) Finite element modeling of nano-indentation technique to characterize thin film coatings. J King Saud Univ Eng Sci 31: 61–69. https://doi.org/10.1016/j.jksues.2017.02.001 doi: 10.1016/j.jksues.2017.02.001
    [17] Jimenez-Pique E, Gonzalez-Garcia L, Gonzalez-Elipe AR, et al. (2014) Nanoindentation of nanocolumnar TiO2 thin films with single and stacked zig-zag layers. Thin Solid Films 550: 444–449. https://doi.org/10.1016/j.tsf.2013.10.022 doi: 10.1016/j.tsf.2013.10.022
    [18] Bressan JD, Tramontin A, Rosa C (2005) Modeling of nanoindentation of bulk and thin film by finite element method. Wear 258: 115–122. https://doi.org/10.1016/j.wear.2004.05.021 doi: 10.1016/j.wear.2004.05.021
    [19] Zhang W (2017) Mechanical characterization of YBCO thin films using nanoindentation and finite element method. Int J Mater Res 108: 732–740. https://doi.org/10.3139/146.111533 doi: 10.3139/146.111533
    [20] Mocko W, Szymanska M, Smietana M, et al. (2014) Simulation of nanoindentation experiments of single-layer and double-layer thin films using finite element method. Surf Interface Anal 46: 1071–1076. https://doi.org/10.1002/sia.5473 doi: 10.1002/sia.5473
    [21] Gupta AK, Porwal D, Dey A, et al. (2016) Evaluation of critical depth ratio for soft V2O5 film on hard Si substrate by finite element modeling of experimentally measured nanoindentation response. J Phys D Appl Phys 49: 155302. DOI 10.1088/0022-3727/49/15/155302 doi: 10.1088/0022-3727/49/15/155302
    [22] Cheng SW, Chen BS, Jian SR, et al. (2022) Finite element analysis of nanoindentation responses in Bi2Se3. Coatings 12: 1554. https://doi.org/10.3390/coatings12101554 doi: 10.3390/coatings12101554
    [23] Lichinchi M, Lenardi C, Haupt J, et al. (1998) Simulation of Berkovich nanoindentation experiments on thin films using finite element method. Thin Solid Films 312: 240–248. https://doi.org/10.1016/S0040-6090(97)00739-6 doi: 10.1016/S0040-6090(97)00739-6
    [24] Wang K, Ma Q, Xu L, et al. (2023) Determining the elastic–plastic properties of materials with residual stress included using nanoindentation experiments and dimensionless functions. Eng Fract Mech 282: 109175. https://doi.org/10.1016/j.engfracmech.2023.109175 doi: 10.1016/j.engfracmech.2023.109175
    [25] Noii N, Aghayan I (2019) Characterization of elastic-plastic coated material properties by indentation techniques using optimisation algorithms and finite element analysis. Int J Mech Sci 152: 465–480. https://doi.org/10.1016/j.ijmecsci.2019.01.010 doi: 10.1016/j.ijmecsci.2019.01.010
    [26] Kang JJ, Becker AA, Sun W (2012) Determining elastic–plastic properties from indentation data obtained from finite element simulations and experimental results. Int J Mech Sci 62: 34–46. https://doi.org/10.1016/j.ijmecsci.2012.05.011 doi: 10.1016/j.ijmecsci.2012.05.011
    [27] Kang JJ, Becker AA, Wen W, et al. (2018) Extracting elastic-plastic properties from experimental loading-unloading indentation curves using different optimization techniques. Int J Mech Sci 144: 102–109. https://doi.org/10.1016/j.ijmecsci.2018.05.043 doi: 10.1016/j.ijmecsci.2018.05.043
    [28] Shaoming M, Sun Y, Wang H, et al. (2017) Effect of a minor Sr modifier on the microstructures and mechanical properties of 7075 T6 Al alloys. Metals 7: 13. https://doi.org/10.3390/met7010013 doi: 10.3390/met7010013
    [29] Zhou B, Lui B, Zhang S, et al. (2021) Microstructure evolution of recycled 7075 aluminum alloy and its mechanical and corrosion properties. J Alloys Compd 879: 160407. https://doi.org/10.1016/j.jallcom.2021.160407 doi: 10.1016/j.jallcom.2021.160407
    [30] Pastor A, Svoboda G (2013) Time-evolution of heat affected zone (HAZ) of friction stir welds of AA7075-T651. J Mater Phys Chem 1: 58–64. https://doi.org/10.12691/jmpc-1-4-1 doi: 10.12691/jmpc-1-4-1
    [31] Krzak-Ros J, Filipiak J, Pezowicz C, et al. (2009) The effect of substrate roughness on the surface structure of TiO2, SiO2, and doped thin films prepared by the sol-gel method. Acta Bioeng Biomech 11: 21–9.
    [32] Chuang LC, Luo CH (2011) Nanomechanical properties of prepared-TiO2 films using nanoindentation technique. Adv Mat Res 214: 388–391. https://doi.org/10.4028/www.scientific.net/AMR.214.388 doi: 10.4028/www.scientific.net/AMR.214.388
    [33] Bastakys L, Marcinauskas L, Milieska M, et al. (2023) Tribological properties of Cr2O3, Cr2O3–SiO2-TiO2 and Cr2O3–SiO2-TiO2-graphite coatings deposited by atmospheric plasma spraying. Coatings 13: 408. https://doi.org/10.3390/coatings13020408 doi: 10.3390/coatings13020408
    [34] Chuang LC, Luo CH, Yang SH (2011) The structure and mechanical properties of thick rutile-TiO2 films using different coating treatments. Appl Surf Sci 258: 297–303. https://doi.org/10.1016/j.apsusc.2011.08.055 doi: 10.1016/j.apsusc.2011.08.055
    [35] Boshkova N, Stambolova I, Stoyanova D, et al. (2023) Protective characteristics of TiO2 sol-gel layer deposited on Zn-Ni or Zn-Co substrates. Coatings 13: 295. https://doi.org/10.3390/coatings13020295 doi: 10.3390/coatings13020295
    [36] Rivero PJ, Maeztu JD, Berlanga C, et al. (2018) Hydrophobic and corrosion behavior of sol-gel hybrid coatings based on the combination of TiO NPs and fluorinated chains for aluminum alloys protection. Metals 8: 1076. https://doi.org/10.3390/met8121076 doi: 10.3390/met8121076
    [37] Shadravan A, Sadeghian Z, Nemati A, et al. (2015) Corrosion protection of 1050 aluminum alloy using a smart self-cleaning TiO2-CNT coating. Surf Coat Technol 275: 224–231. https://doi.org/10.1016/j.surfcoat.2015.05.015 doi: 10.1016/j.surfcoat.2015.05.015
  • This article has been cited by:

    1. Helani Perera, Shalinda Fernando, Miyuru B. Gunathilake, T. A. J. G. Sirisena, Upaka Rathnayake, Antonio Donateo, Evaluation of Satellite Rainfall Products over the Mahaweli River Basin in Sri Lanka, 2022, 2022, 1687-9317, 1, 10.1155/2022/1926854
    2. Miyuru B. Gunathilake, M. N. M. Zamri, Tharaka P. Alagiyawanna, Jayanga T. Samarasinghe, Pavithra K. Baddewela, Mukand S. Babel, Manoj K. Jha, Upaka S. Rathnayake, Hydrologic Utility of Satellite-Based and Gauge-Based Gridded Precipitation Products in the Huai Bang Sai Watershed of Northeastern Thailand, 2021, 8, 2306-5338, 165, 10.3390/hydrology8040165
    3. Yan Li, Teng Ma, Yan Wang, 2022, Chapter 107, 978-3-031-05483-9, 822, 10.1007/978-3-031-05484-6_107
    4. Jingjing Huang, Rahim Khan, Personalized College English Learning Based on Deep Learning under the Background of Big Data, 2022, 2022, 1687-5273, 1, 10.1155/2022/7361746
    5. Juhar Mohammed, Yenesew Mengiste, Vijay P. Singh, Improving spatio-temporal precipitation estimates in data scarce river basins: an application of machine learning-based multi-source data merging, 2022, 1436-3240, 10.1007/s00477-022-02346-4
    6. Sanjana De Zoysa, Helani Perera, Miyuru B. Gunathilake, Upaka Rathnayake, Development of rainfall intensity-duration-frequency curves for the Fiji Islands: integrating TRMM-3B42 and measured gauge data with future projections, 2023, 35, 2766-9645, 361, 10.1080/27669645.2023.2278827
    7. Paul Muñoz, David F. Muñoz, Johanna Orellana-Alvear, Rolando Célleri, Enhancing runoff forecasting through the integration of satellite precipitation data and hydrological knowledge into machine learning models, 2024, 0921-030X, 10.1007/s11069-024-06939-w
    8. Fuzeng Wang, Xuejiao An, Qiusong Wang, Zixin Li, Lin Han, Debin Su, Research on a Rainfall Prediction Model in Guizhou Based on Raindrop Spectra, 2024, 15, 2073-4433, 495, 10.3390/atmos15040495
    9. Sanjana De Zoysa, Jeewanthi Sirisena, Helani Perera, Shalinda Fernando, Miyuru Gunathilake, Upaka Rathnayake, Development of intensity-duration-frequency curves for Sri Lanka using satellite-based precipitation products – Understanding environmental conditions and concerns, 2024, 9, 26660164, 100713, 10.1016/j.cscee.2024.100713
    10. Gilbert Hinge, Swati Sirsant, Amandeep Kumar, Ruchika Gupta, Mohamed A. Hamouda, Enhancing flood prediction in Southern West Bengal, India using ensemble machine learning models optimized with symbiotic organisms search algorithm, 2024, 1436-3240, 10.1007/s00477-024-02712-4
    11. Lorenza Ceferino-Hernández, Francisco Magaña-Hernández, Enrique Campos-Campos, Gabriela Adina Morosanu, Carlos E. Torres-Aguilar, René Sebastián Mora-Ortiz, Sergio A. Díaz, Assessment of PERSIANN Satellite Products over the Tulijá River Basin, Mexico, 2024, 16, 2072-4292, 2596, 10.3390/rs16142596
  • Reader Comments
  • © 2024 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(1662) PDF downloads(149) Cited by(0)

Figures and Tables

Figures(10)  /  Tables(3)

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog