Research article Special Issues

A multiscale perspective on the kinetics of solid state transformations with application to bainite formation

  • Received: 29 May 2015 Accepted: 27 August 2015 Published: 17 September 2015
  • We give an excerpt of recent developments in the experimentally benchmarked modeling of bainite formation in the press hardening process. As the press hardening process poses a heavily multi-parameter dependent modeling challenge, we focus on three main branches which complement each other. We emphasise the combination of basic sharp interface and phase field models with pragmatically adapted multi phase field models and experimentally parametrized implementations of the Johnson-Mehl-Avrami model. In the basic thermodynamic modeling part, we review fundamental aspects of displacive and diffusional-displacive transformations to predict dominant transformation morphologies. These results provide a link to multi-phase-field implementations which allow to simulate isothermal bainitic transformations, supported by available material data from thermodynamic databases. Excellent agreement with experiments, e.g. scanning electron microscopy for the transformed bainite in the high-carbon steel 100Cr6 shows the value of these model implementations. The further connection to Johnson-Mehl-Avrami models offers to extend the understanding to transformation plasticity for the press hardening steel 22MnB5.

    Citation: Claas Hüter, Mingxuan Lin, Diego Schicchi, Martin Hunkel, Ulrich Prahl, Robert Spatschek. A multiscale perspective on the kinetics of solid state transformations with application to bainite formation[J]. AIMS Materials Science, 2015, 2(4): 319-345. doi: 10.3934/matersci.2015.4.319

    Related Papers:

    [1] Minseok Kim, Yeongjong Kim, Yeoneung Kim . Physics-informed neural networks for optimal vaccination plan in SIR epidemic models. Mathematical Biosciences and Engineering, 2025, 22(7): 1598-1633. doi: 10.3934/mbe.2025059
    [2] Yoon-gu Hwang, Hee-Dae Kwon, Jeehyun Lee . Feedback control problem of an SIR epidemic model based on the Hamilton-Jacobi-Bellman equation. Mathematical Biosciences and Engineering, 2020, 17(3): 2284-2301. doi: 10.3934/mbe.2020121
    [3] Dan Shi, Shuo Ma, Qimin Zhang . Sliding mode dynamics and optimal control for HIV model. Mathematical Biosciences and Engineering, 2023, 20(4): 7273-7297. doi: 10.3934/mbe.2023315
    [4] Kento Okuwa, Hisashi Inaba, Toshikazu Kuniya . Mathematical analysis for an age-structured SIRS epidemic model. Mathematical Biosciences and Engineering, 2019, 16(5): 6071-6102. doi: 10.3934/mbe.2019304
    [5] Holly Gaff, Elsa Schaefer . Optimal control applied to vaccination and treatment strategies for various epidemiological models. Mathematical Biosciences and Engineering, 2009, 6(3): 469-492. doi: 10.3934/mbe.2009.6.469
    [6] Toshikazu Kuniya, Yoshiaki Muroya, Yoichi Enatsu . Threshold dynamics of an SIR epidemic model with hybrid of multigroup and patch structures. Mathematical Biosciences and Engineering, 2014, 11(6): 1375-1393. doi: 10.3934/mbe.2014.11.1375
    [7] Yoichi Enatsu, Yukihiko Nakata, Yoshiaki Muroya . Global stability for a class of discrete SIR epidemic models. Mathematical Biosciences and Engineering, 2010, 7(2): 347-361. doi: 10.3934/mbe.2010.7.347
    [8] Abdennasser Chekroun, Mohammed Nor Frioui, Toshikazu Kuniya, Tarik Mohammed Touaoula . Global stability of an age-structured epidemic model with general Lyapunov functional. Mathematical Biosciences and Engineering, 2019, 16(3): 1525-1553. doi: 10.3934/mbe.2019073
    [9] Rinaldo M. Colombo, Mauro Garavello . Optimizing vaccination strategies in an age structured SIR model. Mathematical Biosciences and Engineering, 2020, 17(2): 1074-1089. doi: 10.3934/mbe.2020057
    [10] Mostafa Adimy, Abdennasser Chekroun, Claudia Pio Ferreira . Global dynamics of a differential-difference system: a case of Kermack-McKendrick SIR model with age-structured protection phase. Mathematical Biosciences and Engineering, 2020, 17(2): 1329-1354. doi: 10.3934/mbe.2020067
  • We give an excerpt of recent developments in the experimentally benchmarked modeling of bainite formation in the press hardening process. As the press hardening process poses a heavily multi-parameter dependent modeling challenge, we focus on three main branches which complement each other. We emphasise the combination of basic sharp interface and phase field models with pragmatically adapted multi phase field models and experimentally parametrized implementations of the Johnson-Mehl-Avrami model. In the basic thermodynamic modeling part, we review fundamental aspects of displacive and diffusional-displacive transformations to predict dominant transformation morphologies. These results provide a link to multi-phase-field implementations which allow to simulate isothermal bainitic transformations, supported by available material data from thermodynamic databases. Excellent agreement with experiments, e.g. scanning electron microscopy for the transformed bainite in the high-carbon steel 100Cr6 shows the value of these model implementations. The further connection to Johnson-Mehl-Avrami models offers to extend the understanding to transformation plasticity for the press hardening steel 22MnB5.


    The susceptible-infective-removed-susceptible (SIRS) model is one of the most popular epidemic models, in which total host population is divide into three classes called $ susceptible\; (S) $, $ infective\; (I) $ and $ removed\; (R) $. Since the Kermack and McKendrick pioneering work [1], a large number of studies have been carried out on the model of continuous infectious diseases attempting to gain a better understanding of diseases transmission, especially for the control policies and dynamics [2,3,4,5]. Recent years, many researchers have studied the epidemic with SIRS models such as COVID-19(2019), hand-foot-and-mouth diseases(HFMD) and other infectious diseases recently. Although scholars have flocked to study many properties of the SIRS models, there are still many aspects that deserve further study, such as the threshold dynamics of SIRS model with vaccination and optimal control using dynamic programming method.

    In studying the dynamics of the disease, one of the most important concepts is the basic reproduction number $ R_0 $. It is epidemiologically defined as the expected number of secondary cases produced by a typical infectious individual during its entire infectious period in a completely susceptible host population. To investigate the global behavior of the prevalence of infectious diseases, the stability analysis of equilibrium for epidemic models have been carried out (see [6,7] and the references therein). For ordinary differential SIRS models, Mena-Lorca and Hethcote considered several kinds of SIRS epidemic models and a threshold parameters were also found in [6] to determine whether the disease dies out or approaches to an endemic equilibrium, we also note that the global stability of SIRS epidemic model have already been investigated and the threshold theorems are well obtained. In the subsequent studies, it was generally believed that the introduction of age-structure into epidemic models is necessary and reasonable. Since the fact that the age structure of a population affects the dynamics of disease transmission was recognized, various age-structured epidemic models have been investigated by many authors. Some theoretics on the transmission dynamics of age-structured epidemic models have been developed in [8,9,10,11]. He et al. [12] studied the optimal birth control of age-dependent competitive species. Yang et al. [13] proposed epidemic model with age-since-infection and diffusion, the next generation operator $ R $ is also given as the basic regeneration numbers. In paper [14], Yang et al. also studied the threshold dynamics by the method of Lyapunov functionals for an age-of-infection epidemiological model with nonlinear incidence rate. Chekroun et al. [15] showed the stability of the equilibrium by using Lyapunov functions. Lan et al. [16] just studied the impact of hospital resources and environmental perturbation to the dynamics of SIRS model. But the effect of vaccination on threshold dynamics was not considered in these articles.

    In fact, we notice that vaccination has a significant impact on the thresholds dynamics of many infectious diseases, such as COVID-19, HPV, etc. Obviously, vaccination of susceptible individuals can greatly reduce the probability that susceptible individuals become infected. Thus, decision makers hope to achieve the optimal proportion of vaccines with the lowest cost. Optimal control is a promising strategy to control disease outbreaks. The main purpose of using the control in some diseases is to search for the most effective one that reduces the infection to a minimum level while minimizing the cost of applying control measures. In papers [7,9,17], the strategies such as health promotion campaigns or lockdown policies can contribute to reducing the disease transmission. In paper [18], Mu and Zhang investigated the near-optimal control of SIRS multi-strain epidemic model with age-structure using Pontryagin maximum principle. Bolzoni et al. [19] investigated optimal control of epidemic size and duration with limited resources. Zhou et al. [20] researched optimal isolation strategies of emerging infectious diseases with limited resources. In papers [21,22,23], some measures correspond to screening and quarantining of infected were implemented to achieve the purpose of controlling the disease.

    However, we note that the most of these articles by using maximum principle to show optimal controls and only present necessary conditions. Defectively, the method of maximum principle only addresses the optimal control problem with initial value case of $ t_0 = 0 $. In fact, vaccination can be administered at any time $ t = t_0 $ in the control of the disease. Surpassingly, the dynamic programming method can give the necessary and sufficient conditions for the existence of optimal control via the HJB equation at any initial times $ t = t_0 $ and initial states. In other words, the dynamic programming method is an extension of the maximum principle, namely, the dynamic programming method is equivalent to maximum principle method when $ t_0 = 0 $.

    Overall, the threshold dynamics results about age-structured SIRS epidemic models with the proportion of the vaccine injections are comparatively scarce, as well as few articles show that dynamic programming approach to solve the optimal control problem. Thus, in this paper, taking into account for the effect of the proportion of the vaccine injections, we incorporate the vaccine injection rate in solving the threshold of age-structured SIRS epidemic models. And vaccination optimal control is also obtained through dynamic programming approach.

    The layout of this paper is as follows: in section 2, an age-structure SIRS epidemic model is reproduced and new age-structure SIRS epidemic model with vaccination is given; in section 3, we give the basic regeneration number and discuss the threshold dynamics; in section 4, we formulate and solve the corresponding optimal control problem. The existence of optimal control is proved by HJB equation and the expression of optimal control is present. in section 5, some numerical simulations are performed to demonstrate the theoretical results. Brief conclusions are given in section 6.

    We reproduce the age-structured SIRS epidemic model based on paper [24] as follows

    $ {dS(a,t)=[S(a,t)a+λ(a)μ(a)S(a,t)β(a)S(a,t)I(a,t)+ρ(a)η(a)I(a,t)]dt,dI(a,t)=[I(a,t)a+β(a)S(a,t)I(a,t)δ(a)I(a,t))]dt,dR(a,t)=[R(a,t)a+(1ρ(a))η(a)I(a,t)(μ(a)+γ(a))R(a,t)]dt,S(a,0)=S0(a),I(a,0)=I0(a),R(a,0)=R0(a),
    $
    (2.1)

    where $ S(a, 0) = S_0(a), I(a, 0) = I_0(a), R(a, 0) = R_0(a) $ are initial conditions, with the following conditions:

    $ {S(0,t)=A0β(a)S(a,t)da,I(0,t)=k1A0β(a)I(a,t)da,R(0,t)=k2A0β(a)R(a,t)da,
    $
    (2.2)

    where $ k_1, \; k_2\in[0, 1] $ present the weight of the infected and recovered, respectively. All the parameters are related to age $ a $, positive and bounded, These parameters and their meaning are listed in the Table 1, where $ \eta (a) = \sigma (a)+\alpha(a) $, and $ \delta(a) = \mu(a)+\phi(a)+ \sigma (a)+\alpha(a) $.

    Table 1.  Parameter values used in numerical simulations.
    Notation Biological meanings
    $ S(a, t) $ Densities of susceptible individuals of age $ a $ at time $ t $
    $ I(a, t) $ Densities of infective individuals of age $ a $ at time $ t $
    $ R(a, t) $ Densities of recovery individuals of age $ a $ at time $ t $
    $ \lambda (a) $ Recruitment rate of age $ a $
    $ \mu (a) $ Mortality rate of age $ a $
    $ \beta (a) $ Transmission coefficient of age $ a $
    $ \gamma (a) $ Immunity loss rate of age $ a $
    $ \delta (a) $ Infected removal rate of age $ a $
    $ \eta (a) $ Infected recovery rate of age $ a $
    $ \rho (a) $ Switch between $ SIS \leftrightarrow SIR(S) $ of age $ a $
    $ \alpha (a) $ Treatment rate of age $ a $
    $ \sigma (a) $ Nature cure rate of age $ a $
    $ \phi (a) $ Infected-induced mortality rate of age $ a $

     | Show Table
    DownLoad: CSV

    As we all know that the age density distribution function of the total population reaches a stable state, i.e

    $ S(a,t)+I(a,t)+R(a,t)=P(a).
    $
    (2.3)

    We focus on threshold dynamics on an age-structured SIRS epidemic model with vaccination. The spread of many epidemics are age-related, such as Human Papilloma Virus (HPV) spread across the population after ten years olds usually. And the World Health Organization believes that the most appropriate age for HPV vaccination was 11-12 years old. Therefore, it makes sense to study the vaccination is only given for a certain age to age-structured SIRS epidemic model.

    The introduced vaccination rate $ u(A_0) $ satisfies $ 0\leq u(A_0)\leq 1 $ and is only given for a certain age and susceptible individuals of age $ A_0 $ is inoculated at some point. The disease does not spread among people younger than age $ A_0 $, that is, there are no patients before the vaccination age $ A_0 $. Therefore, we can get that the system satisfied by the susceptible, infective, and recovering individuals based on model (2.1) as following

    $ {S(a,t)=P(a),I(a,t)=0,R(a,t)=0,0<a<A0,dS=(Sa+λ(a,t)μ(a)Sϱ(a,t)S+ρ(a)η(a)I)dt,A0<a<A,dI=(Ia+ϱ(a,t)Sδ(a)I)dt,A0<a<A,dR=(Ra+(1ρ(a))η(a)Iμ(a)Rdt,A0<a<A,S(A0,t)=(1u(A0))S(A0,t),I(A0,t)=I(A0,t)=0,R(A0,t)=u(A0)S(A0,t),S(a,0)=S0(a),I(a,0)=I0(a),R(a.0)=R0(a),ϱ(a,t)=k(a)AA0I(ˆa,t)dˆa,
    $
    (2.4)

    we call infectivity

    $ ϱ(a,t)=β(a)I(a,t)=k(a)AA0I(ˆa,t)dˆa<,
    $

    where $ S(A_0^-, t) $ and $ I(A_0^-, t) $ are the densities of susceptible and recovering individuals in the left neighborhood of age $ A_0 $, $ k(a) $ is a positive continuous function. The vaccination rate $ u(A_0) $ is a constant and $ k(a) $ is a nonnegative continuous function defined on $ 0\leq a \leq A $, where $ A $ is the maximum age.

    In this paper, unless otherwise specified. The $ \mathbb{H} $ represents the Hilbert space, with inner product

    $ φ1,ϕ1=A0φ1(a)ϕ1(a)da,φ1,ϕ1H,
    $
    (2.5)

    and we consider the state equations on space $ L^1(0, A; \mathbb{H}^3) $, when $ \varphi = (\varphi_1, \varphi_2, \varphi_3) $ and $ \phi = (\phi_1, \phi_2, \phi_3) $, $ \varphi_i, \; \phi_i \in \mathbb{H} $, the definition of inner product on $ L^1(0, A; \mathbb{H}^3) $ as

    $ φ,ϕ=A0(φ1(a)ϕ1(a)+φ2(a)ϕ2(a)+φ3(a)ϕ3(a))da=A0φ(a)ϕ(a)da,φ,ϕH3.
    $
    (2.6)

    The existence of time-independent equilibrium solutions for the model (2.4) is discussed. Obviously, the model (2.4) has disease-free equilibrium

    $ S(a,t)=P(a),I(a,t)=0,R(a,t)=0,0<a<A0,S(a,t)=(1u(A0))P(a),I(a,t)=0,R(a,t)=u(A0)P(a),A0<a<A.
    $
    (3.1)

    In order to obtain the endemic disease equilibrium solution, $ S(a, t) = S(a), \; I(a, t) = I(a) $ and $ R(a, t) = R(a) $ are substituted into the model (2.4) obtained

    $ R(a)=0,0<a<A0,R(a)=u(A0)P(a),A0<a<A.
    $
    (3.2)

    And we know when $ R_0(a) = u(A_0)P(a) $, the $ R(a, t) = u(A_0)P(a) $, and when the population reaches a dynamic equilibrium $ P(a) $, Combining the equations (2.3) and (3.2), and substitute into (2.4) for the equation that satisfies $ I(a, t) $. From (2.3), it shows that

    $ I(a)a+I(a)t=δ(a)I(a)+k(a)(P(a)I(a)R(a))D(t),
    $
    (3.3)

    where $ D (t) = \int _{A_0}^A I(\hat{a}, t)d\hat{a} $. Let $ \zeta (a, D(t)) = \delta(a)+ k(a)D(t) $, in order to solve the (3.3), when $ A_0 < a < A $, the expression of $ I(a, t) $ is

    $ I(a,t)=I0(at)exp(t0ζ(at+θ,D(θ))dθ)+[1u(A0)]t0exp(tτζ(at+θ,D(θ))dθ)×k(at+θ)P(at+θ)D(τ)dτ,a>t+A0,I(a,t)=[1u(A0)]aA0exp(aτζ(θ,D(ta+θ))dθ)×k(τ)P(τ)D(ta+τ)dτ,at+A0.
    $
    (3.4)

    we can substitute (3.4) into the $ D(t) $ yields when $ A_0 < a < A $

    $ D(t)=[1u(A0)]AA0aA0k(τ)P(τ)D(tτ)×exp(a+τaζ(θ,D(taτ+θ))dθ)dadτ+aA0I0(a)exp(t0ζ(a+θ,D(θ))dθ)da.
    $
    (3.5)

    Since

    $ \lim\limits_{t\rightarrow \infty} \int _{A_0}^a I_0(a) exp{(-\int_0^t \zeta (a+\theta , D(\theta))d \theta )}da \leq \int _{A_0}^a I_0(a) exp{(-\int_0^t \mu (a+\theta , N) d\theta )}da = 0, $

    obviously, we have $ D(t) = D_\infty $ if $ I(a, t) = I_{\infty} $ as $ t\rightarrow \infty $ is equilibrium of model, then

    $ D=[1u(A0)]AA0aA0k(τ)P(τ)exp(a+τaζ(θ,D)dθ)Ddadτ+aA0I0(a)exp(t0ζ(a+θ,D)dθ)da,
    $
    (3.6)

    then

    $ D=D[1u(A0)]AA0aA0k(τ)P(τ)exp(a+τaδ(θ)dθ)×exp(Da+τak(θ)dθ)dτda.
    $
    (3.7)

    Then it expresses $ I_{\infty} = 0 $ or

    $ 1=[1u(A0)]AA0aA0k(τ)P(τ)exp(a+τaδ(θ)dθ)×exp(Da+τak(θ)dθ)dτda.
    $
    (3.8)

    It is quite obvious the right-hand side of the equation (3.8) is the minus function about $ D_{\infty} $, thus, we can define

    $ R0=(1u(A0))A0aA0k(τ)p(τ)exp(aτδ(θ)dθ)dτda.
    $
    (3.9)

    When $ R_0 < 1 $, (3.7) has only a unique zero solution and no positive root, when $ R_0 > 1 $, the (3.7) formula has zero root and the only positive root $ I^* > 0 $ satisfies (3.8) formula, thus there are only positive root $ E^* = (S^*, I^*, R^*) $ of (2.4). For ease of analysis, we converse $ i(a, t) = \frac{I(a, t)}{(1-u(A_0))P} $, equation (2.4) about $ I(a, t) $ can be transformed into

    $ {i(a,t)t+i(a,t)a=[1i(a,t)][1u(A0)]k(a)AA0P(ˆa)i(ˆa,t)dˆaδ(a)i(a,t),a>A0,i(A0,t)=0,t0,i(a,0)=i0(a)=I0(a)(1u(A0))P,aA0.
    $
    (3.10)

    From the transformation, we can see that the existence and stability of the positive equilibrium of model (2.4) is equivalent to positive equilibrium of (3.10) model. $ R_0 $ can be used as a threshold for the existence or absence of endemic diseases. It is also a threshold for the disappearance of diseases. Let $ \hbar = 1-\frac{1}{R_0} $ and we call $ \hbar $ vaccination elimination point.

    Similar to the [25] method, the following conclusion can be obtained. Before discussing the global stability of equilibriums, we first propose the lemma 3.1, and the similar proof process have been presented in [26].

    Lemma 3.1. Let $ i_1(a, t) $ and $ i_2(a, t) $ are the solutions of the $(3.10)$ equations, satisfying $ i_1(a, 0) = i_{10}(a) \leq i_2(a, 0) = i_{20}(a) $, where $ i_{10}(a), \; i_{20}(a)\in K $, then there are following properties

    (1) $ i_j(a, t)\in K, \; j = 1, \; 2; $

    (2) $ \; i_1(a, t) \leq i_2(a, t), \; t \geq 0 $;

    (3) When $ 0 < \xi < 1 $, the solution $ i(a, t) $ of $(3.10)$ that the initial value satisfies the $ i(a, 0) = i_{10}(a)\xi $, then it must satisfy the $ i_1(a, t)\xi \leq i(a, t) $,

    where the $ K $ is denoted

    $ K={f(a,t)|f(a,t)C,f(0,t)=0,0f(a,t)1}.
    $

    Theorem 3.2. When $ R_0 < 1 $, the model $(3.10)$ has only a disease-free equilibrium, which is globally asymptotically stable, while when $ R_0 > 1 $, the model $(3.10)$ also has a unique endemic equilibrium $ E^* $, which is globally asymptotically stable and the disease-free equilibrium is unstable.

    Proof. We only consider the asymptotic behavior of solutions where the initial values of the model (3.10) satisfy $ 0\leq i(a, 0) \leq 1 $ in $ (A_0, A) $. The theorem of the global stability is given for solutions.

    It is not difficult to see if the model (3.10) has a positive equilibrium equivalent to the (3.5) has a positive root. Since the right of the (3.5) is reduced about $ D(t) $ and tends to zero when the $ D(t) $ tends to positive infinity. There is no positive real root to the (3.5) when $ R_0 < 1 $, and there is a unique positive real root when $ R_0 > 1 $, i.e. the model (2.4) has only a unique positive equilibrium when $ R_0 > 1 $ and, instead, no positive equilibrium when $ R_0 < 1 $.

    When $ R_0 < 1 $, in order to study the global stability of the disease-free equilibrium, expression of the solution of (3.10) about $ I(a, t) $ is obtained by using the characteristic line method

    $ i(a,t)=i0(at)exp(t0[δ(at+τ)+[1u(A0)]k(at+τ)AA0P(ˆa)i(ˆa,τ)dˆa]dτ)+t0[1u(A0)]k(at+τ)AA0P(ˆa)i(ˆa,τ)dˆa×exp(tτ[δ(at+θ)+[1u(A0)]k(at+θ)AA0P(ˆa)i(ˆa,θ)dˆa]dθ)dτ,at+A0,
    $
    (3.11)
    $ i(a,t)=aA0[1u(A0)]k(τ)AA0P(ˆa)i(ˆa,ta+τ)dˆa×exp(aτ[δ(θ)+[1u(A0)]k(θ)AA0P(ˆa)i(ˆa,ta+θ)]dθ)dτ,a<t+A0.
    $
    (3.12)

    When $ t > A $, we obtain

    $ i(a,t)=aA0[1u(A0)]k(τ)χ(ta+τ)×exp(aτ[δ(θ)+[1u(A0)]k(θ)χ(ta+θ)]dθ)dτ,a<t+A0,
    $
    (3.13)

    where $ \chi(t) = \int_{A_0}^A P(\hat{a}) i(\hat{a}, t)d\hat{a} $, when $ t < A $, the $ \chi(t) $ satisfy is

    $ χ(t)=[1u(A0)]A0P(ˆa)ˆaA0k(τ)χ(tˆa+τ)×exp(ˆaτ[δ(θ)+[1u(A0)]k(θ)χ(tˆa+θ)]dθ)dτdˆa,a<t+A0.
    $
    (3.14)

    For natural numbers $ n > 1 $, let $ \chi_n = \max \limits_{nA\leq t \leq (n+1)A}\chi(t), \; L_n = \max \{\chi_{n-1}, \chi_n \} $. According to the (3.14) when $ R_0 < 1 $, we get

    $ χ(t)[1u(A0)]A0P(ˆa)ˆaA0χ(tˆa+τ)exp(ˆaτδ(θ)dθ)dτdˆaR0maxtAτtχ(τ).
    $
    (3.15)

    So we obtain $ \chi(t) \leq R_0 \max \{ \chi_{n-1}, \; \chi_n \} $ when $ nA \leq t \leq (n+1)A $. Because of $ R_0 < 1 $, we also obtain $ \chi_n < L_n = \chi_{n-1} $, i.e. $ \chi_n < \chi_{n-1} $, $ \chi_n $ is monotonous and reduced and $ \chi_n < \chi_{n-1}R_0 $, from this delivery $ \lim\limits_{n\rightarrow +\infty} \chi_n = 0 $. Considering above conditions and (3.13) can also derive $ \lim\limits_{t\rightarrow +\infty} i(a, t) = 0 $. Finally, it is concluded from (3.13) that the disease-free equilibrium of the model (3.10) is stable $ R_0 < 1 $.

    When $ R_0 > 1 $, substituting $ i(a, t) = i(a) $ into the (3.10), the only endemic equilibrium of (3.10) is given by

    $ i(a)=[1u(A0)]aA0k(a)AA0P(ˆa)i(ˆa)dˆa×exp(aτδ(θ)+[1u(A0)]k(τ)AA0P(ˆa)i(ˆa)dˆadθ)dτ.
    $
    (3.16)

    When $ A_0\leq a < A $, the $ \forall\; i(a, t) $ satisfy $ i(a, 0) > 0 $, and $ \chi(t) > 0 $ can be derived from $ (3.11) $ and $ (3.14) $ when $ 0 \leq t\leq A $, and we also know that $ i(a, A) > 0\; (A_0 < a < A) $ by (3.13). Repeat this process then the $ \chi(t) > 0 $ is obtained when $ 0 \leq t\leq 3A $, let $ \chi^0 = \min\limits_{A \leq t\leq 2A} \chi(t) $, $ \chi^* = \max\limits_{A \leq t\leq 2A}\chi(t) $, the following is got as

    $ i(a,2A)[1u(A0)]aA0χ0k(τ)×exp(aτ(δ(θ)+[1u(A0)]χk(θ))dθ)dτ.
    $
    (3.17)

    Comparing the forms of $ (3.16) $ and $ (3.17) $, we can find a positive number $ 0 < \xi < 1 $ such that $ i(a, 2A) \geq \xi i(a) $, namely, the following formula holds for the solution $ i(a, t) $ of $ (2.4) $

    $ ξi(a)i(a,2A)1.
    $

    The Lemma 3.1 shows that the solutions $ i_{\xi}(a, t) $ and $ i_1(a, t) $ with initial values of $ \xi \cdot i(a, 0) $ and $ 1 $ respectively when $ t = 0 $, must have

    $ iξ(a,t)i(a,t+2A)i1(a,t).
    $

    Based on the fact that $ i(a) $ is the equilibrium of $ (3.10) $, then the $ \xi i(a) \leq i_{\xi}(a, t) $, i.e. $ i_{\xi}(a, t) $ is a monotone increasing function about $ t $. Since $ i_1(a, t) $ is a monotone subtraction function about $ t $, and the corresponding $ \chi_{\xi}(t) $ and $ \chi_1(t) $ are monotone increasing and monotone decreasing, respectively, and there are $ \lim\limits_ {t\rightarrow \infty} \chi_{\xi}(t) = \chi_{\xi}^0 \leq \lim\limits_ {t\rightarrow \infty} \chi_{1}(t) = \chi_{1}^0 $. From (3.14), the equation has only unique normal number solutions when $ R_0 > 1 $, i.e. $ \chi_{\xi}^0 = \chi_{1}^0 = \xi^{*} $, the $ \chi (t) $ corresponding to the $ i(a, t) $ also obtain $ \lim\limits_ {t\rightarrow \infty} \chi (t) = \xi^{*} $. Reuse (3.13) when $ t > A_0 $, can be obtained $ \lim\limits_ {t\rightarrow \infty} i(a, t) = i(a) $, i.e. the local disease equilibrium of (3.10) is globally attractive. The stability of the local disease equilibrium can be obtained by linearization. Therefore, the local disease equilibrium is globally asymptotically stable.

    Remark 3.3. 3.2 shows that the epidemic disease will die out if $ R_0 < 1 $, while the disease will be prevalent if $ R_0 > 1 $. This implies that $ R_0 $ is the threshold of model $(2.4)$.

    In section 3, we consider the effect of a age determined vaccine immunization of the system on threshold dynamics, which is one of the cases in which vaccines are added to all susceptible persons in the system. In this section, we focus on the optimal control of (4.1) on the basis of age-structured SIRS by implementing proportion of vaccination $ u(a, t) $ as control strategy at all ages of system, and vaccination rate $ u(a, t) $ satisfies $ 0 \leq u(a, t) \leq 1 $. We aim to minimize the infected population, while keeping the cost of applying such control strategy at the minimum level.

    Let $ ((a, t_0), x(a, t_0))\in [0, T)\times(\mathbb{H}\times\mathbb{H}\times\mathbb{H}) $, and consider the following control system over on $ (a, t)\in (0, A)\times[t_0, T] $:

    $ {S(a,t)t+S(a,t)a=λ(a)μ(a)S(a,t)β(a)S(a,t)I(a,t)+ρ(a)η(a)I(a,t)+γ(a)R(a,t)u(a,t)S(a,t),I(a,t)t+I(a,t)a=β(a)S(a,t)I(a,t)δ(a)I(a,t)),R(a,t)t+R(a,t)a=(1ρ(a))η(a)I(a,t)(μ(a)+γ(a))R(a,t)+u(a,t)S(a,t),S(a,t0)=St0(a),I(a,t0)=It0(a),R(a,t0)=Rt0(a),
    $
    (4.1)

    where $ S(a, t_0) = S_{t_0}(a), I(a, t_0) = I_{t_0}(a), R(a, t_0) = R_{t_0}(a) $ is varying initial time value and denote the initial state $ x(a, t_0) $ as $ x_{t_0}(a) $. Here the control $ \{u(\cdot, \cdot):(0, A)\times[t_0, T]\rightarrow \mathcal {U}\} $, and $ u(\cdot, \cdot) $ is measurable. We obtain the schematic diagram of model (4.1) with introducing vaccination control (see Figure 1).

    Figure 1.  Schematic diagram of the model (4.1).

    The object is to design the optimal controller $ u^*(a, t), \; (a, t)\in (0, A)\times[t_0, T] $ for the system (4.1) which minimizes the following performance index during the finite time interval $ [t_0, T] $ [27,28,29,30]. Our goal is to minimize the total number of the infected and susceptible individuals by using minimal control efforts. The cost function as follows

    $ J(t0,xt0(a),u(a,t))=Tt0A0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))dadt+A0h(I(a,T))da,
    $
    (4.2)

    where $ I(a, t): (a, t)\in [0, A]\times[t_0, T] $ and we denote $ x = (S, I, R), x_{t_0} \in (\mathbb{H}\times\mathbb{H}\times\mathbb{H}) $, the terminal value is $ V(x, T) = \int_0^Ah(I(a, T))da $. Let $ L(t, x, u): = \tau_1I(a, t)+\tau_2u(a, t)+\frac{1}{2}\tau_3u^2 (a, t) $, where $ \tau_1 $, $ \tau_2 $ and $ \tau_3 $ are the positive weighting factors, representing the cost per unit time of the components $ I(a, t) $, $ u(a, t) $ and $ u^2(a, t) $, respectively. In particular, $ \int_{t_0}^T\int_0^A\tau_1I(a, t)dadt $ is the cost that infected individual creates for the society due to standard medical care, not including the vaccination treatment $ u(a, t) $. $ \int_{t_0}^T \int_0^A(\tau_2u(a, t)+\frac{1}{2}\tau_3u^2 (a, t))dadt $ is the cost of treating infected individuals. $ h:\mathbb{H}\times\mathbb{H}\times\mathbb{H} \rightarrow \mathbb{H} $ are measurable and $ h(I(a, T)) $ denotes the total number disease in population depend on the age at the time $ T $. The value function as

    $ {V(t0,xt0(a))=infu()UJ(t0,xt0(a),u(a,t)),t0[0,T),V(T,xt0(a))=A0h(I(a,T))da.
    $
    (4.3)

    The function $ V(\cdot, \cdot) $ will play an important role in obtaining optimal control, thus, we would like to study $ V(\cdot, \cdot) $ in great detail. Let $ (x^*, u^*) $ is optimal pair, where $ x^* = (S^*, I^*, R^*) $. Noting that the initial time $ (t = 0) $ and the initial state $ (x(a, 0) = x_0) $ are fixed in the model [18]. The $ h(I(a, T)) $ satisfies

    $ A0h(I(a,T))h(˜I(a,T))daCI(a,T)˜I(a,T)C,
    $
    (4.4)

    where $ C $ is a positive constant, $ \|\cdot\| $ represents norm on space $ \mathbb H $ and $ I(a, T) $, $ \tilde{I}(a, T)\in \mathbb H $. The following hypothesis are put forward and need to be satisfied.

    D1 $ (\mathcal{U}, d) $ is a separable metric space and $ T > 0 $.

    D2 The control set $ \mathcal{U} $ is bounded convex

    D3 Different controls correspond to same terminal value.

    Then we study dynamic programming method, to give optimal strategy of model (4.1). The basic idea of this approach applied to optimal controls is to consider a family of optimal control problems with different initial times $ (t = t_0\geq 0) $ and initial states $ (x(a, t_0) = x_{t_0}(a)) $ [see [31,32,33]], to establish relationships among these problems through the HJB equation. Firstly, we present the Bellman's principle of optimality as follows.

    Theorem 4.1. Let (D1)-(D3) hold. Then for any $ (t_0, x_{t_0}(a))\in [0, T)\times \mathbb{H}\times\mathbb{H}\times\mathbb{H} $,

    $ V(t0,xt0(a))=infu()U[t0,T]{ˆt0t0A0L(t,x(t;t0,xt0(a),u(,)),u(a,t))dt+V(ˆt0,x(ˆt0;t0,xt0(a),u(,)))},0t0ˆt0<T.
    $
    (4.5)

    where $ L(t, x(t; t_0, x_{t_0}(a), u(\cdot, \cdot)), u(a, t)) = \tau_1 I_{t_0}(a)+\tau_2u(a, t)+\frac{1}{2}\tau_3u^2 (a, t) $.

    Proof. Let us denote the right-hand side of (4.5) by $ \overline{V}(t_0, x_{t_0}(a)) $. By (4.3), we have

    $ V(t0,xt0(a))J(t0,xt0(a),u(a,t))=ˆt0t0A0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))dadt+J(ˆt0,x(a,ˆt0),u(,)),u(,)U.
    $

    Thus, taking the infimum over $ u(\cdot, \cdot) \in \mathcal {U} $ we get

    $ V(t0,xt0(a))¯V(t0,xt0(a)).
    $
    (4.6)

    Conversely, for any $ \varepsilon > 0 $, there exists a $ u_{\varepsilon}(\cdot, \cdot) \in \mathcal {U} $ such that

    $ V(t0,xt0(a))+εJ(t0,xt0(a),uε(,))ˆt0t0A0(τ1Iε(a,t)+τ2uε(a,t)+12τ3u2ε(a,t))dadt+V(ˆt0,xε(a,ˆt0))¯V(t0,xt0(a)).
    $
    (4.7)

    Combing (4.6) and (4.7), we obtain (4.5).

    The equation (4.5) is dynamic programming equation, which gives a relationship among globally and locally optimal via value function.

    We would like to further study the (4.5), trying to obtain an equation for $ V(t_0, x_{t_0}(a)) $ simpler form. The following results give a partial differential equation that a continuously differentiable value function ought to satisfy. We know the $ V(t_0, x_{t_0}(a)) $ is a continuously differentiable value function via 4.1 and basic theorem of calculus.

    Define the Lyapunov equation (LE)

    $ \langle V_x, (f(a, t, x)+g(a, t, x)u(a, t))\rangle+\int_{0}^{A}(\tau_1I(a, t)+\tau_2u(a, t)+\frac{1}{2}\tau_3u^2 (a, t))da = 0. $

    That is

    $ 0=VS,(Sa+λ(a)μ(a)Sβ(a)SI+ρ(a)η(a)I+γ(a)Ru(a,t)S)+VR,(Ra+(1ρ(a))η(a)I(μ(a)+γ(a))R+u(a,t)S)+VI,(Ia+β(a)SIδ(a)I)+A0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))da=A0VS(Sa+λ(a)μ(a)Sβ(a)SI+ρ(a)η(a)I+γ(a)Ru(a,t)S)da+A0VR(Ra+(1ρ(a))η(a)I(μ(a)+γ(a))R+u(a,t)S)da+A0VI(Ia+β(a)SIδ(a)I)da+A0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))da,
    $
    (4.8)

    where $ V_{x} = (\frac{\partial V}{\partial S}, \; \frac{\partial V}{\partial I}, \; \frac{\partial V}{\partial R})^{\top} $, $ V_{t} = \frac{\partial V}{\partial t} $. The equation LE is the HJB equation when $ u(a, t) $ is optimal control of (4.1).

    $ {0=A0VS(Sa+λ(a)μ(a)Sβ(a)SI+ρ(a)η(a)I+γ(a)Ru(a,t)S)da+A0VR(Ra+(1ρ(a))η(a)I(μ(a)+γ(a))R+u(a,t)S)da+A0VI(Ia+β(a)SIδ(a)I)da+A0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))da,Vt=T=h(x),xH×H×H.
    $
    (4.9)

    In this subpart, we write the HJB equation of model (4.1) as well as its viscosity solution is discussed as following.

    Theorem 4.2. Let (D1)-(D3) and $ V(t, x)\in C^{1}([0, T]\times \mathbb{H}\times\mathbb{H}\times\mathbb{H}) $. Then $ V(t, x) $ is a solution of the following terminal value problem of HJB equation:

    $ {Vt+infuUH(t,x,u,Vx)=0,t[0,T),V|t=T=h(x),(t,x,u)[0,T]×H×H×H×U.
    $
    (4.10)

    which can be written as $(4.9)$.

    Proof. Fix $ (t_0, x_{t_0}(a))\in [0, T)\times \mathbb{H}\times\mathbb{H}\times\mathbb{H} $ and $ u\in \mathcal{U} $. Let $ x(\cdot, \cdot) $ be the state trajectory corresponding to the control $ u(\cdot, \cdot)\in \mathcal{U}^{\omega}[s, T] $ with $ u(a, t)\equiv u $. On the one hand, we can write according to (4.5)

    $ 0{V(ˆt0,x(ˆt0,a))V(t0,xt0(a))}ˆt0t0A0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))dadt.
    $
    (4.11)

    Rearranging the equation (4.11) and dividing the both side of the equation by $ \hat{t}_0-t_0 $ yield

    $ 0{V(ˆt0,x(a,ˆt0))V(t0,xt0(a))}ˆt0t01ˆt0t0ˆt0t0A0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))dadt={V(ˆt0,x(a,ˆt0))V(ˆt0,xt0(a))+V(ˆt0,x(a,ˆt0))V(t0,x(a,ˆt0))}ˆt0t01ˆt0t0ˆt0t0A0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))dadt.
    $
    (4.12)

    By (4.5) with $ \hat{t}_0\downarrow t_0 $, we obtain

    $ 0Vt(t0,xt0(a))Vx(t0,xt0(a)),dxdtA0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))da=A0{VS(t0,xt0(a))(S(a,t)a+λ(a)μ(a)S(a,t)β(a)S(a,t)I(a,t)+ρ(a)η(a)I(a,t)+γ(a)R(a,t)u(a,t)S(a,t))VI(t0,xt0(a))(I(a,t)a+β(a)S(a,t)I(a,t)δ(a)I(a,t))VR(t0,xt0(a))(R(a,t)a+(1ρ(a))η(a)I(a,t)(μ(a)+γ(a))R(a,t)+u(a,t)S(a,t))}daA0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))da.
    $
    (4.13)

    This results in

    $ 0Vt(t,x(a,t))+infuUA0{VS(t,x(a,t))(S(a,t)aμ(a)S(a,t)β(a)S(a,t)I(a,t)+ρ(a)η(a)I(a,t)+λ(a)+γ(a)R(a,t))+VI(t,x(a,t))(I(a,t)a+β(a)S(a,t)I(a,t)δI(a,t))+VR(t,x(a,t))(R(a,t)a+(1ρ(a))η(a)I(a,t)(μ(a)+γ(a))R(a,t))}da+A0(τ1I(a,t)τS(VS(t,x(a,t))VR(t,x(a,t)))+12τ3(S(VS(t,x(a,t))VR(t,x(a,t))))2)da.
    $
    (4.14)

    On the other hand, for any $ \varepsilon > 0 $, $ \forall0\leq t_0\leq \hat{t_0} < T $ with $ \hat{t}_0-t_0 > 0 $ small enough, there exists a $ u(\cdot, \cdot)\equiv u_{\varepsilon, \hat{t}_0}(\cdot, \cdot)\in \mathcal{U} $ such that

    $ V(t0,xt0(a))+ε(ˆt0t0)V(ˆt0,x(ˆt0;t0,xt0(a),u(a,t)))+ˆt0t0A0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))dadt.
    $
    (4.15)

    Thus, it follows that as (noting $ V\in C^1([0, T]\times \mathbb{H}\times\mathbb{H}\times\mathbb{H}) $)

    $ ε{V(ˆt0,x(ˆt0))V(t0,x(t0))}ˆt0t01ˆt0t0ˆt0t0A0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))dadt1ˆt0t0ˆt0t0[Vt(t0,xt0(a))Vx(t0,xt0(a)),dxdtA0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))da]dt1ˆt0t0ˆt0t0[Vt(t0,xt0(a))Vx(t0,xt0(a)),dxdtinfuUA0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))da]dt1ˆt0t0ˆt0t0Vt(t0,xt0(a))+infuUH(t,x(a,t),u(a,t),Vx(t,x(a,t)))dt.
    $
    (4.16)

    Taking limit on (4.16) with respect to $ \hat{t}_0\downarrow t_0 $, we have

    $ 0Vt(t,x(a,t))+infuUA0{VS(t,x(a,t))(S(a,t)aμ(a)S(a,t)β(a)S(a,t)I(a,t)+ρ(a)η(a)I(a,t)+λ(a)+γ(a)R(a,t))+VI(t,x(a,t))(I(a,t)a+β(a)S(a,t)I(a,t)δI(a,t))+VR(t,x(a,t))(R(a,t)a+(1ρ(a))η(a)I(a,t)(μ(a)+γ(a))R(a,t))}da+A0(τ1I(a,t)τS(VS(t,x(a,t))VR(t,x(a,t)))+12τ3(S(VS(t,x(a,t))VR(t,x(a,t))))2)da.
    $
    (4.17)

    Combining (4.14) and (4.17), we obtain our results.

    Definition 4.3. [34] A function $ v \in C([0, T]\times \mathbb{H}\times\mathbb{H}\times\mathbb{H}) $ is called a $ viscosity \; solution $ of $(4.9)$ if

    $ v(T,x(a,T))h(x(a,T)),x(a,T)H×H×H,
    $
    (4.18)

    and for any $ \psi \in C^1([0, T]\times \mathbb{R}) $, whenever $ v-\psi $ attains a local maximum at $ (t, x)\in [0, T)\times \mathbb{H}\times\mathbb{H}\times\mathbb{H} $, we have

    $ ψt(t,x)+supuUH(t,x,u,ψx(t,x))0.
    $
    (4.19)

    The function $ v \in C([0, T]\times \mathbb{H}\times\mathbb{H}\times\mathbb{H}) $ is called a $ viscosity\; solution $ of $(4.9)$ if in $(4.18)-(4.19)$ the inequalities $ \leq $ are changed to $ \geq $ and $ local\; maximum $ is changed to $ local\; minimum $. In the case that $ v $ is both a viscosity subsolution and supersolution of $(4.9)$, it is called a viscosity solution of $(4.9)$.

    Theorem 4.4. Let (D1)-(D3) hold. Then the value function of $ V(t, x) $ satisfies

    $ |V(t0,xt0(a))V(ˆt0,xˆt0(a))|K{|t0ˆt0|},(t0,xt0(a)),(ˆt0,xˆt0(a))[0,T]×H×H×H,
    $
    (4.20)

    for some $ K > 0 $. Moreover, $ V(t, x(a, t)) $ is the only solution of $(4.9)$ in the class $ C([0, T]\times\mathbb{H}\times\mathbb{H}\times\mathbb{H}) $.

    Proof. Through 4.2 and $ t_0\leq \hat{t}_0 $ we obtain

    $ |V(t0,xt0(a))V(ˆt0,xˆt0(a))|infuU{ˆt0t0A0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))dadt}ˆt0t0A0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))dadt|t0ˆt0|maxt0tˆt0A0(τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t))daK|t0ˆt0|.
    $
    (4.21)

    Thus, $ V(t, x(a, t)) $ satisfies (4.21). We can know that $ V(t, x(a, t)) $ is unique viscosity of the HJB equation via the same method as [[34], Theorem 2.5].

    In this subsection, we study the existence of the optimal control pair in finite time for the system (4.1) and construct the Hamiltonian of the optimal control problem.

    Theorem 4.5. There exists an optimal control $ u^* \in \mathcal{U} $ and a corresponding optimal state $ S^*, \; I^*, \; R^* $ such that

    $ V(t0,xt0(a))=infu()UJ(t0,xt0(a),u(a,t)),(t0,xt0(a))[0,T)×H×H×H,
    $

    subject to the control system (4.1).

    Proof. To show the existence of the optimal control for the problem under consideration, we use the result in [35,36]. We notice that the state and control variables are non-negative, and the control set $ \mathcal{U} $, by definition, is closed and bounded. The optimal system is bounded, which ensures the compactness needed for the existence of an optimal control. Further, the integrand $ \int_{t_0}^T \int_0^A(\tau_1I(a, t)+\tau_2u(a, t)+\frac{1}{2}\tau_3u^2(a, t))dadt $ is convex on the control set $ \mathcal{U} $ due to the biquadratic and quadratic nature of control variable $ u(a, t) $, which represents vaccination proportion. In addition, there exists a constant $ \nu > 1 $ and positive numbers $ \kappa_1 $, $ \kappa_2 $ such that

    $ τ1I(a,t)+τ2u(a,t)+12τ3u2(a,t)κ1|u(a,t)|νκ2.
    $

    We complete the existence of the optimal control of state variables.

    Theorem 4.6. Let $ u^*(a, t) $ be optimal control variable, $ S^*(a, t) $, $ I^*(a, t) $ and $ R^*(a, t) $ be corresponding optimal state variable of model $(4.1)$. The corresponding optimal control is given as follows:

    $ u^*(a, t) = \min\{\max[0, -\frac{1}{\tau_3}S(a, t)^*(V_R(x, t)-V_S(x, t))+\tau], 1\} $

    where $ \frac{\tau_2}{\tau_3} = \tau $.

    Proof. We define the Hamiltonian function $ H(t, x, u(a, t), V_x), \; t\in[0, T], \; u(a, t)\in \mathcal{U}_{ad}, \; x, \; V_x(x, t) \in \mathbb{H}\times\mathbb{H}\times\mathbb{H} $ and by constructing the Hamiltonian function and HJB equation to characterize this optimal control $ u(a, t) $ as follows:

    $ H=(Sa+λ(a)μ(a)Sβ(a)SI+ρ(a)η(a)I+γ(a)Ru(a,t)S)VS(x,t)+(Ia+β(a)SIδ(a)I)VI(x,t)τ1Iτ2u(a,t)12τ3u2(a,t)+(Ra+(1ρ(a))η(a)I(μ(a)+γ(a))R+u(a,t)S)VR(x,t).
    $
    (4.22)

    Let $ u^*(a, t) $ be given optimal control and $ S^*(a, t) $, $ I^*(a, t) $ and $ R^*(a, t) $ be corresponding optimal state variable of model (4.1). By the Hamiltonian function (4.22), optimal control $ u(a, t) = u^*(a, t) $ is obtained

    $ u(a,t)=1τ3S(a,t)(VR(x,t)VS(x,t))+τ.
    $
    (4.23)

    From the properties of the control set with the findings, then there is exists optimal control

    $ u(a,t)={0,1τ3S(VRVS)+τ<0,1τ3S(VR(x,t)VS(x,t))+τ,01τ3S(VRVS)+τ<1,1,11τ3S(VRVS)+τ.
    $
    (4.24)

    Remark 4.7. If we could obtain the value function $ V(t, x(a, t)) $ by solving the HJB equation, then we could, at least formally, construct an optimal pair for each of optimal control problem for different initial values. The principle involves the following steps to solve optimal control problem

    $ Step.1. $ Solve the HJB equation $(4.10)$ to find the $ V(T, x(a, t)) $;

    $ Step.2. $ Find $ u^*(a, t) $ through Hamiltonian function $(4.22)$;

    $ Step.3. $ Combine the $ u^*(a, t) $ to solve model $(4.1)$ to get the optimal pair $ (x^*(a, t), u^*(a, t)) $.

    This section aims to illustrate the effectiveness of our theoretical results obtained in previous sections.

    In the numerical simulation, the parameters are given in Table 2. We assume that $ (a, t) \in [0, 80]\times [0,500], \; \Delta t = 0.1, \; \Delta a = 1 $, take $ p_{\infty} = 1.0, \; \varrho = 0.2, \; \delta = 0.1827, \; k = 0.2 $, then $ R_0\approx 2.69 > 1 $. The parameters that we partially need to assume are reasonable according to the study of Yang et al. [14] has shown that the range of the basic reproduction number of the pandemic influenza in $ 1.0-4.0 $. We see from 3.2 the $ E^* $ is globally asymptotically stable. Actually, as showing in the following Figure 2 (a), the density of the infected individual $ I(a, t) $ tends to be a positive constant over time.

    Table 2.  Parameter values used in numerical simulations.
    Parameter Value Source of data
    $ \lambda $ 0.20548 [16], Lan et al.
    $ \beta $ $ 1.6805\times 10^{-4} $ [16], Lan et al.
    $ \mu $ $ 3.4246\times 10^{-5} $ [18], Mu $ \& $ Zhang
    $ \eta $ 0.02 [18], Mu $ \& $ Zhang
    $ \rho $ 0.05 [18], Mu $ \& $ Zhang
    $ \gamma $ 0.027 [16], Lan et al.
    $ \delta $ 0.2227 Assumed
    $ \varrho $ 0.2 Assumed
    $ u $ 0.35 Assumed
    $ A $ 80 years [16], Lan et al.
    $ T $ 500 days Assumed

     | Show Table
    DownLoad: CSV
    Figure 2.  The evolution of $ I(a, t) $ of (2.4) for $ R_{0}\approx 2.69 > 1 $ (a) and $ R_0\approx 0.78 < 1 $ (b).

    Take $ \varrho = 0.1, \; \delta = 0.1827, \; k = 0.1 $, then $ R_0\approx 0.78 < 1 $, thus we can know the $ E^0 $ is globally asymptotically stable, as showing in the following Figure 2 (b), the density of the infected individual $ I(a, t) $ tends to zero over time. These are the same conclusions as 3.2.

    We choose $ \delta = 0.1827 $, $ 0.2027 $ and $ 0.2227 $ to study the effect on the infected population and all the other parameters are fixed as in Table 2. We calculate the corresponding basic reproduction number $ R_0 \approx 0.9427 $, $ 0.8979 $ and $ 0.7802 $ when $ R_0 < 1 $, respectively. In Figure 3 (a), we can see that the extinction time of decreases with the increase of infected removal rate $ \delta $. In order to simulate the cases of $ R_0 > 1 $, we calculate the corresponding basic reproduction number $ R_0 \approx 1.3468 $, $ 2.2441 $ and $ 2.6930 $, respectively. The simulation result is presented in Figure 3 (b). Synthetically, it is obvious to see from Figure 3 that the peak of $ I(a, t) $ decreases significantly as infected removal rate $ \delta $ of age a goes up.

    Figure 3.  (a): the path of $ I(a, t) $ under different $ \delta $ when $ R_0 < 1 $ at age $ a = 30 $; (b): the path of $ I(a, t) $ under different $ \delta $ when $ R_0 > 1 $ at age $ a = 30 $.

    In this section, we show the result of optimal control using a numerical example. The parameters are given in Table 2. $ (a, t) \in [0, 80]\times [0,500], \; \Delta t = 0.1, \; \Delta a = 1 $ by performing some numerical simulations. We take $ \tau_1 = 11.7 $, $ \tau_2 = 3.4 $ and $ \tau_3 = 12 $. The corresponding paths of control and infective populations $ I(a, t) $ are plotted in Figure 4, Figure 5 (a), Figure 6 and Figure 7 (a) show the section view at age $ a = 20 $ and $ a = 30 $ of the image variation tendency of $ I(a, t) $ with and without control and the tendency of corresponding optimal control $ u(a, t) $. We observe that the comprehensive use of control $ u(a, t) $ works better than no control applied during the course of disease. The effect of applying control reduces the number of infective populations $ I(a, t) $ dramatically. In short, control variable $ u(a, t) $ play an important role in the control of the disease. At the same time, we can also observe a reduced epidemic time of population infection.

    Figure 4.  When $ R_0>1 $, (a): the paths of $ I(a, t) $ with and without control $ u(a, t) $ at age $ a = 20 $; (b): The paths of $ I(a, t) $ with and without control $ u(a, t) $ at age $ a = 30 $.
    Figure 5.  (a): when $ R_0>1 $, the paths of control at age $ a = 20 $ and $ a = 30 $; (b): the path of control at $ t = 500 $.
    Figure 6.  When $ R_0<1 $, (a): the paths of $ I(a, t) $ with and without control $ u(a, t) $ at age $ a = 20 $; (b): The paths of $ I(a, t) $ with and without control $ u(a, t) $ at age $ a = 30 $.
    Figure 7.  (a): when $ R_0 < 1 $, the paths of control at age $ a = 20 $ and $ a = 30 $; (b): the paths of control at $ t = 500 $.

    We note that a significant decline in the number of infected persons $ I(a, t) $ during vaccine control of susceptible at age $ a = 20 $ and $ a = 30 $ from these following figures.

    When $ R_0 > 1 $, $ I(a, t) $ trends is regard to the optimal controls, where the upper-level image is the evolution of $ I(a, t) $ without control and the lower-level image is the evolution of $ I(a, t) $ with control on in Figure 4. In Figure 5 (a) shows the optimal control pathes at ages $ a = 20 $ and $ a = 30 $, repectively. Eventually the disease tends to extinction with control. Compared with the path of non-vaccinated $ I(a, t) $, the epidemic time of $ I(a, t) $ with vaccination is significantly shortened, as well as the density of infected population verge to zero eventually.

    When $ R_0 < 1 $, $ I(a, t) $ trends is regard to the optimal controls, where the upper-level image is the evolution of $ I(a, t) $ without control and the lower-level image is the evolution of $ I(a, t) $ with control in Figure 6. In Figure 7 (a) shows the optimal control pathes at ages $ a = 20 $ and $ a = 30 $, repectively. Eventually the disease tends to extinction with control. Compared with Compared with the path of non-vaccinated $ I(a, t) $, the epidemic time of $ I(a, t) $ with vaccination is significantly shortened. Combined with Figure 4 and the corresponding conclusions, it is clear that vaccination have a certain degree of hindrance to the transmission of epidemics. We also observe the Figure 5 (b) and Figure 7 (b) that vaccine control varies for different ages, mainly adolescents and middle age, and this control strategy is consistent with the actual control of many epidemics such as pneumonia, influenza and so on. The pneumococcal vaccine, COVID-19 vaccine and influenza vaccine are not suitable for the elderly. Overall, population age has a significant effect on the implementation of vaccine control strategies for many common diseases.

    In this paper, we recurred an age-structured SIRS epidemic model (2.1) with vaccination control, and obtained the condition for disease extinction and persistence. The research results showed that when $ R_0 < 0 $, the disease extincted; when $ R_0 > 0 $, the disease persisted. The conditions of disease extinction and permanence were given in 3.2. And we verified the resluts by numerical simulation. However, since the bifurcation phenomenon will occur when $ R_0 = 0 $, this part will also be an important point for our further study of the age-structured SIRS epidemic system dynamics. Besides, we also studied the control problem of SIRS epidemic model with age structure. By proving the existence of viscosity solution, we obtain the existence of optimal control. Then, we obtained the implicit expression for optimal control and HJB equation. Finally, by using Milstein method, the optimal proportion of vaccination is obtained. In addition, we plan that in the next step, traveling wave solutions of an age-structured SIRS epidemic model will become one of our main research contents. Wu et al. [37] studied the existence and non-existence of traveling wave solutions for a diffusive age-structured SIR epidemic model. Their work provides some insights on how to deal with the high dimensional age-structured epidemic models.

    The authors thank the editor and referees for their careful reading and valuable comments. The research was supported in part by the Ningxia Natural Science Foundation (No. 2020AAC03067).

    The authors have no conflict of interest.

    [1] Smith CS (1960) A History of Metallography, University of Chicago Press.
    [2] Fielding L (2013) The bainite controversy. Mater Sci Techn 29:383–399. doi: 10.1179/1743284712Y.0000000157
    [3] van Bohemen S, Siemtsma J (2008) Modeling of isothermal bainite formation based on nucleation kinetics. Int J Mat Res 99:739–747. doi: 10.3139/146.101695
    [4] Gaude-Fugarolas D, Jacques PJ (2006) A new physical model for the kinetics of the bainite transformation. ISIJ Inter 46:712–717. doi: 10.2355/isijinternational.46.712
    [5] Rees GI, Bhadeshia H (1992) Bainite transformation kinetics Part 1: Modified model. Mater Sci Techn 8:985–993. doi: 10.1179/mst.1992.8.11.985
    [6] Tzeng TC (2000) Autocatalysis in bainite transformations. Mater Sci Eng A 293:185–190. doi: 10.1016/S0921-5093(00)01221-1
    [7] Zolotoresvky N, Nesterova E, Titovets Y, et al. (2013) Modeling the effect of austenite deformation on the bainite structure parameters in low carbon microalloyed steels. Int J Mat Res 104:337–343. doi: 10.3139/146.110872
    [8] Hunkel M, Lübben T, Hoffmann F, et al. (1999) Modellierung der bainitischen und perlitischen Umwandlung bei Stählen. HTM Härterei-Techn Mitt 54:365–373.
    [9] Quidort D, Brechet YJM (2002) A model of isothermal and non isothermal transformation kinetics of bainite in 0.5% C steels. ISIJ Inter 42:1010–1017.
    [10] Maier H-J, Ahrens U (2002) Isothermal bainitic transformation on low alloy steels: factors limiting prediction of the resulting material’s properties. Z Metallk 93:712–718. doi: 10.3139/146.020712
    [11] Freiwillig R, Kudrman J, Chraska P (1976) Bainite transformation in deformed austenite. Metall Trans A 7:1091–1097. doi: 10.1007/BF02656591
    [12] Holzweissig M, Canadinc D, Maier H-J (2012) In-situ characterization of transformation plasticity during an isothermal austenite-to-bainite phase transformation. Mater Char 65:100–108. doi: 10.1016/j.matchar.2012.01.007
    [13] Su T, Veaus M, Aeby-Gautier E, et al. (2003) Effect of tensile stresses on bainitic isothermal transformation. J Phys IV France 112:293–296. doi: 10.1051/jp4:2003886
    [14] Su T, Aeby-Gautier E, Denis S (2006) Morphology changes in bainite formed under stress. Scripta Mater 54:2185–2189. doi: 10.1016/j.scriptamat.2006.02.031
    [15] Hase K, Garcia-Mateo C, Bhadeshia H (2004) Bainite Formation influenced by large stress. Mater Sci Techn 20:1499–1505. doi: 10.1179/026708304X6130
    [16] Kundu S, Hase K, Bhadeshia S (2007) Crystallographic texture of stress affected bainite. Proc Royal Soc A 463:2309–2328. doi: 10.1098/rspa.2007.1881
    [17] Fuijiwara K, Okaguchi S, Ohtani H (1995) Effect of hot deformation on bainite structure in low carbon steels. ISIJ Inter 15:1006–1012.
    [18] Min J, Lin J, Min Y, et al. (2012) On the ferrite and bainite transformation in isothermally deformed 22MnB5 steels. Mater Sci Eng A 550:375–387. doi: 10.1016/j.msea.2012.04.091
    [19] Nikravesh M, Naderi M, Akbari GH (2012) Influence of hot plastic deformation and cooling rate on martensite and bainite start temperatures in 22MnB5 steel. Mater Sci Eng A A 540:24–29. doi: 10.1016/j.msea.2012.01.018
    [20] Karbasian H, Tekkaya AE (2010) A review on hot stamping. J Mat Proc Tech 210:2103. doi: 10.1016/j.jmatprotec.2010.07.019
    [21] Feuser P, Schweiker T, Merklein M (2011) Partially hot-formed parts from 22MnB5 - process window. ICTP Aachen 10:408.
    [22] Chen L (2002) Phase fields models for microstructure evolution. Annu Rev Mater Res 32:113. doi: 10.1146/annurev.matsci.32.112001.132041
    [23] Wang Y, Jin YM, Khachaturyan AG (2004) The effects of free surfaces on martensite microstructures: 3D phase field microelasticity simulation study. Acta Mat 52:1039. doi: 10.1016/j.actamat.2003.10.037
    [24] Micress microstructure evolution simulation software, www.micress.de.
    [25] Loginova I, Agren J, Amberg G (2004) On the formation ofWidmanstätten ferrite in a binary Fe-C phase-field approach. Acta Materialia 52:4055–4063. doi: 10.1016/j.actamat.2004.05.033
    [26] Song W, Prahl U, Bleck W, et al. (2011) Phase field simulations of bainitic phase transformation in 100Cr6. Supplemental proceedings: Materials Fabrication, Properties, Characterization, and Modeling 2:417–425.
    [27] Song W (2014) Characterization and simulation of bainite transformation in high carbon bearing steel 100Cr6, PhD thesis RWTH Aachen University.
    [28] Arif T, Qin R (2013) A phase field model for bainitic transformation. Computational Materials Science 77:230–235. doi: 10.1016/j.commatsci.2013.04.044
    [29] Arif T, Qin R (2014) A phase field model for the formation of martensite and bainite. Advanced materials research 922:31–36. doi: 10.4028/www.scientific.net/AMR.922.31
    [30] Qin R, Bhadeshia H (2009) Phase field model to study the effect of interface anisotropy on the crystal morphological evolution of cubic metals. Acta Materialia 57:2210–2216. doi: 10.1016/j.actamat.2009.01.024
    [31] Bhadeshia H (1987)Worked examples in the geometry of crystals, Institute of metals, London and Brookfield.
    [32] Bouville M, Ahluwalia R (2006) Interplay between diffusive and displacive phase transformations: Time-Temperature-Transformation diagrams and microstructures. Phys Rev Lett 97:055701. doi: 10.1103/PhysRevLett.97.055701
    [33] Kundin J, Raabe D, Emmerich H (2011) A phase field model for incoherent martensitic transformations including plastic accommodation processes in the austenite. Journal of the mechanics and physics of solids 59:2082–2102. doi: 10.1016/j.jmps.2011.07.001
    [34] Kundin J, Pogorelov E, Emmerich H (2015) Numerical investigation of the interaction between the martensitic transformation front and the plastic strain in austenite. Journal of the mechanics and physics of solids 76:65–83. doi: 10.1016/j.jmps.2014.12.007
    [35] Levitas V, Javanbakht M (2013) Phase field approach to interaction of phase transformation and dislocation evolution. Applied Physics Letters 102:251904. doi: 10.1063/1.4812488
    [36] Roters F, Eisenlohr P, Hantcherli L, et al. (2010) Overview of constitutive laws, kinematics, homogenization and multiscale methods in crystal plasticity finite element modeling: Theory, experiments, applications. Acta Materialia 58:2210–2216.
    [37] Johnson W, Mehl R (1939) Reaction kinetics in process of nucleation and growth. Trans AIME 135:416–458.
    [38] Avrami M (1941) Kinetics of phase change III: granulation, phase change and microstructure. J Chem Phys 9:177–184. doi: 10.1063/1.1750872
    [39] Lee G, Kim S, Han H (2009) Finite element investigations for the role of transformation plasticity on springback in hot press forming process. Comp Mater Sci 47:556567.
    [40] HunkelM(2012) Anisotropic transformation strain and transformation plasticity: two corresponding effects. Mat -wiss u Werkstofftech. 43:150–157.
    [41] Lütjens J, Hunkel M (2013) The influence of the transformation plasticity effect on the simulation of partial press-hardening. Proc 4th Int Conf CHS2 319–327.
    [42] Brener EA, Marchenko VI, Spatschek R (2007) Influence of strain on the kinetics of phase transitions in solids. Phys Rev E 75:041604. doi: 10.1103/PhysRevE.75.041604
    [43] Fratzl P, Penrose O, Lebowitz JL (1999) Modelling of Phase Separation in Alloys with Coherent Elastic Misfit. J Stat Phys 95:1429. doi: 10.1023/A:1004587425006
    [44] Freund L (1998) Dynamic fracture mechanics, Cambridge University Press.
    [45] Spatschek R, Brener E, Karma A (2011) Phase field modeling of crack propagation. Phil Mag 91:75. doi: 10.1080/14786431003773015
    [46] Chien FR, Clifton RJ, Nutt SR (1995) Stress-induced phase transformation in single crystal titanium carbide. J Am Ceram Soc 78:1537. doi: 10.1111/j.1151-2916.1995.tb08849.x
    [47] Spatschek R, Müller-Gugenberger C, Brener E A, et al. (2007) Phase field modelling of fracture and stress-induced phase transitions. Phys Rev E 75:066111. doi: 10.1103/PhysRevE.75.066111
    [48] Spatschek R, Eidel B (2013) Driving forces for interface kinetics and phase field models. Int J Solid and Structures 50:2424. doi: 10.1016/j.ijsolstr.2013.03.016
    [49] Steinbach I (2011) Phase field models in materials science. Modelling and Simulation in Materials Science and Engineering 17:073001.
    [50] Steinbach I, Shchyglo O (2011) Phase field modelling of microstructure evolution in solids: Perspectives and challenges. Current opinion in solid state and materials science 15:87. doi: 10.1016/j.cossms.2011.01.001
    [51] Rao M, Sengupta S (2003) Nucleation of solids in solids: ferrite and martensite. Phys Rev Lett 91:045502. doi: 10.1103/PhysRevLett.91.045502
    [52] Brener EA, Iiordanskii SV, Marchenko VI (1999) Elastic effects on the kinetics of a phase transition. Phys Rev Lett 82:1506. doi: 10.1103/PhysRevLett.82.1506
    [53] Brener EA, Boussinot G, Hüter C, et al. (2009) Pattern formation during diffusional transformations in the presence of triple junctions and elastic effects. J Phys Cond Mat 21:464106. doi: 10.1088/0953-8984/21/46/464106
    [54] Pilipenko D, Brener EA, Hüter C (2008) Theory of dendritic growth in the presence of lattice strain. Phys Rev E 78:060603. doi: 10.1103/PhysRevE.78.060603
    [55] Ivantsov GP (1947) PhD thesis Akad. Nauk. SSSR.
    [56] Steinbach I, Pezzolla F (1999) A generalized field method for multiphase transformations using interface fields. Physica D: Nonlinear Phenomena 134:385–393. doi: 10.1016/S0167-2789(99)00129-3
    [57] Song W, von Appen J, Choi P, et al. (2013) Atomic-scale investigation of epsilon and theta precipitates in bainite in 100Cr6 bearing steel by atom probe tomography and ab initio calculations. Acta Materialia 61(20):7582–7590.
    [58] Eiken J, Boettger B, Steinbach I (2006) Multi-phase-field approach for multi-component alloys with extrapolation scheme for numerical application. Phys Rev E 73:066122. doi: 10.1103/PhysRevE.73.066122
    [59] Steinbach I, Pezzolla F, Prieler R (1995) Grain selection in faceted crystal growth using the phase field theory. In: The 7th conference on on modeling of casting, welding and advanced solidification processes.
    [60] Lukas H, Fries S, Sundman B (2007) Computational thermodynamics: The CALPHAD method, Cambridge University Press.
    [61] Rees GI, Shipway PH (1997) Modelling transformation plasticity during the growth of bainite under stress. Materials Science and Engineering A 223:168–178. doi: 10.1016/S0921-5093(96)10478-0
    [62] Wolff M, Böhm M, Dalgic M, et al. (2008) Evaluation of models for TRIP and stress-dependent transformation behavior for the martensitic transformation of the steel 100Cr6. Comput Mater Sci 43:108114.
    [63] Denis S (1997) Considering stress-phase transformation interaction in the calculation of heat treatment residual stresses. Series: International Centre for Mechanical Sciences 368:293–317.
    [64] Leblond J, Deveaux JC (1989) Mathematical modelling of transformation plasticity in steels I: Case of ideal-plastic phases. Int J Plasticity 5:551–572. doi: 10.1016/0749-6419(89)90001-6
    [65] Fisher FD, Sun QP, Tanaka K (1996) Transformation-induced plasticity. Appl Mech Rev 49:317364.
    [66] Leblond JB, Deveaux J (1984) A new kinetic model for anisothermal metallurgical transformation in steels including effect of austenite grain size. Acta Metallurgica 32:137146.
    [67] ASTM International standard test methods for tension testing of metallic materials (2011). Available from: www.astm.org.
    [68] HunkelM(2009) Anisotropic transformation strain and transformation plasticity: two corresponding effects. Mat -wiss u Werkstofftech. 40(5-6):466–472.
    [69] Devaux J, Leblond JB, Bergheau JM (2000) Numerical study of the plastic behaviour of a low alloy steel during phase transformation. Journal of Shanghai Jiaotong University 3:206–212.
    [70] Zwigl P, Dunand DC (1997) A non-linear model for internal stress superplasticity. Acta Materialia 45(12):5285–5294.
    [71] Schicchi DS, Hunkel M (2015) Transformation plasticity effect during bainite transformation on a 22MnB5 Steel Grade. IDE 2015, Bremen, Germany.
  • Reader Comments
  • © 2015 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(8094) PDF downloads(1501) Cited by(8)

Figures and Tables

Figures(13)

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog