Applied Mathematical Modelling 37 (2013) 3777–3788
Contents lists available at SciVerse ScienceDirect
Applied Mathematical Modelling journal homepage: www.elsevier.com/locate/apm
Analysis and approximation of a reliable model Houbao Xu a,⇑, Weiwei Hu b a b
Department of Mathematics, Beijing Institute of Technology, Haidian District, Beijing 100081, China ICAM, Department of Math., Virginia Polytechnic Institute and State University, Blacksburg, VA 24061, United States
a r t i c l e
i n f o
Article history: Received 19 November 2011 Received in revised form 4 June 2012 Accepted 6 July 2012 Available online 8 September 2012 Keywords: Trotter–Kato theorem Convergence Dynamical solution
a b s t r a c t In this paper we consider the problem of approximating the dynamical system that models reliability of a system consisting of two machines separated by a finite storage buffer. The system is described as a distributed parameter system defined by a coupled partial and ordinary differential equations and formulated as an abstract Cauchy problem. To derive the dynamical solution and some instantaneous indexes of the model, we present a simple finite difference scheme and establish the convergence of this scheme by employing Trotter–Kato Theorem. Numerical results are given to illustrate the effectiveness of the scheme. Ó 2012 Elsevier Inc. All rights reserved.
1. Introduction Consider two machines separated by a storage buffer as illustrated in Fig. 1. Both machines have to process certain objects, but have different service rates. It is assumed that the pieces to be processed first arrive at machine 1 which process the pieces one by one. When a piece is finished by machine 1, it either passes from machine 1 to the buffer with probability g1 ð0 < g1 < 1Þ, or it is reprocessed by machine 1 immediately with probability q1 ðq1 ¼ 1 g1 Þ. Machine 2 takes the pieces from the buffer and also processes them one by one. When a piece is finished by machine 2, it either exits the system with probability g2 ð0 < g2 < 1Þ or it is reprocessed by machine 2 immediately with probability q2 ðq2 ¼ 1 g2 Þ. The buffer can store at most N pieces. If the buffer is full, then machine 1 rests until a piece leaves the buffer. This scenario is typical in computer integrated manufacturing system (CIMS). The formulation of the CIMS was proposed by Shu [1]. Reference [2] developed a mathematical model consisting of ordinary and partial differential equations, and studied the steady-state solution of the model. The time parameter t 2 ½0; 1Þ refers to the time evolution of the whole system, whereas x 2 ½0; 1Þ gives the elapsed service time of the current service. Let p0 ðtÞ denote the probability that at time t machine 1 is processing a piece and there are no pieces in the buffer and there is no piece in machine 2. Let pm ðx; tÞ; 1 6 m 6 N þ 1, denote the state that both machines are busy and there are m 1 pieces in the buffer. Here x is the elapsed service time of the piece which is currently processed by machine 2. The probability denoted by pm ðtÞ for this state at time t is given by
pm ðtÞ ¼
Z 0
1
pm ðx; tÞdx:
Finally, pNþ2 ðx; tÞ corresponds to the situation that machine 2 is processing a piece, and the elapsed service time for this piece is x. Meanwhile, the buffer is full, one piece is waiting to be stored in the buffer, and machine 1 is blocked. The probability denoted by pNþ2 ðtÞ for this state at time t is given by ⇑ Corresponding author. Tel.: +86 13520096549; fax: +86 1068948372. E-mail address:
[email protected] (H. Xu). 0307-904X/$ - see front matter Ó 2012 Elsevier Inc. All rights reserved. http://dx.doi.org/10.1016/j.apm.2012.07.056
3778
H. Xu, W. Hu / Applied Mathematical Modelling 37 (2013) 3777–3788
Fig. 1. Behavior of the system.
pNþ2 ðtÞ ¼
Z
1 0
pNþ2 ðx; tÞdx:
In reference [2], a system consisting of a storage buffer and two machines is formulated as the mathematical model defined by
dp0 ðtÞ ¼ kg1 p0 ðtÞ þ g2 dt
Z
1
lðxÞp1 ðx; tÞdx;
ð1Þ
@p1 ðx; tÞ @p1 ðx; tÞ þ ¼ ðkg1 þ lðxÞÞp1 ðx; tÞ; @t @x
ð2Þ
0
@pm ðx; tÞ @pm ðx; tÞ þ ¼ ðkg1 þ lðxÞÞpm ðx; tÞ þ kg1 pm1 ðx; tÞ; @t @x
2 6 m 6 N þ 1;
@pNþ2 ðx; tÞ @pNþ2 ðx; tÞ þ ¼ lðxÞpNþ2 ðx; tÞ þ kg1 pNþ1 ðx; tÞ; @t @x
ð3Þ
ð4Þ
with boundary conditions,
p1 ð0; tÞ ¼ kg1 p0 ðtÞ þ q2
pm ð0; tÞ ¼ q2
Z
0
1
lðxÞp1 ðx; tÞdx þ g2
1
lðxÞpm ðx; tÞdx þ g2
0
pNþ2 ð0; tÞ ¼ q2
Z
Z 0
Z 0
Z
1
lðxÞp2 ðx; tÞdx;
ð5Þ
lðxÞpmþ1 ðx; tÞdx; 1 6 m 6 N þ 1;
ð6Þ
0
1
1
lðxÞpNþ2 ðx; tÞdx;
ð7Þ
and initial conditions,
p0 ð0Þ ¼ 1;
pm ðx; 0Þ ¼ 0;
1 6 m 6 N þ 2:
ð8Þ
Here, k represents the precessing rate of machine 1, and lðxÞ represents the hazard rate of machine 2’s service time. Reference [3] converted the model (1)–(8) into an abstract Cauchy problem in the Banach space B,
(
_ PðtÞ ¼ APðtÞ; Pð0Þ ¼ ð1; 0; . . . ; 0Þ;
ð9Þ
P B ¼ R ðL1 ½0; 1ÞÞNþ2 , with k kB ¼ j j þ Nþ2 i¼1 k kL1 . Reference [3] also proved the existence of a unique positive solution by employing C 0 -semigroup theory. Furthermore, the C 0 -semigroup SðtÞ generated by A is contractive and the system is conservative, i.e. kPðtÞkB ¼ kSðtÞPð0ÞkB ¼ 1. In reference [4], it was proved that 0 is a simple eigenvalue of the system operator and also the unique spectral point on the imaginary axis. Moreover, using the approach in [5] one can prove that there is a continuous spectrum in a left-half part rc ðAÞ ¼ fr 2 CjRer 6 lg, and there exist at most a finite number of isolated point spectrum with finite algebraic multiplicity in the strip fr 2 Cj l < Rer 6 0g. All other points in this strip belong to the resolvent set qðAÞ. Again, using semigroup theory, [6] proved the system (1)–(8) is asymptotically stable. Much work has been done on the qualitative analysis of this system, however, until now the issue of numerical approximation has not been studied. Although there are considerable literatures on functional differential equations in Lp space, most of the computational work is restricted to L2 space (see [7,8] and etc.). Here, L1 is the natural space for this problem and the computational challenge is more difficult since L1 is not reflexive. The remainder of this paper is organized as follows. Section 2 reviews the Trotter–Kato Theorem. Section 3 presents an approximation scheme and its convergence is established by applying Trotter–Kato Theorem. Numerical results are given in Section 4. Section 5 concludes the paper.
H. Xu, W. Hu / Applied Mathematical Modelling 37 (2013) 3777–3788
3779
2. The Trotter–Kato theorem Let X and Xn be Banach spaces with norm k k; k kn ; n ¼ 1; 2; . . ., respectively, and assume that TðtÞ is a C 0 -semigroup on X. A 2 GðM; x; XÞ; M P 1; x 2 R, means that A is the infinitesimal generator of a C 0 -semigroup TðtÞ; t P 0, satisfying kTðtÞk 6 Mext . Also assume that, for each n ¼ 1; 2; . . ., there exist bounded linear operators P n : X ! Xn and En : Xn ! X satisfying. (A1) kPn k 6 M 1 ; kEn k 6 M2 , where M 1 and M 2 are independent of n; (A2) kEn P n x xk ! 0 as n ! 1, for all x 2 X; (A3) Pn En ¼ In , where In is the identity operator on Xn .
Theorem 2.1 (Trotter–Kato). Assume that (A1) and (A3) are satisfied. Let A resp. An be in GðM; x; XÞ resp. in GðM; x; Xn Þ and let TðtÞ and T n ðtÞ be the semigroups generated by A and An on X and Xn , respectively. Then the following statements are equivalent: (a) There exists a k0 2 qðAÞ \
T1
n¼1
qðAn Þ such that, for all x 2 X,
kEn ðk0 In An Þ1 Pn x ðk0 In AÞ1 xk ! 0 as n ! 1: (b) for every x 2 X and t P 0,
kEn T n ðtÞPn x TðtÞxk ! 0 as n ! 1: uniformly on bounded t-intervals. If (a) or (b) is true, then (a) holds for all k with Rek > x. The assumption that, An 2 GðM; x; Xn Þ; n ¼ 1; 2; . . ., or equivalently, that kT n ðtÞkn 6 Mext is called the stability property of the approximation. Condition (a) is called the consistency property and condition (b) is called the convergency property of the approximation. Then, the Theorem 2.1 essentially states that, under the assumption of stability, consistency is equivalent to convergence. This is the functional analysis form of the Lax equivalent theorem (see [9]). In order to apply the Theorem 2.1 to a special problem, it is usually difficult to verify the resolvent convergence. To overcome this difficulty, reference [10] pointed out that condition (a) can be verified by the following theorem. Theorem 2.2. Let the assumptions of Theorem 2.1 be satisfied. Then statement (a) of the Theorem 2.1 is equivalent to (A2) and the following two statements: (C1) There exists a subset D dom A such that D ¼ X and ðk0 I AÞD ¼ X for a k0 > x; (C2) For all x 2 D there exist a sequence ðxn Þn2N with ðxn Þn2N 2 dom An such that
lim En xn ¼ x; and lim En An xn ¼ Ax:
n!1
n!1
Condition (C1) is the statement that D is the core of the operator A. Thus, condition (C2) is the statement that An converge pointwise on a core of A (see in [11, p. 88]). The goal of this paper is to develop a convergent approximation scheme for the system (1)–(8). We shall apply the Trotter–Kato theorem to obtain convergence. Thus, we formulate this system as an abstract system in an appropriate state space and construct a numerical scheme by employing a finite difference approximation of operator A. We consider the case with finite service time as that 0 6 x 6 M. Particular, the basic model is described by the equations,
dp0 ðtÞ ¼ kg1 p0 ðtÞ þ g2 dt
Z
M
lðxÞp1 ðx; tÞdx;
ð10Þ
0
@p1 ðx; tÞ @p1 ðx; tÞ þ ¼ ðkg1 þ lðxÞÞp1 ðx; tÞ; x 6 M; @t @x
ð11Þ
@pm ðx; tÞ @pm ðx; tÞ þ ¼ ðkg1 þ lðxÞÞpm ðx; tÞ þ kg1 pm1 ðx; tÞ; @t @x @pNþ2 ðx; tÞ @pNþ2 ðx; tÞ þ ¼ lðxÞpNþ2 ðx; tÞ þ kg1 pNþ1 ðx; tÞ; @t @x
x 6 M; 2 6 m 6 N þ 1;
x 6 M;
ð12Þ
ð13Þ
with boundary conditions,
p1 ð0; tÞ ¼ kg1 p0 ðtÞ þ q2
Z 0
M
lðxÞp1 ðx; tÞdx þ g2
Z 0
M
lðxÞp2 ðx; tÞdx;
ð14Þ
3780
H. Xu, W. Hu / Applied Mathematical Modelling 37 (2013) 3777–3788
Z
pm ð0; tÞ ¼ q2
M
lðxÞpm ðx; tÞdx þ g2
0
pNþ2 ð0; tÞ ¼ q2
Z
Z
M
lðxÞpmþ1 ðx; tÞdx; 2 6 m 6 N þ 1;
ð15Þ
0
M
lðxÞpNþ2 ðx; tÞdx;
ð16Þ
0
and initial conditions,
p0 ð0Þ ¼ 1;
pm ðx; 0Þ ¼ 0;
1 6 m 6 N þ 2:
ð17Þ T
Note that since p0 ðtÞ and pm ð; tÞ are probabilities, 0 6 p0 ðtÞ; pm ðx; tÞ 6 1, let ZðtÞ , ðp0 ðtÞ; p1 ð; tÞ; p2 ð; tÞ; . . . ; pNþ2 ð; tÞÞ , and define A on the domain
(
T
DðAÞ ¼
P ¼ ðp0 ; p1 ðxÞ; . . . ; pNþ2 ðxÞÞ 2 Xj
pm ðxÞ 2 W 1;1 ½0; M; and ðp1 ð0Þ; . . . ; pNþ2 ð0ÞÞT ¼ UP
)
:
3 2 RM 3 kg1 p0 þ g2 0 lðxÞp1 ðxÞdx p0 d 7 6 7 6 p1 ðxÞ 7 6 dx þ KðxÞ p1 ðxÞ 7 7 6 6 6 p ðxÞ 7 6 kg p ðxÞ d þ KðxÞ p ðxÞ 7 7 7 6 6 2 2 1 1 dx 7; 7¼6 by A6 .. 7 7 6 6 .. 7 7 6 6 . . 7 7 6 6 4 pNþ1 ðxÞ 5 6 kg p ðxÞ d þ KðxÞ p ðxÞ 7 5 4 Nþ1 1 N dx pNþ2 ðxÞ kg p ðxÞ d þ lðxÞ p ðxÞ 2
1 Nþ1
Nþ2
dx
where,
2
kg1 6 0 6 6 0 U :¼ 6 6 6 4 0
q2 w g2 w 0 0 q2 w g2 w 0 0 q2 w .. . 0 0 0
0 0 0
0 0 0
3 7 7 7 7; 7 7 5
and wðf Þ :¼
Z
M
lðxÞf ðxÞdx:
0
0 q2 w
Then system (10)–(16) can be rewritten in state space X ¼ R ðL1 ½0; MÞNþ2 as
(
_ ZðtÞ ¼ AZðtÞ;
ð18Þ
Zð0Þ ¼ ð1; 0; . . . ; 0Þ:
Therefore, using the approach given in [3,6,4], together with the reasonable assumption 0 6 lðxÞ 6 hM ðx 2 ½0; MÞ, we can establish the follow result, Corollary 2.3. A generates a C 0 -semigroup TðtÞ on state space X ¼ R ðL1 ½0; MÞNþ2 . Proof. For self contained and conciseness, we only sketch the proof of Corollary 2.3. Firstly, we divided operator A into two parts A ¼ A0 þ AB , where
2
kg1 6 0 6 6 0 A0 ¼ 6 6 6 4 0
3 0 0 0 d dx 0 0 7 7 d 0 7 0 dx 7; 7 7 .. 5 . d 0 0 dx
3 R 3 2 p0 g2 0M lðxÞp1 ðxÞdx 7 6 p1 ðxÞ 7 6 7 KðxÞp1 ðxÞ 7 6 6 7 6 p ðxÞ 7 6 7 g p ðxÞ K ðxÞp ðxÞ k 7 6 6 2 2 1 1 7 6 7¼6 AB 6 7: .. .. 7 6 6 7 7 6 6 . . 7 7 6 6 4 pNþ1 ðxÞ 5 4 kg p ðxÞ KðxÞp ðxÞ 7 5 Nþ1 1 N pNþ2 ðxÞ kg1 pNþ1 ðxÞ lðxÞpNþ2 ðxÞ 2
1 We can prove that kðcI A0 Þ1 k < ch when c > hM , where hM satisfies 0 6 lðxÞ 6 hM . Moreover, it is easy to prove that M
DðA0 Þ=DðAÞ is dense in Banach space X. Hille–Yosida theorem (see [11]) yields that operator A0 generates a C 0 -smeigroup. Since AB is bounded operator, with perturbation theory of semigroup, we know that A generates a C 0 -semigroup. Thus we complete the proof of Corollary 2.3. In the following parts of the paper, we denote the C 0 -semigroup generated by A as TðtÞ. h 3. An approximation scheme In this section we introduce a finite difference method and employ the Trotter–Kato Theorem to the abstract version (18) to obtain convergence. Divide the interval ½0; M into n equal subintervals to get ½xk1 ; xk , and set Dx ¼ Mn, so that xk ¼ kDx; k ¼ 1; 2; . . . ; n. Then 0 ¼ x0 < x1 < < xk < < xn ¼ M.
H. Xu, W. Hu / Applied Mathematical Modelling 37 (2013) 3777–3788
3781
Applying a forward difference in space yields the approximating system n X dp0 ðtÞ ¼ kg1 p0 ðtÞ þ g2 Dx lðxk Þp1 ðxk ; tÞ; dt k¼1
ð19Þ
dp1 ðxk ; tÞ p1 ðxk1 ; tÞ p1 ðxk ; tÞ ¼ ðkg1 þ lðxk ÞÞp1 ðxk ; tÞ; dt Dx
ð20Þ
dpm ðxk ; tÞ pm ðxk1 ; tÞ pm ðxk ; tÞ ¼ ðkg1 þ lðxk ÞÞpm ðxk ; tÞ þ kg1 pm1 ðxk ; tÞ; dt Dx
ð21Þ
dpNþ2 ðxk ; tÞ pNþ2 ðxk1 ; tÞ pNþ2 ðxk ; tÞ ¼ lðxk ÞpNþ2 ðxk ; tÞ þ kg1 pNþ1 ðxk ; tÞ; p1 ð0; tÞ ¼ kg1 dt Dx " # n n X X p0 ðtÞ þ Dx q2 lðxk Þp1 ðxk ; tÞ þ g2 lðxk Þp2 ðxk ; tÞ ; k¼1
" pm ð0; tÞ ¼ Dx q2
n X
n X
k¼1
k¼1
lðxk Þpm ðxk ; tÞ þ g2
#
lðxk Þpmþ1 ðxk ; tÞ ;
ð24Þ
n X
lðxk ÞpNþ2 ðt; xk Þ;
ð25Þ
k¼1
where 2 6 m 6 N þ 1; pi ðt; xk Þ is the value of pi ðt; xÞ at the kth nodal point xk ¼ kDx, and so the same as lðxk Þ. From Eqs. (19)–(25), it is clear that the approximating generators An on Xn ¼ R ðRn ÞNþ2 have the matrix form,
An ¼ ½Anij ðNþ3ÞðNþ3Þ ; where
An12 ¼ Dxg2 ½lðx1 Þ; . . . ; lðxn Þ1n ;
An13 ¼ An14 ¼ ¼ An1ðNþ3Þ ¼ ½01n ; An21 ¼
T kg1 ; ; 0; . . . ; 0 Dx 1n 2 1
An22
6 6 6 6 ¼6 6 6 4
Dx
An24 ¼ An25 ¼ ¼ An2ðNþ3Þ ¼ ½0nn ;
Kðx1 Þ þ q2 lðx1 Þ 1 Dx
q2 lðx1 Þ 1 Dx
q2 lðx2 Þ
q2 lðxn Þ
0
0
Kðx3 Þ .. .
0
Kðx2 Þ
0
1 Dx
0
0
1 Dx
0
1 Dx
3 7 7 7 7 7; 7 7 5
Kðxn Þ
here; Kðxi Þ ¼ kg1 þ lðxi Þ; 2 6 6 An23 ¼ 6 6 4
g2 lðx1 Þ g2 lðx2 Þ g2 lðx3 Þ g2 lðxn Þ 0
0
0
.. .
0
0
0
0
0
Ann1 ¼ ½0n1 ; Annðnþ1Þ ¼ An23 ;
Ann2 ¼ ¼ Annðn2Þ ¼ ½0nn ;
3 7 7 7; 7 5
Annðn1Þ ¼ diagðkg1 Þnn ;
Annðnþ2Þ ¼ ¼ AnðNþ3Þ ¼ ½0nn ;
AnðNþ3Þ1 ¼ ð0Þn1 ;
ð23Þ
k¼1
pNþ2 ð0; tÞ ¼ q2 Dx
An11 ¼ ½kg1 11 ;
ð22Þ
n ¼ 3; 4; . . . ; N þ 2;
AnðNþ3Þ2 ¼ ¼ AnðNþ3ÞðNþ1Þ ¼ ½0nn ;
Annn ¼ An22 ;
3782
H. Xu, W. Hu / Applied Mathematical Modelling 37 (2013) 3777–3788
AnðNþ3Þ;ðNþ2Þ ¼ diagðkg1 Þnn ; and; 2 1 AnðNþ3Þ;ðNþ3Þ
6 6 6 6 ¼6 6 6 4
Dx
lðx1 Þ þ q2 lðx1 Þ 1 Dx
q2 lðx1 Þ 1 Dx
q2 lðx2 Þ
q2 lðxn Þ
0
0
lðx3 Þ .. .
0
lðx2 Þ
0
1 Dx
0
0
1 Dx
0
1 Dx
3 7 7 7 7 7: 7 7 5
lðxn Þ
Since An is the infinitesimal generator of a C 0 -semigroup of contraction T n ðtÞ on Xn , let Pn ; En ; k kn be defined as
En ðyÞ ¼
p0 ;
n X
p1 ðxk Þvðxk1 ;xk ; . . . ;
k¼1
! n X pNþ2 ðxk Þvðxk1 ;xk ; k¼1
Z x1 Z x1 Z xn Z xn 1 1 1 1 p1 ðxÞdx; . . . ; p1 ðxÞdx ; . . . ; pNþ2 ðxÞdx ; . . . ; pNþ2 ðxÞdxÞ ; ðPn /Þk ¼ p0 ; Dx x0 Dx xn1 Dx x 0 Dx xn1 n n X X kykn ¼ jp0 j þ Dx jp1 ðxk Þj þ þ Dx jpNþ2 ðxk Þj; k¼1
k¼1
where
y ¼ ðp0 ; ðp1 ðx1 Þ; . . . ; p1 ðxn ÞÞ; . . . ; ðpNþ2 ðx1 Þ; . . . ; pNþ2 ðxn ÞÞÞT 2 Xn ; / ¼ ðp0 ; p1 ðxÞ; . . . ; pNþ2 ðxÞÞT 2 X: Then, we have the following theorem, which is also the main theorem of the paper. Theorem 3.1. Let En ; Pn ; A and An be defined as above, and let TðtÞ and T n ðtÞ be the semigroup generated by A and An on X and Xn , respectively, then for every x 2 X and t P 0; kEn T n ðtÞPn x TðtÞxk ! 0 as n ! 1 uniformly on bounded t-intervals. To prove Theorem 3.1, we only need to verify that En and P n satisfy (A1)–(A3), and prove that both (C1) and (C2) in Theorem 2.2 hold.Firstly, we prove that En and Pn satisfy (A1)–(A3) 8y ¼ ðp0 ; ðp1 ðx1 Þ; . . . ; p1 ðxn ÞÞ; . . . ; ðpNþ2 ðx1 Þ; . . . ; pNþ2 ðxn ÞÞÞT 2 Xn , since x0 ¼ 0 and xn ¼ M
kEn yk ¼ jp0 j þ
Nþ2 Z X
M
0
i¼1
X n Nþ2 X n Z xk Nþ2 X X
jpi ðxk Þjdx 6 jp0 j þ kpi ðxÞkn ¼ kykn ;
pi ðxk Þvðxk1 ;xk dx 6 jp0 j þ
k¼1 i¼1 k¼1 xk1 i¼1
which yields kEn k 6 1. 8/ðxÞ ¼ ðp0 ; p1 ðxÞ; . . . ; pNþ2 ðxÞÞT 2 X,
kPn /kn ¼ jp0 j þ
Nþ2 X
Dx
i¼1
!
Z !
n
Nþ2 X n Z xk xk X X
1 pi ðxÞdx 6 jp0 j þ jpi ðxÞjdx 6 k/ðxÞk:
Dx xk1 i¼1 k¼1 k¼1 xk1
This implies kPn k 6 1. Therefore, (A1) is satisfied.
kEn Pn / /k ¼ ¼ ¼
XNþ2 Z i¼1
XNþ2 i¼1
XNþ2 i¼1
X n
ðP p ðxÞÞ v p ðxÞ
dx
n i i k ðxk1 ;xk i¼1
0 0 k¼1
!
!
Z
n Z xk X n n Z xk
xk X X X Nþ2
1 dx ¼ ðP p ðxÞÞ v p ðxÞ p ðxÞdx p ðxÞ
dx
n i i i i k ðxk1 ;xk i¼1
D x x x x k1 k¼1 k1 k1 k¼1 k¼1 ! ! Z Z n n x x XNþ2 X k X k jpi ðnk Þ pi ðxÞjdx ; nk 2 ðxk1 ; xk Þ 6 Dx jp0i ðgk Þjdx ; i¼1 M
jEn Pn pi ðxÞ pi ðxÞjdx ¼
k¼1
xk1
XNþ2 Z
M
k¼1
xk1
gk 2 ðxk1 ; xk Þ ! 0; as n ! 1. This follows since p0i ðxÞ 2 L1 ½0; MÞ. Thus, (A2) is satisfied. Pn En ¼ In is obvious. Therefore, A(3) follows. Secondly, we will verify the consistency property, that is, we will prove that both (C1) and (C2) hold. To this end, we choose D ¼ dom A, which establishes condition (C1) in Theorem 2.2 with x ¼ 0, the proof is rather straightforward and will therefore be omitted. 8/ ¼ ðp0 ; p1 ðxÞ; . . . ; pNþ2 ðxÞÞT 2 dom A X, we define /n 2 Xn by
3783
H. Xu, W. Hu / Applied Mathematical Modelling 37 (2013) 3777–3788
/n ¼ ðp0 ; ðp1 ðx1 Þ; . . . ; p1 ðxn ÞÞ; . . . ; ðpNþ2 ðx1 Þ; . . . ; pNþ2 ðxn ÞÞÞT ; kEn /n /kL1 ¼
n Z X
xk
jEn /n /jdx ¼
xk1
k¼1
n Z X
xk
j/ðxk Þ /jdx ¼
n Z X
xk1
k¼1
xk
j/0 ðnðxÞÞðxk xÞjdx 6 Dx
xk1
k¼1
n Z X k¼1
xk
j/0 ðnðxÞÞjdx
xk1
6 Dxk/0 kL1 ; which proves limn!1 En /n ¼ /. To investigate kEn An /n A/kX , we first consider the expression of An /n . By the expression of An , we know that n A /n =½ai ðnðNþ2Þþ1Þ1 , and
a1 ¼ kg1 p0 þ Dxg2
n X
lðxk Þp1 ðxk Þ;
k¼1
a2 ¼
n n X X kg1 1 q2 lðxk Þp1 ðxk Þ þ g2 lðxk Þp2 ðxk Þ; p0 þ Kðx1 Þ p1 ðx1 Þ þ Dx Dx k¼1 k¼1
aiþ1 ¼
1 1 p1 ðxði1Þ Þ þ Kðxi Þ p1 ðxi Þ; Dx Dx
amnþ2 ¼ kg1 pm ðx1 Þ
amnþiþ1 ¼ kg1 pm ðxi Þ þ
i ¼ 2; 3; . . . ; n;
n n X X 1 q2 lðxk Þpmþ1 ðxk Þ þ g2 lðxk Þpmþ2 ðxk Þ; m ¼ 1; 2; . . . ; N; þ Kðx1 Þ pmþ1 ðx1 Þ þ Dx k¼1 k¼1
1 1 pmþ1 ðxi1 Þ þ Kðxi Þ pmþ1 ðxi Þ; Dx Dx
aðNþ1Þnþ2 ¼ kg1 pNþ1 ðx1 Þ
i ¼ 2; . . . ; n; m ¼ 1; 2; . . . ; N;
aðNþ1Þnþiþ1 ¼ kg1 pNþ1 ðxi Þ þ
n X 1 q2 lðxk ÞpNþ2 ðxk Þ; þ lðx1 Þ pNþ2 ðx1 Þ þ Dx k¼1
1 1 pNþ2 ðxi1 Þ þ lðxi Þ pNþ2 ðxi Þ; Dx Dx
i ¼ 2; . . . ; n:
Then En An /n ¼ ðc1 ; c2 ðxÞ; . . . ; cNþ3 ðxÞÞT can be expressed as
a1 ;
n X
akþ1 vðxk1 ;xk ;
k¼1
n n X X a1nþkþ1 vðxk1 ;xk ; . . . ; aðNþ1Þnþkþ1 vðxk1 ;xk k¼1
!T :
k¼1
Together with definition of operator A, it is ready to see that A/ can be expressed as A/ ¼ ðbj ÞðNþ3Þ1 ¼ ðb1 ; b2 ðxÞ; . . . ; bNþ3 ðxÞÞT .
b1 ¼ kg1 p0 þ g2
Z
M
lðxÞp1 ðxÞdx ¼ kg1 p0 þ
0
n X
g2
k¼1
Z
xk
lðxÞp1 ðxÞdx;
xk1
d b2 ðxÞ ¼ þ KðxÞ p1 ðxÞ; dx bi ðxÞ ¼ kg1 pi2 ðxÞ
d þ KðxÞ pi1 ðxÞ; dx
bNþ3 ðxÞ ¼ kg1 pNþ1 ðxÞ
ð3 6 i 6 N þ 2Þ
d þ lðxÞ pNþ2 ðxÞ: dx
The initial value of each bi ðxÞ is as follows
3784
H. Xu, W. Hu / Applied Mathematical Modelling 37 (2013) 3777–3788
b2 ð0Þ ¼ kg1 p0 ðtÞ þ q2 ¼ kg1 p0 ðtÞ þ q2
Z
M
lðxÞp1 ðxÞdx þ g2
0 n Z xk X
xk1
k¼1
bi ð0Þ ¼ q2
Z
bNþ3 ð0Þ ¼ q2
Z
n Z X
0
xk
lðxÞp2 ðxÞdx;
xk1
M
lðxÞpi ðxÞdx; ð3 6 i 6 N þ 2Þ ¼ q2
M
0
lðxÞp2 ðxÞdx
0
k¼1
lðxÞpi1 ðxÞdx þ g2 Z
M
lðxÞp1 ðxÞdx þ g2
M
0
Z
lðxÞpNþ2 ðx; tÞdx ¼ q2
n Z X k¼1
n Z X k¼1
xk xk1
xk
xk1
lðxÞpi1 ðxÞdx þ g2
n Z X k¼1
xk
xk1
lðxÞpi ðxÞdx;
lðxÞpNþ1 ðxÞdx:
Therefore, we have
kEn An /n A/kXM ¼ jc1 b1 j þ
Nþ2 X kci ðxÞ bi ðxÞk:
ð26Þ
i¼2
We can prove that
"Z " #
#
Z xk
n n
X
X xk
jc1 b1 j ¼ g2
Dxlðxk Þp1 ðxk Þ lðxÞp1 ðxÞdx ¼ g2
ðlðxk Þp1 ðxk Þ lðxÞp1 ðxÞÞdx
k¼1
k¼1 xk1 xk1
" #
Z
X n xk
6 g2
sup ðlðx þ zÞp1 ðx þ zÞ lðxÞp1 ðxÞÞdx ; jzj 6 Dx
i¼1 z xk1
" #
Z M
n Z xk X
ðlðx þ zÞp1 ðx þ zÞ lðxÞp1 ðxÞÞdx 6 g2 sup ðlðx þ zÞp1 ðx þ zÞ lðxÞp1 ðxÞÞdx
6 g2 sup
z i¼1 xk1 z 0 6 g2 xðlðxÞp1 ðxÞ; DxÞ ! 0; as Dx ! 0;
RM where xðlðxÞpj ðxÞ; DxÞ ¼ supf 0 jlðx þ zÞpj ðx þ zÞ lðxÞpj ðxÞjdx : kzk 6 Dxg is called modulus of continuity for
kc2 ðxÞ b2 ðxÞk ¼
Z
M
jc2 ðxÞ b2 ðxÞjdx ¼
0
n Z X k¼1
xk
xk1
jc2 ðxÞ b2 ðxÞjdx ¼
n Z X k¼1
xk
xk1
lðxÞpj ðxÞ.
jakþ1 vðxk1 ;xk ðxÞ þ p01 ðxÞ þ KðxÞp1 ðxÞjdx
kg n X 1
1 0 ðq2 lðxk Þp1 ðxk Þ þ g2 lðxk Þp2 ðxk ÞÞ þ p1 ðxÞ þ KðxÞp1 ðxÞ dx p0 þ Kðx1 Þ p1 ðx1 Þ þ 6
D x D x x0 k¼1
Z n X xk 1
p1 ðxk1 Þ 1 þ Kðxk Þ p1 ðxk Þ þ p0 ðxÞ þ KðxÞp1 ðxÞ dx þj 1
Dx
D x k¼2 xk1 Z x1 Z x1 1 6 ½jDxp01 ðxÞ p01 ðsÞdsj þ DxjKðx1 ÞÞp1 ðx1 Þ KðxÞp1 ðxÞjdx Dx x 0 x0
# "
Z xk
n Z xk
1X
0 0 þ p1 ðsÞdsj þ DxjKðxk Þp1 ðxk Þ KðxÞp1 ðxÞ dx
Dxp1 ðxÞ
Dx k¼2 xk1
xk1
# " Z Z
n xk
x 1X
k 0 6 ðp1 ðxÞ p01 ðsÞÞdsj þ DxjKðxk Þp1 ðxk Þ KðxÞp1 ðxÞ dx
Dx k¼1 xk1 xk1 Z n xk X 6 jðp01 ðxÞ p01 ðnk ÞÞjdx þ xðKðxÞp1 ðxÞ; DxÞ Z
x1
k¼1
6 sup z
xk1
n Z X k¼1
xk
xk1
jðp01 ðx þ zÞ p01 ðxÞÞjdx þ xðKðxÞp1 ðxÞ; DxÞ 6 xðp01 ðxÞ; DxÞ þ xðKðxÞp1 ðxÞ; DxÞ ! 0; as Dx ! 0;
RM where xðp01 ðxÞ; DxÞ ¼ supf 0 jp01 ðx þ zÞ p01 ðxÞjdx : kzk 6 Dxg is called modulus of continuity for p01 ðxÞ. Again, we can prove that kci ðxÞ bi ðxÞk; i ¼ 3; . . . ; N þ 3 tends to 0 as Dx ! 0, therefore, kEn An /n A/kX ! 0 as Dx ! 0. Thus, (C2) in Theorem 2.2 holds. With Theorem 2.1, Theorem 2.2 and the fact that operator A generates a C 0 -semigroup TðtÞ of contraction in X, and operator An generates a C 0 -semigroup T n ðtÞ of contraction in Xn , we know that Theorem 3.1 holds. Combining with Theorem 3.1, we have following corollary. Corollary 3.2. Since the dynamical solution of system (10)–(16) can be expressed as p ¼ TðtÞu0 , here u0 is the initial value of the system. Then, we have following relationship,
kEn T n ðtÞPn u0 TðtÞu0 k ! 0 as n ! 1;
H. Xu, W. Hu / Applied Mathematical Modelling 37 (2013) 3777–3788
Fig. 2. Dynamical solution of p0 ðtÞ with different g2 .
Fig. 3. Dynamical solution of pNþ2 ðtÞ with different g2 .
Fig. 4. Dynamical solution of pNþ2 ðx; tÞ with exponential distribution.
3785
3786
H. Xu, W. Hu / Applied Mathematical Modelling 37 (2013) 3777–3788
Fig. 5. Dynamical solution of pNþ2 ðtÞ.
Fig. 6. Dynamical solution of pNþ2 ðx; tÞ with Weibull distribution.
n
where T n ðtÞ is semi-group with generator An , and T n ðtÞ ¼ eA t . Proof. In Corollary 2.3, we have verified that A generates C 0 -semigroup TðtÞ. With Eq. (18), we know that if u0 is the initial value of the system (10)–(16), the dynamical solution could be expressed as p ¼ TðtÞu0 (see [11]). Moreover, An is bounded operator, which generates a uniformly continuous semigroup. There for, the conditions of Theorem 3.1 are satisfied, and
kEn T n ðtÞPn u0 TðtÞu0 k ! 0 as n ! 1: holds. Corollary 3.2 provides us with a useful tool to calculate the dynamical solution and instantaneous indexed of the model. h
H. Xu, W. Hu / Applied Mathematical Modelling 37 (2013) 3777–3788
3787
4. Numerical results In this section, some numerical experiments are presented to show the dynamical solutions of the system. To simplicity the simulation, here we assume that N ¼ 1. Although this is not the most general form of the system, the framework of simulating dynamical solution shown in Section 2 is sufficient to illustrate the issues raised in this paper. Here, we assume g2 is a variety which can be adjusted in the process. Fig. 2 shows the dynamical behavior of p0 ðtÞ with different g2 . p0 ðtÞ means the probability that at time t the machine 1 is processing a piece and there is no piece in the buffer and in the machine 2. With Fig. 2, we know that the value of p0 ðtÞ is decreasing to a steady value, but the steady value will vary with different g2 . Moreover, if g2 decrease, the corresponding steady value will decrease too. Fig. 3 describes the dynamR1 ical behavior of pNþ2 ðtÞ ¼ 0 pNþ2 ðx; tÞdx with different g2 . From the introduction section of this paper, we know that pNþ2 ðtÞ means the probability that at time t the buffer is full. From Fig. 3, we can find that each of the dynamical solution of pNþ2 ðtÞ converges and the steady value will increase as g2 decreasing. In terms of the method of finite differential scheme, we can present the solution pNþ2 ðx; tÞ as a 3D figure, just as shown in Fig. 4. The parameters in Fig. 4 are, k ¼ 0:4; g1 ¼ 0:6; g2 ¼ 0:8 and lðxÞ ¼ 0:6. Here lðxÞ ¼ 0:6 means that the service time of machine 2 follows exponential distribution with parameter equalling to 0.6. Other situations, like g2 ¼ 0:2, or g2 ¼ 0:4, can also be simulated. Although, in Figs. 2–4, lðxÞ is assumed to be constant, we can simulate the dynamical solution if lðxÞ is function of x. For example, let service time of machine 2 follows Weibull distribution with hazard rate lðxÞ ¼ lbðg2 xÞb1 , here l ¼ 0:6 and shape parameter b ¼ 2. With the numerical method presented in this paper, we can simulate the dynamical behavior of pNþ2 ðtÞ. The corresponding parameters used in Fig. 5 and Fig. 6 are, k ¼ 0:4; g1 ¼ 0:6; g2 ¼ 0:8, and lðxÞ ¼ 1:2 ð0:6 xÞ. Obviously, lðxÞ is bounded in interval ½0; M. The physical meaning of pNþ2 ðtÞ is the probability that the buffer of the queueing system is full at time t. From Fig. 5, we can find that pNþ2 ðtÞ increases and converges to a certain value with time t grows. Moreover, Fig. 6 just reflects the dynamical behavior of solution pNþ2 ðx; tÞ with Weibull distribution in a 3-D space. 5. Conclusion There are several papers concerning the well-posedness problem of reliable system and repairable system (see [12,13] and reference therein). However, none of these papers deal with numerical approximations. We formulated a finite difference scheme to simulate the CIMS system and prove convergence by Trotter–Kato theorem. The same method can also be applied to other reliable systems. From the viewpoint of practice, every dynamical solution has its own physical means, for example, p0 ðtÞ means the probability the buffer is empty and machine 2 is idle, pNþ2 ðtÞ means the probability that buffer is full. Moreover, many instantaneous index of system can be deduced directly from these dynamical solutions, like instantaneous availability, instantaneous queue length and so on. Therefore, it is an important research work to study the dynamical solution of the system. If we denote the solution of Eqs. (10)–(17) as p ¼ ðp0 ðtÞ; p1 ðx; tÞ; . . . ; pNþ2 ðx; tÞÞT and extend p by letting pi ðx; tÞ ¼ 0, when x > M; i ¼ 1; 2; . . . ; N þ 2. Then, it can be verified that
kp pkB ! 0; as M ! 1; where p ¼ ðp0 ðtÞ; p1 ðx; tÞ; . . . ; pNþ2 ðx; tÞÞT is the dynamical solution of Eqs. (1)–(8). Acknowledgements The authors thank Professor John Burns for his valuable help and advice on this paper. This work is supported by NSFC (11001013) and supported by the Interdisciplinary Center for Applied Mathematics during the first author’s visit at Virginia Tech. This work is also supported by the Air Force Office of Scientific Research Grants FA9550-07-1-0273 and FA9550-10-10201 and by the Environmental Security Technology Certification Program (ESTCP) under a subcontract from United Technologies. The authors would also like to thank the referees for their helpful comments and suggestions. Appendix A. Supplementary data Supplementary data associated with this article can be found, in the online version, at http://dx.doi.org/10.1016/ j.apm.2012.07.056. References [1] S.G. Shu, An analysis of the repairable CIMS with buffers and a study of system reliability, Acta Autom. Sin. 18 (1992) 15–22. [2] W. Li, J.H. Cao, Some performance measures of transfer line consisting of two unreliable machines with reprocess rule, J. Syst. Sci. Syst. Eng. 7 (1998) 283–292. [3] G. Gupur, Well-posedness of a reliability model, Acta Anal. Func. Appl. 5 (2003) 193–209. [4] H.B. Xu, W.H. Guo, J.Y. Yu, G.T. Zhu, The asymptotic stability of series CIMS with a single storage buffer, Syst. Eng. Theory Practice 5 (2004) 91–96. [5] W.W. Hu, H.B. Xu, J.Y. Yu, G.T. Zhu, Exponential stability of a reparable multi-state device, J. Syst. Sci. Complexity 20 (2007) 437–443. [6] A. Haji, A. Radl, A semigroup approach to queueing system, Semigroup Forum 75 (2007) 609–623.
3788 [7] [8] [9] [10] [11] [12] [13]
H. Xu, W. Hu / Applied Mathematical Modelling 37 (2013) 3777–3788
H.T. Banks, J.A. Burns, Hereditary control problems: numerical methods based on average approximations, SIAM J. Control Optim. 16 (1978) 169–208. H.T. Banks, F. Kappel, Spline approximations for functional differential equations, J. Differ. Equ. 34 (1979) 496–522. P.D. Lax, R.D. Richtmyer, Survey of the stability of linear finite differential equations, Commun. Pure Appl. Math. IX (1956) 267–293. K. Ito, K. Kappel, The Trotter–Kato theorem and approximation of PDEs, Math. Comput. 69 (1998) 21–44. A. Pazy, Semigroups of Linear Operator and Applications to Partial Differential Equations, Springer Press, 1983. H.B. Xu, J.Y. Yu, G.T. Zhu, Asymptotic property of a reparable multi-state device shape Quarterly of Applied Mathematics, LXIII 2005, 779–789. H.B. Xu, W.W. Hu, Availability optimisation of repairable system with preventive maintenance policy, Int. J. Syst. Sci. 39 (2008) 655–664.