Mokhtar Kirane
email: mokhtar.kirane@univ-lr.fr

Salman A Malik
email: salman.malik@univ-lr.fr

Profile of blowing-up solutions to a nonlinear system of fractional differential equations

Keywords: Nonlinear system, Fractional derivative, Integral equations, Blow up time, Laplace transform. 2010 MSC: 26A33, 45J05, 34K37

We investigate the profile of the blowing up solutions to the nonlinear nonlocal system (FDS)

where

and D β 0 + , 0 < α < 1, 0 < β < 1 stand for the Riemann-Liouville fractional derivatives. Our method of proof relies on comparisons of the solution to the (FDS) with solutions of the subsystems obtained from (FDS) by dropping either the usual derivatives or the fractional derivatives.

Introduction

In this paper, we are concerned with the profile of the blowing up solutions to the nonlinear system of fractional differential equations (FDS)

u (t) + D α 0 + (u -u 0 )(t) = |v(t)| q , t > 0, v (t) + D β 0 + (v -v 0 )(t) = |u(t)| p , t > 0,
subject to the initial conditions

u(0) = u 0 , v(0) = v 0 , (1)
where p > 1, q > 1, u(0) = u 0 > 0, v(0) = v 0 > 0 are constants, D α 0 + and D β 0 + stand for the Riemann-Liouville fractional derivatives of orders 0 < α < 1 and 0 < β < 1, respectively. It has been proven by Furati and Kirane in [START_REF] Furati | Necessary conditions for the existence of global solutions to systems of fractional differential equations[END_REF] that, under the constraint 1 -1/pq ≤ α + β/p or 1 -1/pq ≤ β + α/q, any solution emerging from (u 0 , v 0 ), (u 0 > 0, v 0 > 0) blows up in a finite time. Moreover, they gave an upper bound for the blow up time. Their method of proof relies on the weak formulation of the system (FDS) with a judicious choice of the test function. Here, we present estimates from above and below of the blowing up solutions by comparing them with the solutions of the particular systems obtained from (FDS):

u (t) = λ|v(t)| q , t > 0, p > 1, v (t) = λ|u(t)| p , t > 0, q > 1, (2)
the system of ordinary differential equations (ODS) with either λ = 1 or λ = 1/2, and

D α 0 + (u -u 0 )(t) = µ|v(t)| q , t > 0, D β 0 + (v -v 0 )(t) = µ|u(t)| p , t > 0, (3)
the system of fractional differential equations (PFDS), with either µ = 1 or µ = 1/2. As the solutions to (ODS) are easy to obtain, we will investigate in particular the system (PFDS).

Our paper is organized as follows: some definitions and known results needed in the sequel are stated in the next section. Section 3 is devoted to our main results for the investigation of the profile of the blowing up solutions of system (FDS). In the last section, we provide numerical approximations of the solutions, and the profiles are given for particular examples.

Preliminaries

In this section, the necessary definitions and notations from fractional calculus needed for the paper are recalled for the sake of the reader.

The Riemann-Liouville fractional integral of order 0 < α < 1 is

J α 0 + f (t) := 1 Γ(α) t 0 f (τ) (t -τ) 1-α dτ, t > 0, (4)
where Γ(α) is the Euler Gamma function and f : R + → R is a real valued integrable function. The integral (4) can be written as a convolution

J α 0 + f (t) = (φ α f )(t), (5)
where

φ α := t α-1 Γ(α) , t ≥ 0, 0 t < 0. ( 6
)
The Laplace transform of the Riemann-Liouville integral of order 0 < α < 1 is

L{J α 0 + f (t) : s} = L{ f (t) : s}/s α .
The Riemann-Liouville fractional derivative of order 0 < α < 1 is

D α 0 + f (t) := d dt J 1-α 0 + f (t) = 1 Γ(1 -α) d dt t 0 f (τ) (t -τ) α dτ. ( 7
)
In particular, when α = 0, D 0 0 + f (t) = f (t); note that the Riemann-Liouville fractional derivative of a constant is not equal to zero.

The Caputo fractional derivative of order 0 < α < 1 is defined by

C D α 0 + f (t) := J 1-α 0 + d dt f (t) = 1 Γ(1 -α) t 0 f (τ) (t -τ) α dτ. ( 8
)
The relation between the Riemann-Liouville and the Caputo fractional derivatives is

C D α 0 + f (t) = D α 0 + f (t) -f (0) , 0 < α < 1. (9)
The formula for the integration by parts in [0, T ] is given by

T 0 f (t)D α 0 + g(t)dt = T 0 D α T -f (t)g(t)dt, (10)
see [START_REF] Kilbas | Theory and Applications of Fractional Differential Equations[END_REF], where D α Tis the right sided fractional derivative defined by

D α T -f (t) := - 1 Γ(1 -α) d dt T t f (τ) (τ -t) α dτ, 0 < α < 1. ( 11
)
The Laplace transforms of the Riemann-Liouville derivative D α 0 + f (t) and the Caputo derivative C D α 0 + f (t) are

L{D α 0 + f (t) : s} = s α L{ f (t) : s} -J 1-α 0 + f (0), 0 < α < 1, ( 12
)
L{ C D α 0 + f (t) : s} = s α L{ f (t) : s} -s α-1 f (0), 0 < α < 1. ( 13
)
At this stage let us recall the test function considered in [START_REF] Furati | Necessary conditions for the existence of global solutions to systems of fractional differential equations[END_REF]:

ϕ(t) = T -λ (T -t) λ , t ∈ [0, T ], 0 t > T ; (14) it satisfies T 0 D α T -ϕ(t)dt = C α,λ T 1-α , C α,λ = Γ(λ + 1) Γ(λ -α + 2) , (15)
and for p > 1, λ > αp -

1 T 0 ϕ 1-p (t) D α T -ϕ(t) p dt = C p,α T 1-αp , (16)
where

C p,α = 1 λ -pα + 1 Γ(λ + 1) Γ(λ + 1 -α) p . ( 17
)
The main result in [START_REF] Furati | Necessary conditions for the existence of global solutions to systems of fractional differential equations[END_REF] is the Theorem 2.1 ( [START_REF] Furati | Necessary conditions for the existence of global solutions to systems of fractional differential equations[END_REF]). Suppose that 0 < α, β < 1, p, q > 1 and u 0 > 0, v 0 > 0, then solutions to system (FDS) subject to (1) blow up in a finite time, whenever

1 - 1 pq ≤ α + β p , or 1 - 1 pq ≤ β + α q .
Let us recall that the blowing up solution of the system (ODS) with λ = 1

u (t) = |v(t)| q , t > 0, p > 1, v (t) = |u(t)| p , t > 0, q > 1, (18)
subject to the initial conditions (u(0

) = u 0 > 0, v(0) = v 0 > 0) is u od (t) = C 1 T max -t -q+1 pq-1 , v od (t) = C 2 T max -t -p+1 pq-1 , (19)
where T max , C 1 and C 2 are given by

T max = C 2 u 0 C 1 v 0 pq-1 p-q , C 1 = (p + 1) q (q + 1) (pq -1) q+1 1 pq-1 , C 2 = (p + 1)(q + 1) p (pq -1) p+1 1 pq-1 . (20)
Alike, the blowing up solution of the system (ODS) with λ = 1/2

u (t) = 1 2 |v(t)| q , t > 0, u(0) = u 0 > 0, v (t) = 1 2 |u(t)| p , t > 0, v(0) = v 0 > 0, (21)
for p, q > 1, is given by

u * od (t) = C * 1 T * max -t -q+1 pq-1 , v * od (t) = C * 2 T * max -t -p+1 pq-1 ,
where T * max , C * 1 and C * 2 are given by

T * max = C * 2 u 0 C * 1 v 0 pq-1 p-q , C * 1 = 2 q+1 pq-1 (p + 1) q (q + 1) (pq -1) q+1 1 pq-1 , C * 2 = 2 p+1 pq-1 (p + 1)(q + 1) p (pq -1) p+1 1 pq-1 .
Notice that the solution of the system (18) blows up earlier than the solution of the system (21) since T * max = 2T max .

Main Results

In this section, we consider the nonlinear system of fractional differential equations (PFDS)

D α 0 + (u -u 0 )(t) = µ|v(t)| q , t > 0, D β 0 + (v -v 0 )(t) = µ|u(t)| p , t > 0,
for 0 < α, β < 1, with either µ = 1 or µ = 1/2, subject to the initial conditions [START_REF] Furati | Necessary conditions for the existence of global solutions to systems of fractional differential equations[END_REF]. For µ = 1 we will give the necessary conditions for solutions to blow up in a finite time; moreover an upper bound on the blow up time is provided. Let us set

s = pα + β 1 -pq , s = α + qβ 1 -pq , C v = K 2 (K 1 ) p pq-1 , C u = K * 2 (K * 1 ) q pq-1 ,
where K 1 , K 2 , K * 1 and K * 2 are:

K 1 = C 1 p p ,α C 1 pq q ,β , K 2 = C 1 q q ,β C -1 β,λ , K * 1 = C 1 q q ,β C 1 qp p ,α , K * 2 = C 1 p p ,α C -1
α,λ , and the constants {C α,λ , C β,λ } and {C p ,α , C q ,β } are given by ( 15) and ( 17) respectively.

The first result is given by

Theorem 3.1. Let p, q > 1 and u 0 > 0, v 0 > 0. If 1 - 1 pq ≤ α + β p , or 1 - 1 pq ≤ β + α q ,
then any solution to the system (PFDS) with µ = 1 subject to the initial conditions (1) blows up in a finite time T max . Furthermore, an upper bound on the blow up time T max is given by min{T u , T v } where

T v = v 0 C v 1/s T u = u 0 C u 1/ s.
Proof. The proof is by contradiction; suppose (u, v) is a global solution of the system (PFDS). Multiplying the equations of system (PFDS) by ϕ and integrating over [0, T ], we, by virtue of (10), obtain

T 0 uD α T -ϕ = u 0 T 0 D α T -ϕ + T 0 |v| q ϕ, ( 22
) T 0 vD β T -ϕ = v 0 T 0 D β T -ϕ + T 0 |u| p ϕ, (23)
for any function ϕ such that D α Tϕ exists and ϕ(T ) = 0. Using Hölder's inequality, we obtain

T 0 uD α T -ϕ ≤ T 0 |u| p ϕ 1/p T 0 D α T -ϕ p ϕ -p /p 1/p , ( 24
) T 0 vD β T -ϕ ≤ T 0 |v| q ϕ 1/q T 0 D β T -ϕ q ϕ -q /q 1/q , ( 25
)
where p, p , q, q satisfy p + p = pp and q + q = qq . If we set

I := T 0 |u| p ϕ, J := T 0 |v| q ϕ, A := T 0 D α T -ϕ p ϕ -p /p , B := T 0 D β T -ϕ q ϕ -q /q , (26)
then ( 22), ( 23), ( 24) and (25) lead to:

J ≤ I 1/p A 1/p and I ≤ J 1/q B 1/q .
As u 0 > 0, v 0 > 0 and

T 0 D α T -ϕ > 0 when ϕ(t) is selected as in (14), it is clear that J 1-1 pq ≤ A 1/p B 1/pq , (27)
I 1-1 pq ≤ A 1/qp B 1/q . (28)
Using ( 16), with λ ≥ max{αp -1, βq -1}, we have

A = C p ,α T 1-αp , B = C q ,β T 1-βq .
Consequently, (27) takes the form

J 1-1 pq ≤ K 1 T s 1 , (29)
where

s 1 = 1 -α - β p - 1 pq .
Thus J is bounded under the constraint s 1 ≤ 0 which is equivalent to

1 - 1 pq ≤ α + β p .
The equality (23) and the estimates (25) allow us to write

v 0 T 0 D β T -ϕ ≤ J 1/q B 1/q , ( 30
)
which in virtue of (30) becomes

v 0 ≤ K 2 J 1/q T r 1 , (31)
where

r 1 = -1/q.
Letting T → ∞, in (31) we obtain the contradiction 0 < v 0 ≤ 0. Similarly, an analysis could be done along I; it leads to 1 -

1 pq ≤ β + α q .
To obtain an estimation on the blow up time, we estimate J in (31) by ( 29) to obtain

v 0 ≤ K 2 (K 1 ) p pq-1 T ps 1 pq-1 +r 1 = C v T s
where s < 0. Whereupon a bound on the blow-up time is given by

T max ≤ T v = v 0 C v 1/s .
A similar bound on T max can be obtained in terms of u 0 if we use ( 22) and (28); it reads

T max ≤ T u = u 0 C u 1/ s.
Finally, T max ≤ Min{T u , T v }.

Now, we present estimates from below and above of the blowing up solutions to the system (FDS) via the solutions of systems (PFDS) and (ODS) for µ = 1, 1/2 and λ = 1, 1/2. First observe that system (PFDS) is equivalent to the nonlinear system of Volterra equations

                 u(t) = u 0 + µ Γ(α) t 0 |v(τ)| q (t-τ) 1-α dτ, t > 0, v(t) = v 0 + µ Γ(β) t 0 |u(τ)| p (t-τ) 1-β dτ, t > 0. (32)
Due to the positivity of initial conditions (1), the solution of the system (32) is always greater or equal than the solution of the system

                 ũ(t) = µ Γ(α) t 0 |ṽ(τ)| q (t-τ) 1-α dτ, t > 0, ṽ(t) = µ Γ(β) t 0 |ũ(τ)| p (t-τ) 1-β dτ, t > 0, (33)
with µ = 1 or µ = 1/2. Following [START_REF] Olmstead | Coupled Volterra equations with blow-up solutions[END_REF], growth rates near the blow-up time for blowing up solutions to the system (33) are available. They are stated in Theorem 3.2. The profiles of the components of the solution (ũ(t), ṽ(t)) for the system (33) with µ = 1 are given by

ũ(t) ∼ C 1,α,β (T α,β max -t) -δ , ṽ(t) ∼ C 2,α,β (T α,β max -t) -ξ , as t → T α,β max , ( 34
)
where T α,β max is the blow-up time and

C 1,α,β = Γ(pδ) Γ(pδ -β) q pq-1 Γ(qξ) Γ(qξ -α) 1 pq-1 , C 2,α,β = Γ(qξ) Γ(qξ -α) p pq-1 Γ(pδ) Γ(pδ -β) 1 pq-1 , with δ = α + qβ pq -1 , ξ = β + pα pq -1 . (35)
Alike, the profiles of u(t), v(t) for the system (33) with µ = 1/2 are given by

u(t) ∼ C * 1,α,β (τ α,β max -t) -δ , v(t) ∼ C * 2,α,β (τ α,β max -t) -ξ , as t → τ α,β max , (36)
where

C * 1,α,β = Γ(pδ) 2Γ(pδ -β) q pq-1 Γ(qξ) Γ(qξ -α) 1 pq-1 , C * 2,α,β = Γ(qξ) 2Γ(qξ -α) p pq-1 Γ(pδ) Γ(pδ -β) 1 pq-1 ,
where δ and ξ are defined in (35).

Let us notice at this stage that the analysis in [START_REF] Olmstead | Coupled Volterra equations with blow-up solutions[END_REF] gives no information on T α,β max or τ α,β max . Now, let e α (t) := E α,1 (-t α ) where E α,β (t) is the Mittag-Leffler function of two parameters α and β.

Lemma 3.3. The solution (u, v) to the system (FDS) satisfies

u(t) = u 0 + t 0 e 1-α (t -τ)|v(τ)| q dτ, ( 37
)
v(t) = v 0 + t 0 e 1-β (t -τ)|u(τ)| p dτ. ( 38
)
Proof. Taking the Laplace transform of both sides of the first equation of the system (FDS) and using ( 9) and (13), we obtain

L{u(t) : s} = u 0 s + L{|v(t)| q : s} (s + s α ) ,
which, via the inverse Laplace transform, leads to:

u(t) = u 0 + L -1 L{|v(t)| q : s} (s + s α ) .
Notice that L{|v(t)| q : s} (s

+ s α ) = s -α (s 1-α + 1) L{|v(t)| q : s};
by means of the inverse Laplace transform of convolution of functions, equation ( 37) is obtained. Following the same steps, we obtain equation (38).

Lemma 3.4. For the system (FDS), the functions u (t) and v (t) satisfy

u (t) = |v(t)| q + t 0 e 1-α (t -τ)|v(τ)| q dτ, ( 39
)
v (t) = |u(t)| p + t 0 e 1-β (t -τ)|u(τ)| p dτ. ( 40
)
Proof. By virtue of the definition of the Caputo derivative, for 0 < α < 1, the equation corresponding to u(t) in system (FDS) takes the form

u (t) + J 1-α 0 + u (t) = |v(t)| q , t > 0, (41)
which can be written as

u (t) = (1 + J 1-α 0 + ) -1 |v(t)| q , t > 0. (42)
Since

(1 + J 1-α 0 + ) -1 = 1 + +∞ k=1 (-1) k J (1-α)k ,
we have

u (t) = |v(t)| q + +∞ k=1 (-1) k t (1-α)k-1 Γ((1 -α)k) |v(t)| q , = |v(t)| q + t 0 e 1-α (t -τ)|v(τ)| q dτ,
where relations (5) and ( 6) have been used. Let us mention that the same result can be obtained by using the Laplace transform technique: taking the Laplace transform of both sides of equation ( 41), it follows

L{u (t) : s} = s 1-α L{|v(t)| q : s} s 1-α + 1 . ( 43
)
There are different ways of writing the solution from equation ( 43); for example we can write

s 1-α L{|v(t)| q : s} s 1-α + 1 = s s -α L{|v(t)| q : s} s 1-α + 1 ,
and we obtain

u (t) = d dt e 1-α |v(t)| q .
Another way of writing equation ( 43) is

L{u (t)} = s -α s 1-α + 1 sL{|v(t)| q : s} -|v(0)| q + s 1-α s -α + 1 |v(0)| q ;
taking the inverse Laplace transform, we get

u (t) = t 0 d dt |v(t -τ)| q e α (τ)dτ + |v(0)| q e α (t).
Noting that e 1-α (0) = E 1-α,1 (0) = 1, equation (43) can be written as

L{u (t) : s} = s 1-α s 1-α + 1 -1 |v(s)| q + |v(s)| q , hence u (t) = t 0 |v(τ)| q e 1-α (t -τ)dτ + |v(t)| q .
The equation (40) can be obtained similarly.

Analysis of the results

It is a well known fact [START_REF] Pollard | The completely monotonic character of the Mittag-Leffler function E a (-x)[END_REF] that the Mittag-Leffler function E 1-α,1 (-t 1-α ) =: e 1-α (t) is completely monotone for t > 0, 0 < α < 1; so, in particular e 1-α (t) ≤ 0 for t > 0. Then equations (39) and (40) allow us to write

u (t) ≤ |v(t)| q , v (t) ≤ |u(t)| p , t > 0. (44)
The differential inequalities (44) lead to the estimates from above of the solution (u, v) for 0 < t < T max < +∞

u(t) ≤ C 1 T max -t -q+1 pq-1 , v(t) ≤ C 2 T max -t -p+1 pq-1 ,
where T max , C 1 and C 2 are given in (20). To obtain estimates from below for u and v, we notice that system (FDS) can be written as:

D α 0 + (u -u 0 )(t) = |v(t)| q -u (t), D β 0 + (v -v 0 )(t) = |u(t)| p -v (t), t > 0. ( 45
)
Using inequalities (44) we obtain

D α 0 + (u -u 0 )(t) ≥ 0, D β 0 + (v -v 0 )(t) ≥ 0. Lemma 3.5. For any u 0 > 0, v 0 > 0, it holds: u > 0, v > 0. Proof. Let u -= max(0, -u ), u + = max(0, u ), v -= max(0, -v ), v + = max(0, v ), so, u = u + -u -and v = v + -v -.
We are intended to show that u -= 0 and v -= 0.

From system (FDS) and using DD α 0 + = D α 0 + D (where D is the first order integer derivative) we have

u + D α 0 + u = qv q-1 v , (46)
v + D β 0 + v = pu p-1 u . (47)
Multiplying both sides of the equation ( 46) with u -and using the definitions of u + and u -, we get the equation

(u -)(u -) + (u -)D α 0 + (u -) = -qv q-1 v (u -). ( 48
)
Taking the integral over [0, T ] and using the inequality

T 0 w(t)D α 0 + w ≥ 0, w(0) = 0, 0 < α < 1 [5], the following inequality is obtained T 0 (u -)(u -) ≤ -q T 0 v q-1 v (u -), = -q T 0 v q-1 v + (u -) + q T 0 v q-1 v -(u -), ≤ q T 0 v q-1 v -(u -),
Using Hölder's inequality, we obtain

d dt 1 2 T 0 (u -) 2 ≤ qC(T ) T 0 (v -) 2 + qC(T )C * ( ) T 0 (u -) 2 , ≤ C 1 T 0 (v -) 2 + T 0 (u -) 2 , ( 49
)
where C 1 = max{qC(T ) , qC(T )C * ( )}, since v is locally bounded. Likewise, we obtain

d dt 1 2 T 0 (v -) 2 ≤ C 2 T 0 (v -) 2 + T 0 (u -) 2 . (50)
Adding ( 49) and (50), we obtain

d dt T 0 (u -) 2 + (v -) 2 ≤ C T 0 (v -) 2 + T 0 (u -) 2 , (51)
where C = max(2C 1 , 2C 2 ). Using Gronwall's inequality and the fact that u -(0) = v -(0) = 0, we conclude that u -= 0 and v -= 0; hence the result.

At this stage, the following configuration is obtained

D α 0 + (u -u 0 )(t) ≤ |v(t)| q , D β 0 + (v -v 0 )(t) ≤ |u(t)| p , u (t) ≤ |v(t)| q , v (t) ≤ |u(t)| p . (52)
Furthermore, the system (FDS) along with u > 0, v > 0 and

D α 0 + (u -u 0 )(t) ≥ 0, D β 0 + (v -v 0 )(t) ≥ 0 allow us to write either u (t) ≥ 1 2 |v(t)| q or D α 0 + (u -u 0 )(t) ≥ 1 2 |v(t)| q , and either v (t) ≥ 1 2 |u(t)| p or D β 0 + (u -u 0 )(t) ≥ 1 2 |u(t)| p .
The analysis of the above inequalities gives us the precise profile of the blowing up solutions of system (FDS).

Remark 3.6. The lack of estimates on the blow up time for (PFDS) prevent us to have bilateral estimates on the blow up time for the solutions of (FDS).

Numerical implementation

For the numerical treatment of the system (FDS), we will approximate the solution (u, v) through equations ( 37) and (38). The approximation of solutions to the system (PFDS) is obtained via its equivalent form given by the system of Volterra equations in (32).

The numerical approximations of the solutions to both systems (FDS) and (PFDS) require the approximations of the convolution integrals; this has been studied to a reasonable extent in the literature. One of the techniques is the convolution quadrature method in which the quadrature weights are determined by the Laplace transform of the convolution kernel and a linear multistep method (see [6] and references therein).

In this work, the backward Euler convolution quadratures are used for approximating the convolution integral Example 3:

For figure 3, we set p = 1.1, q = 1.8, α = 0.25, β = 0.4; the initial conditions are u 0 = 5, v 0 = 1. The profile of the solution (u, v), for the system (FDS) is given by the system (ODS) with λ = 1 and the system (PFDS) with µ = 1/2. As expected, our simulations show a dependance between the parameters (α, β, u 0 , v 0 , p, q) and the calculated blow up.

t 0 K

 0 (t -τ)g(τ)dτ. (53) Let > 0 be the time step of the discretization, thent 0 K(t -τ)g(τ)dτ = t-τ) L{K(τ) : s}ds g(τ)dτ = 1 2πi γ y(s, t)L{K(τ) : s}ds,where γ is a suitable path connecting -i∞ to +i∞ andy(s, t) = t 0 e s(t-τ) g(τ)dτ u 0 = 3, v 0 = 2.In fig(1) (a) the solution curves u(t) corresponding to the systems (FDS), (PFDS) and (ODS) are plotted. The dotted curve is the solution of the system (ODS) with λ = 1. As expected, the middle solid curve corresponds to the solution of the system (FDS). Finally, the dash followed by a dot curve is the solution of the system (PFDS) with µ = 1/2. Likewise, in fig(1) (b) the solution curves v(t) for the systems are plotted. It could be seen that the solution curves of the system (FDS) are between the solution curves of (ODS) and (PFDS).

Figure 1 :

 1 Figure 1: Solution curves for p = 1.5, q = 2, α = 0.75, β = 0.5, u 0 = 3, v 0 = 2 (a) Solution curves u(t) for the systems (FDS), (PFDS) and (ODS) (b) Solution curves v(t) for the systems (FDS), (PFDS) and (ODS).

Example 2 :

 2 In the second example, we take the parameters p = 1.1, q = 1.4, α = 0.5 = β; the initial conditions u 0 = 3 = v 0 . In fig(2) (a) the solution curves u(t) corresponding to the systems (FDS), (PFDS) and (ODS) are plotted. The dotted curve (which serve as lower bound for the solution curve of (FDS)) is the solution of the system (ODS) with λ = 1/2. The solid curve is the solution of the system (FDS) and the dash followed by dot curve is the solution of the system (PFDS) with µ = 1.In fig(2) (b) the corresponding solution curves v(t) for the systems are plotted. The profiles of the solution curves for the system (FDS) match the analysis of section 3.

Figure 2 :

 2 Figure 2: Solution curves for p = 1.1, q = 1.4, α = 0.5, β = 0.5, u 0 = 3, v 0 = 3 (a) Solution curves u(t) for the systems (FDS), (PFDS) and (ODS) (b) Solution curves v(t) for the systems (FDS), (PFDS) and (ODS).

is the solution of the initial value problem y = sy + g, y(0) = 0, 0 ≤ t ≤ T.

(54)

The backward Euler convolution quadrature is obtained as

where t n = n , and y n (s) stands for the approximation of y(s, t n ) reached by the backward Euler method applied to (54). Therefore, if δ(σ) is the generating polynomial for the backward Euler method that is δ(σ) = 1 -σ, then

the convolution quadrature weights q α j are given by the coefficients of

For more details on the convolution quadrature methods one is referred to [6] and references therein. Let for n ≥ 0 the approximation of u(t n ) be u n , then the discritization of the system (FDS) via equations (37), (38) by the backward Euler convolution quadrature method gives

for n = 1, 2, 3, ..., where the quadrature coefficients q α j , q β j are determined from equation (56) with K(t) = e 1-α and K(t) = e 1-β , respectively. Similarly, a numerical scheme for the system (PFDS) can be obtained.

Example 1:

For the figure 1 we set p = 1.5, q = 2, α = 0.75, β = 0.5; the initial conditions are


Ionut Mugurel

Andreica

Angela Andreica

An Algorithmic Perspective on Some Network Design, Construction and Analysis Problems

Keywords: network design, network analysis, network construction

Efficient network design, construction and analysis are important topics, considering the highly dynamic environment in which data communication occurs nowadays. In this paper we address several problems concerning these topics from an algorithmic point of view.

Introduction

Efficient network design and analysis techniques are strongly required in the highly dynamic environment in which data transfers occur nowadays. In this paper we focus on several algorithmic aspects of network design and analysis. We consider several constrained network design and (re)construction problems (Section 2), and several edge and vertex classification problems (Section 3). In Section 4 we discuss related work and we conclude.

Optimal Network Design and (Re)construction Problems

We consider N≥3 points in the plane, whose coordinates (x(i),y(i)) are not given. Every point i has a weight w(i)>0 (1≤i≤N). The arrangement of points is considered balanced if for every 3 distinct points A, B and C, we have w(A)+w(B)+w(C)=AreaTri(A,B,C) (AreaTri(A,B,C) is the area of the triangle ABC). We want to compute the coordinates of the points, such that the resulting arrangement is balanced. Moreover, we have w(N-1)=w(N). We first notice that there is no solution for N≥5. Thus, we only need to handle the cases N=3 and N=4. For N=3 we can choose: x(2)=y(2)=0, x(3)=2, y(3)=0, x(1)=0, y(1)=w( 1)+w [START_REF] Lin | The Most Vital Edges in the Minimum Spanning Tree Problem[END_REF]+w [START_REF] Ahuja | Network Flows: Theory, Algorithms, and Applications[END_REF]. The triangle has a right angle at (0,0). The two perpendicular sides have lengths C 1 =w(1)+w( 2)+w [START_REF] Ahuja | Network Flows: Theory, Algorithms, and Applications[END_REF] and C 2 =2. Thus, the triangle's area is (C 1 •C 2 )/2=w(1)+w( 2)+w [START_REF] Ahuja | Network Flows: Theory, Algorithms, and Applications[END_REF]. For N=4 we have 2 subcases. In the first subcase we have w( 1)=w [START_REF] Lin | The Most Vital Edges in the Minimum Spanning Tree Problem[END_REF] (apart from w(3)=w( 4)). In this subcase we will choose the coordinates of the points 3 and 4 at (0,0) and (2,0). The coordinates of point 1 will be (0, w(1)+w( 3)+w( 4)), and those of point 2 will be (2•(w(1)+w( 2)+w( 4))/(w(1)+w( 3)+w( 4)), w(1)+w( 3)+w( 4)). The 4 points form a trapezoid whose parallel sides are the segments connecting the points 1-2 and 3-4. Triangles 1-3-4 and 2-3-4 have heights of (w(1)+w( 3)+w( 4)), perpendicular on sides of length 2. Triangles 1-2-3 and 1-2-4 have heights equal to (w(1)+w( 3)+w( 4)) and perpendicular on a side with length 2•(w( 1)+w( 2)+w( 4))/(w(1)+w( 3)+w( 4)); thus, their area is w(1)+w( 2)+w( 4)=w( 1)+w( 2)+w [START_REF] Ahuja | Network Flows: Theory, Algorithms, and Applications[END_REF]. The second subcase for N=4 occurs when w(1)+w(2)=2•w [START_REF] Mccabe | Network Analysis, Architecture, and Design[END_REF]. In this subcase, the coordinates of the points 3 and 4 will be (0,0) and (2,0), and those of the points 1 and 2 will be (1, (w(1)+w( 3)+w( 4))) and (1, -(w( 2)+w( 3)+w( 4))). We notice that the 2 subcases for N=4 correspond to the situations in which the points 1 and 2 are on the same side or on opposite sides of the segment 3-4. We notice that for N≥5 (or N=4 and weights which do not satisfy the conditions of the two subcases) we cannot have a solution, because: if we have 4 points, then either the segment 1-2 is parallel to the segment 3-4, or the segment 1-2 intersects the segment 3-4 at its midpoint. For N≥5, one of these 2 situations should occur for the pairs [START_REF] Young | Faster Parametric Shortest Path and Minimum Balance Algorithms[END_REF][START_REF] Lin | The Most Vital Edges in the Minimum Spanning Tree Problem[END_REF], [START_REF] Young | Faster Parametric Shortest Path and Minimum Balance Algorithms[END_REF][START_REF] Ahuja | Network Flows: Theory, Algorithms, and Applications[END_REF] and (2,3) (relative to the segment 4-5, respectively to the midpoint of the segment 4-5). This indicates that either the area of the triangle 1-2-3 is 0 (but w(1)+w( 2)+w(3)>0), or that one of the points 1, 2 or 3 (let this point be P) is the midpoint of the segment 4-5, in which case the area of the triangle 4-5-P would be 0 (but w(4)+w( 5)+w(P)>0).

The 2 nd problem considers the reconstruction of a polygon with a known number N of vertices, whose edges may self-intersect. Let's consider that vertex i of the polygon has coordinates (x(i), y(i)) (1≤i≤N). On the support line of every edge (i, i+1) of the polygon (where i+1=1, for i=N), a point p(i) is chosen at distance t(i)•D(i,i+1) away from vertex i and towards vertex i+1. D(i,i+1) is the length of the edge (i,i+1) and t(i)≠0. The coordinates of the point p(i) are (xp(i),yp(i)). We want to reconstruct the polygon from the N points p(i). The coordinate xp(i) (yp(i)) of a point p(i) can be written as follows, based on the coordinates of the vertices i and i+1:

xp(i)= x(i)+t(i)•(x(i+1)-x(i))=(1-t(i))•x(i)+t(i)•x(i+1) (yp(i)=(1-t(i))•y(i)+ t(i)•y(i+1)).
Based on these equations we will write the coordinate x (y) of every vertex as a linear function of the coordinate x(1) (y(1)). We have x(i)=ax(i)•x(1)+bx(i) (y(i)=ay(i)•y(1)+by(i)). For i=1 we have ax(1)=1, bx(1)=0, ay(1)=1, and by(1)=0. For 2≤i≤N+1, we will obtain x(i) (y(i)) from the equation of point p(i-1): x(i)=(xp(i-1)-( 1

t(i-1))•x(i-1))/t(i-1)=(t(i-1)-1)/t(i-1)•x(i-1)+xp(i-1)/t(i-1) (the equation for y(i) is similar; we just replace xp(i-1) by yp(i-1), and x(i-1) by y(i-1)). From here we obtain: ax(i)=(t(i-1)-1)/t(i-1)•ax(i-1), bx(i)=(t(i-1)-1)/t(i-1)•bx(i-1)+xp(i-1)/t(i-1) (the equation for y(i)

is similar -we just replace the character "x" by "y" in the previous equation). After computing the values for the vertex N+1 (which is, in fact, vertex 1), we will write the equations: ax(N+1)•x( 1

)+bx(N+1)=x(1) => x(1)•(ax(N+1)-1)=-bx(N+1) (similarly, y(1)•(ay(N+1)-1)=-by(N+1)).

We will determine the coordinates x and y of the polygon vertices independently. We will present the method only for the coordinate x, because it is identical for the coordinate y (we just replace x by y). If (ax(N+1)-1)≠0), then we obtain x(1)=-bx(N+1)/(ax(N+1)-1), and then we will compute the coordinates of the other vertices: 1)+bx(i). If (ax(N+1)-1)=0, then we have 2 subcases. In the first subcase, bx(N+1)=0; thus, we can choose any value for x(1), afterwards obtaining the x coordinates of the other vertices. In the second subcase, bx(N+1)≠0 and the problem admits no solution. If we allow t(q)=0 (i.e. xp(q)=x(q) and yp(q)=y(q)), then we split the points q into maximal chains [i..j], such that t(j)=t(i-1)=0 and t(i), t(i+1), …, t(j-1) are non-zero (with q+1=1, for q=N, and q-1=N, for q=1). For every chain [i..j], we write x(k) (y(k)) as a linear function fx(k) (fy(k)) of x(i) (y(i)) (k belongs to [i..j]), like we did before, when we considered linear functions of x(1) (y(1)). Since x(j)=xp(j) (y(j)=yp(j)), we compute x(i) (y(i)) from fx(j)=xp(j) (fy(j)=yp(j)) and, based on them, all the other values x(k) and y(k) of the points k of the chain [i..j]. The time complexity of the algorithm is O(N) in any case.

x(i)=ax(i)•x(
In the 3 rd problem we consider a triangle ABC, in which the length of the side AB (denoted by lc), the length of the side AC (denoted by lb) and the length of the median AM (denoted by lm) are given. AM is the segment which connects the vertex A of the triangle ABC to the midpoint of the segment BC. We want to compute the coordinates of the 3 vertices. We will binary search the length of the segment BC. Let's assume that we considered this length to be 2•a. We will test if this length is valid. We will consider that the coordinates of the vertex B are (0,0), and those of vertex C are (2•a, 0). This way, the midpoint of the segment BC is at coordinates (a,0). If lm+a<lc or lm+a<lb, then the length 2•a is too small. Otherwise, we will binary search, between 0 and π, the angle AMB (determined by the points A, M and B). Let's assume that, during the binary search, we selected an angle alpha. We compute the corresponding length of the side AB, as being c'=sqrt(lm 2 +a 2 -2•lm•a•cos(alpha)) (we denote by sqrt(x) the square root of x). If c'<lc, then alpha is too small and we will consider a larger angle next; otherwise, we will consider a smaller angle next. The binary search of the angle will end when the length of the search interval becomes smaller than a constant ueps. After finding the angle alpha for which the length of the side AB is equal to lc, we will compute the length of the side AC corresponding to this angle: b'=sqrt(lm 2 +a 2 -2•lm•a•cos(π-alpha)). If b'<lb, then the length 2•a chosen within the binary search for the length of the side BC is too small and we will search for a larger value; otherwise, we will search for a smaller value. We finish the binary search when the length of the search interval becomes smaller than a predefined constant leps. The algorithm contains 2 nested binary searches and its complexity is O(log(LMAX)•log(UMAX)); LMAX is the maximum possible length of BC and UMAX=π. Alternatively, we can solve the system of equations (where we used cos(π-alpha)=-cos(alpha

)): (1) lc 2 =lm 2 +a 2 -2•lm•a•cos(alpha) ; (2) lb 2 =lm 2 +a 2 +2•lm•a•cos(alpha).
By adding (1) and (2) we obtain a=sqrt((lc 2 +lb 2 )/2-lm 2 ). After replacing the value of a in (1), we can also compute cos(alpha).

Classifying Edges and Vertices Relative to Matchings

We consider a bipartite multigraph with n vertices overall on its left and right parts. Every edge (u,v) (u on the left side and v on the right side) has a cost c(u,v)≥0. We want to classify every edge (vertex) of the multigraph in one of the following 3 categories: 1) it belongs to every minimum cost maximum matching ; 2) it belongs to at least one minimum cost maximum matching (but not to all of them) ; 3) it belongs to no minimum cost maximum matching. We start by computing a minimum cost maximum matching M, in O(n 3 ) time. Obviously, the edges (vertices) from category 1 must be among the edges (endpoints of the edges) of M. Then, we will construct a directed multigraph G 1 , as follows. We direct every edge (u,v) which belongs to M from v to u and we assign it a cost c'(v->u)=-c(u,v). We direct every edge (u,v) which does not belong to M from u to v and we assign it a cost c'(u->v)=c(u,v). Then, we add an extra vertex S and all the zero-cost directed edges (S->u) and (v->S) (u and v on the left side), where u is a vertex such that no edges adjacent to it belong to M, and v is a vertex such that there is an edge adjacent to it in M. After constructing this multigraph, we will compute the set E 1 of edges which belong to at least one zero-cost cycle (considering the costs c'), i.e. a cycle whose sum of edge costs is 0. This test can be performed in O(n 3 ) time using a standard method [START_REF] Young | Faster Parametric Shortest Path and Minimum Balance Algorithms[END_REF] for multigraphs with no negative cost cycles (as is our case), i.e. without cycles whose sum of edge costs is negative: we make all the edge costs nonnegative and we consider the directed multigraph G' with the zero cost edges only; the edges of G' which connect two vertices from the same strongly connected component of G' belong to a zero-cost cycle in G 1 . Then, we construct another directed multigraph G 2 , starting from the initial multigraph. We direct the edges and assign them costs just like before. This time we will add a vertex T instead of S, together with all the zero-cost directed edges (v->T) and (T->u) (u and v on the right side), where v is a vertex which has no edge adjacent to it belonging to M, and u is a vertex which has one edge adjacent to it belonging to M. We then compute the set of edges E 2 which belong to at least one zero-cost cycle in G 2 . Let E 1,a be the set of edges (S->v) from E 1 and E 1,b the set of edges (v->S) from E 1 , and let E 2,a be the set of edges (v->T) from E 2 and E 2,b the set of edges (T->v) from E 2 . After this, we set

E 1 =(E 1 \E 1,a )\E

E E ∪

. The edges from category 3 are the other edges of the multigraph (which do not belong to ) ( 21

E E M ∪ ∪
). Let Q(S) be the vertex set containing all the endpoints of the edges in the set of edges S. The set of vertices from category

1 is (Q(M)\Q(E 1,b ))\Q(E 2,b ), from category 2 is } , { \ )) ( ) ( ) ( ) ( ( , 2 , 2 , 1 , 1 T S E Q E Q E Q E Q b a b a ∪ ∪ ∪
; the remaining vertices belong to category 3. Thus, in O(n 3 ) time, we can classify all the edges (vertices) of the multigraph in the 3 categories.

We will now consider a similar problem, but for bipartite multigraphs without edge costs. We want to classify the edges (vertices) of the multigraph as: 1) belonging to every maximum matching ; 2) belonging to at least one maximum matching (but not all of them) ; 3) belonging to no maximum matching. As before, we start by computing a maximum matching M, in O(n 2.5 ) time. Then, we construct the multigraphs G 1 and G 2 like before, except that we do not assign costs to the edges. We compute the strongly connected components of G 1 and G 2 (in O(n 2 ) time) and then we compute the sets of edges E 1 and E 2 , consisting of the edges between two vertices u and v which belong to the same strongly connected component in G 1 and, respectively, G 2 . We also compute QV 1 (QV 2 ) as the set of vertices u such that the edge (u,S) ((u,T)) or (S,u) ((T,u)) exists in E 1 (E 2 ). Then, we remove the edges (S,*) and (*,S) from E 1 , the edges (T,*) and (*,T) from E 2 , and then we interpret the edges in E 1 and E 2 as being undirected. Just like before, the set of edges in category 1 is

) ( \ 2 1 E E M ∪ , the set of edges in category 2 is ) ( 2 1 E E ∪
and the edges in category 3 are those edges which do not belong to ) ( 21

E E M ∪ ∪
. The set of vertices from category 1 is (Q(M)\QV 1 )\ QV 2 and from category 2 is ) (

2 1

QV QV ∪

; the remaining vertices are in category 3. The algorithm is dominated by the computation of a maximum matching, and its complexity is O(n 2.5 ).

A related problem (whose solution was mentioned to us by C. Negruseri) asks for classifying the edges of a connected multigraph with N vertices and M edges as follows: 1) belonging to every minimum spanning tree (MST) ; 2) belonging to at least one MST (but not all) ; 3) belonging to no MST. We will sort the edges in ascending order, according to their cost, and we will then consider groups of edges with equal costs, in increasing order of the costs. We will use the disjoint sets data structure, which provides the functions Find(x) and Union(x,y), for maintaining the connected components of the graph's vertices. Let's assume that we reached the group of all edges with cost C. Let QV(C) be the set of endpoints of the edges in the group. We will construct a multigraph G(C) whose set of vertices is the set {Find(x)|x ∈ QV(C)}; every edge (x,y) from the current group induces an edge (Find(x), Find(y)) (possibly a loop, if Find(x)=Find(y)) in G(C). Find(x) denotes the representative of the connected component containing node x, in the graph G'(C)=the graph containing all the N vertices and all the edges with cost strictly smaller than C. The edges (x,y) inducing critical edges (bridges) in G(C) belong to category 1; the edges inducing non-bridge edges which connect two different nodes in G(C) belong to category 2; the edges inducing loops (i.e. those edges which connect a node to itself) in G(C) belong to category 3. After considering a group of edges with cost C, for every edge (x,y) in the group, we set that x is in the same connected component as y (we call Union(x,y)). overall for computing the critical edges (bridges) for all the groups of equal-cost edges. If we drop the costs of the edges, then every bridge (critical edge) of the multigraph belongs to every spanning tree and every other edge belongs to at least one spanning tree (i.e. for every edge (i,j) there is a spanning tree ST which contains it).

We will now consider a directed flow network with designated source and sink vertices (and possibly parallel edges) and we want to classify the network edges as being upward critical or not (an edge (u->v) is upward critical if by increasing its capacity by 1 unit and leaving the capacities of the other edges the same, the maximum flow in the network increases). First, we compute a maximum flow from the source(s) to the sink(s). Thus, every directed edge (u->v) has a flow f(u->v) on it, such that f(u->v)≤cap(u->v) (where cap(u->v) is the capacity of the edge). We perform a breadth-first search (BFS) starting from the source(s) and considering only those edges (u->v) with f(u->v)<cap(u->v). We mark the reachable vertices i by reachable 1 (i)=true (the other vertices j have reachable 1 (j)=false). Then, we perform a BFS from the sink(s), considering the reverse direction of the edges. before, we only consider those edges (u->v) with f(u->v)<cap(u->v). We mark the reachable vertices i by reachable 2 (i)=true (the other vertices j have reachable 2 (j)=false). A directed edge (u->v) is upward critical if f(u->v)=cap(u->v) and reachable 1 (u)=true and reachable 1 (v)=false and reachable 2 (v)=true and reachable 2 (u)=false. If the multigraph is undirected then: (1) we compute a maximum flow ; (2) we consider only the edges (u,v) with f(u,v)<cap(u,v) when performing the two breadth-first searches ; (3) for every edge (u,v), we consider its two possible orientations, (u->v) and (v->u), and use the criterion described above; (u,v) is upward critical if it is classified as such for one of the two orientations.

  Computing the critical edges (bridges) of a multigraph with n C vertices and m C edges can be done in O(n C +m C ) time (since we consider only those vertices which are adjacent to at least one of the m C edges => n C =O(m C )). The different stages of the algorithm have the following complexities: O(M•log(M)) for sorting the edges (or O(M+CMAX) if the edge costs are integer and are bounded by a small integer value CMAX); O(M•log*(N)) for all the Find(*) and Union(*,*) calls; O(M)

  1,b and E 2 =(E 2 \E 2,a )\E 2,b and we interpret the edges in E 1 and E 2 as being undirected. The edges from category 1 are those from the set

	M	\	(	1 E	∪	2 E	)	; those from category 2 form the set	(	1	2	)

Conclusions and Related Work

In this paper we addressed several network design and analysis problems, for which we presented novel algorithmic solutions. Similar or related problems were considered in [2,4]. Some of our solutions are based on algorithmic techniques presented in [1,3].


Ashwin Rao
email: ashwin.rao@inria.fr

Arnaud Legout
email: arnaud.legout@inria.fr

Walid Dabbous
email: walid.dabbous@inria.fr

BitTorrent Experiments on Testbeds: A Study of the Impact of Network Latencies

In this paper, we study the impact of network latency on the time required to download a file distributed using BitTorrent. This study is essential to understand if testbeds can be used for experimental evaluation of BitTorrent. We observe that the network latency has a marginal impact on the time required to download a file; hence, BitTorrent experiments can performed on testbeds.

I. INTRODUCTION

Testbeds such as PlanetLab and Grid5000 are widely used to study the performance of communication protocols and networking applications. One commonly used practice while performing experiments on such testbeds is to run multiple instances of the application being studied on the same machine. However, one primary shortcoming of this approach is the absence of any network latency between the instances of the application running on the same machine. Further, in experiments involving more than one machine, the latency between the machines present in the same local area network (LAN) is negligible. In this paper we study the impact of network latency on the outcome of experiments that are performed on testbeds to evaluate the performance of BitTorrent.

The BitTorrent Protocol internally uses the Transmission Control Protocol (TCP) while distributing the content [START_REF] Cohen | The BitTorrent Protocol Specification[END_REF]. The steady-state throughput of TCP is function of the roundtrip time (RTT) [START_REF] Mathis | The Macroscopic Behavior of the TCP Congestion Avoidance Algorithm[END_REF]. Further, the slow start and congestion avoidance phase of TCP introduce a ramp up period which is required to attain a throughput equal to the minimum of the network throughput and the rate at which the application is sending data. This ramp up period is a function of the RTT and the rate at which the data is being uploaded. BitTorrent allows the users to limit the rate at which data is uploaded; as the time duration of an upload by a peer is in the order of seconds, we believe that the time required to transfer pieces of a file is not affected by such variations in the TCP throughput. Our experiments show that the RTT (and hence the latency) between the peers in the torrent has a marginal impact (less than 15%) on the time required to download a file.

The details of the methodology and the tools used are presented in Section II. We initially assume the latency between any two peers in the torrent to be the same (homogeneous latency); the impact of homogeneous latency on the time required to download a file are presented in Section III. The results without this assumption are presented in Section IV, followed by the conclusions in Section V.

II. METHODOLOGY

In this paper we use the terminology used by the BitTorrent community. A torrent, also known as a BitTorrent session or a swarm, consists of a set of peers that are interested in having a copy of the given content. A peer in a torrent can be in two states: the leecher state when it is downloading the contents, and the seed state when it has a copy of the content being distributed. A tracker is a server that keeps track of the peers present in the torrent.

A. Experiment Scenarios

We consider a torrent consisting of one tracker and a finite number of peers; a few of these peers are seeds, while the rest of the peers are leechers. We assume that the peers remain in the torrent until all the leechers have finished downloading the file.

The metric used to study the impact of the network latency between the peers is the download completion time, the time required to the download the file distributed using BitTorrent. We use the following network topologies to evaluate the impact of latency on download completion time of a file.

1) Homogeneous Latency. The latency between any two peers in the torrent is the same in this network topology. This topology provides an upper bound on the download completion time when the maximum round trip time between the peers in a torrent is known. Further, this setting was used to give an insight on the threshold of the latency between the peers beyond which the latency affects the download completion time. 2) Heterogeneous Latency. The peers are grouped together to abstract Autonomous Systems (AS). We assume that the latency between any two peers in a given AS is the same and that all ASes are fully meshed. Further, we assume that the inter-AS latency is greater than the intra-AS latency; we also assume symmetric latency in the upload and download links within an AS and between ASes. All the experiments were performed in a private torrent consisting of one tracker, one initial seed (henceforth called as the seed), and 300 leechers; these experiments were carried out on the Grid'5000 experimental testbed [3]. A 50 MB file was distributed in this torrent where the upload rates of the leechers and the seed was varied from 10 kB/s to 100 kB/s. As shown in Figure 1, four machines with Linux as their operating system were used to run the instances of the peers

Node

Node Node

Fig. 1: Topology of the peers in the machines used for the experiment. One machine for the tracker and the initial seed, and three machines each with one hundred leechers. Fig. 3: CDF of TCP Payload Length. The maximum payload length with TSO enabled is much greater than that observed in the WIDE Backbone. Disabling TSO ensures that the maximum payload length is similar to that observed in the WIDE backbone.

in the torrent; one machine was used for the tracker and the seed, and each of the other three machines ran 100 instances of the leechers. A pair of peers in the torrent, including the tracker, the seed, and the leechers, communicate with each other using either the loopback interface or the ethernet interface. The latency between the peers was varied from 0 ms to 500 ms using the Netem module [START_REF] Hemminger | Network emulation with NetEm[END_REF]. This latency represents the one way delay observed by a packet, hence the round-trip time between any two peers is at least twice the latency mentioned. Homogeneous latency was achieved by adding the same latency on the loopback and the ethernet interface. Similarly, heterogeneous latency was achieved by adding a latency on the ethernet interface of a given machine that is greater than the latency added on the loopback interface of the same machine.

The following torrent configurations were used to vary the upload rate of the peers:

1) Seed and Leechers Slow. In this setting, the upload rate of the peers was limited to 10 kB/s and 20 kB/s. 2) Seed and Leechers Fast. In this setting, the upload rate of the peers was limited to 50 kB/s and 100 kB/s. 3) Seed Fast and Leechers Slow. In this setting, the upload rate of the initial seed was limited to 50 kB/s while the upload rate of the leechers was limited to 20 kB/s.

B. Testbed Configuration

The Netem module buffers the TCP frames which are in flight for a time period equal to the latency being emulated. A buffer size of 100000 frames was used in each machine to support up to 1000 frames of each peer to be in flight. To ensure that the machines are capable of running 100 peers uploading at 100 kB/s without affecting the added latency, we varied the number of leechers running on a machine from 4 to 100. The TCP_INFO option for the getsockopt method of the socket library was used to sample the RTT estimated by TCP each time a send system call was issued on a socket. Figure 2 shows the average RTT estimate with the error bars representing the 95th and 5th percentiles of all the peers in five iterations. We observe that the number of leechers running a given machine has a marginal impact on the average RTT estimated by TCP when each of the leechers has its maximum upload rate limited to 100 kB/s.

The Maximum Transmission Unit (MTU) for the loopback interface is typically greater than the MTU for other network interfaces such as Wi-Fi and ethernet. To avoid the impact of large frames being exchanged between the peers we set the MTU of the loopback interface to 1500 bytes (the default value set for the ethernet interface). Figure 3 (top plot) shows that despite setting the MTU to 1500 bytes, a significant number of frames have a size greater than the MTU. Further, we observe that increasing the latency between the peers results in a significant number of TCP segments having a large payload length; TCP segments with large payload lengths were also observed in frames being sent over the ethernet interface. This increase in payload lengths is due to a feature called TCP Segmentation Offloading (TSO) [START_REF] Mogul | Tcp offload is a dumb idea whose time has come[END_REF], which is enabled by default in the 2.6 series (the current series) of the Linux kernel. TSO enables the host machine to offload some of the TCP/IP implementation, such as segmentation, and calculation of IP checksum, to the network device. Further, to enhance the throughput, TSO supports the exchange of data in frames of sizes that can be greater that the underlying MTU size. The increase in the frame size can result in significant improvement in throughput; however, this improvement depends on various factors such as CPU processing power and the amount of data being transferred [START_REF] Freimuth | Server Network Scalability and TCP offload[END_REF]. Figure 3 (middle plot) shows the packet lengths obtained from the publicly available traces of the Internet traffic in the WIDE backbone [7]. The values presented are from the sample taken on the WIDE backbone on November 29, 2009. As hardware support on all the devices in the communication link is essential for handling large segments, the graph shows that the devices on most of links (including end hosts or intermediate nodes) in the Internet do not support TSO. In this paper we restrict ourselves to study the impact of the latency, hence TSO was disabled in the subsequent experiments. Figure 3 (bottom plot) shows that the maximum payload length of the frames is similar to that observed in the Internet when the MTU is set to 1500 bytes and TSO is disabled. The impact of the latency and heterogeneous latency between the peers are presented in the subsequent sections. The plots presented are the outcome of 10 iterations.

III. HOMOGENEOUS LATENCY

We now present the impact of homogeneous network latency between any two peers in a torrent on the download completion time of a file.

A. Presentation of Results

When the maximum upload rate of the seed and leechers is limited to 10 kB/s, Figure 4 shows that an RTT greater than 200 ms results in the average download completion time increasing by at most 15%. Further, when the upload rate of all the peers is limited to 20 kB/s we observe that the download completion time is not a monotonously increasing function of the RTT; peers having an RTT of 1000 ms have a lower download completion time compared to peers having an RTT of 300 ms. However, when the maximum upload rate of the seed is increased to 50 kB/s while that of the leechers is limited to 20 kB/s, we observe that an RTT of even 1000 ms between the peers has a marginal impact on the download completion time of the file. A similar observation is made when the limit on maximum upload rates of all the peers is increased to 50 kB/s; this is also true when the upload rates are limited to 100 kB/s. Figure 5, Figure 6, and Figure 7 show the distribution of the time between successive send system calls while uploading the blocks (and pieces) of the file being distributed. Figure 5a, Figure 6a, and Figure 7 show that the distribution at the seed is similar for an RTT of 0 ms, 400 ms, and 1000 ms; further, from Figure 5b and Figure 6b, we observe that the distribution at the leechers for upload rates of 50 kB/s and 100 kB/s is the similar when the RTT is 0 ms and 400 ms. This shows that the RTT has a marginal impact on the upload process at the peers when their maximum upload rates are limited to 50 kB/s and 100 kB/s. However, for the upload rates of 10 kB/s and 20 kB/s we observe that the leechers tend to have a smaller time between successive send system calls when the RTT is 400 ms (or, 1000 ms) compared to an RTT of 0 ms. As the peers can simultaneously upload to many peers in parallel, the low inter-arrival time implies that the upload capacity is being utilized to upload data to a smaller number of peers.

B. Discussion of Results

For upload rates of 10 kB/s and 20 kB/s, Figure 4, Figure 5, and Figure 6, show that when the time between successive send system calls is less than the RTT, the RTT does not have an impact on the download completion time. Further, we observe that the ramp-up period required to attain a throughput equal to the upload rate of 50 kB/s (or 100 kB/s) does not have an impact on the download completion time. However, we currently do not have an accurate reason for the nonmonotonous increase in the download completion time for upload rates of 20 kB/s.

The above results show that network latency has a negligible impact on the download completion time of a file if the peers are fast (capable of uploading at high rates such as 50 kB/s). However, we observe that the latency affects the download completion time when the peers are slow (upload rates are less than or equal to 20 kB/s). Further, we observe that a single fast seed is capable of mitigating the impact of network latency on a torrent consisting of slow leechers.

IV. HETEROGENEOUS LATENCY

The latency between two peers in a given AS, is usually less than the latency between a peer from the given AS and another peer present in an adjacent AS. We emulated ASes by ensuring that the latency on the loopback interface is less than that on the ethernet interface on each of the machines used; hence, two peers running on a given machine have an RTT less than the RTT between a peer running on the given machine and another peer running on another machine. Further, we assume all the ASes to be fully meshed.

A. Abstraction of ASes

As in the case of homogeneous latency, we consider a torrent consisting of three hundred leechers, one tracker, and one initial seed; we emulated three ASes with one hundred leechers each, while the seed and the tracker were placed in the fourth AS. The four ASes used in these experiments were chosen from a set five ASes; the latency values added on the ingress and egress of the loopback and ethernet device of the machines while emulating these five ASes are given in Table I. III: RTT between the initial seed and the other peers (except the tracker) in the torrent. RTT between the seed in AS ′

1 and a peer in AS 1 is 20 ms.

We now show how Figure 1 and Table I can be used to find the RTT between a pair of peers.

As a machine is used to emulate an AS, a peer in AS 1 uses the ethernet interface to communicate with a peer in AS 2 , the RTT between this pair of peer is therefore 40 ms (5+15+15+5); as a peer in AS 1 uses the loopback interface to communicate with another peer in AS 1 , the RTT between this pair of peers is 8 ms (2+2+2+2). Table II gives the RTT values between all such pairs of peers that are initially leechers.

Similarly, from Figure 1, as the initial seed (henceforth called as the seed) and leechers are placed in different machines, the seed uses the ethernet interface of the machine to communicate with all the leechers present in the torrent. We use AS ′ i to denote that the seed and the tracker are placed in an AS with an inter-AS latency and intra-AS latency equal to that of AS i ; for example, AS ′ 1 implies that the seed and tracker are placed in an AS having the same latency values as AS 1 . Therefore, the RTT between the seed in AS ′ 1 and a leecher in AS 1 is 20 ms (5 + 5 + 5 + 5), while the RTT between the same seed and a leecher in AS 5 is 210 ms. The RTT between the seed and the peers are given in Table III.

From Table II, the RTT between a peer in either AS 1 , AS 2 , or AS 3 , and another peer in either AS 1 , AS 2 , or AS 3 , is less than 100 ms. Further, the RTT between a peer in AS 4 and another peer in AS 4 is 100 ms, while the RTT from this peer to any other peer is greater than 200 ms. Similarly the RTT from a peer in AS 5 to any other peer, irrespective of its AS, is greater than 200 ms.

B. Presentation and Discussion of Results

Figure 8 and Figure 9 show the impact of heterogeneous latency on the download completion time of a 50 MB file when the upload rate of the peers is limited to 20 kB/s and 50 kB/s respectively. The X-axis represents the AS of the leechers present in the torrent, and the Y-axis represents the download completion time in seconds; the error bars indicate the minimum and maximum download completion time of the leechers in 10 iterations.

Figure 8a shows the outcome of three experiments having the leechers placed in AS 1 , AS 2 , and AS 3 ; for the first experiment the seed was placed in AS ′ 1 , for the other two experiments, the seed was placed in AS ′ 2 and AS ′ 3 respectively. The figure shows that the RTT between the peers does not affect the download completion time. According to Table II, the RTT between any two peers in these experiments was less than 120 ms, hence, as in the case of homogeneous latency, we observe that an RTT of less than 120 ms does not affect the download completion time when the upload rates are limited to 20 kB/s. When some of the peers are present in an AS having a large RTT (AS 4 or AS 5 ), and the seed is also present in another AS with a large RTT (AS ′ 4 or AS ′ 5 ), then Figure 8b, and Figure 8c, show that the RTT affects the download completion time. However, we observe that the increase in average download completion time is not more than 15% of the average download completion time when all the peers in a torrent have an RTT less than 120 ms.

Figure 9 shows the impact of heterogeneous latency on the download completion time when the maximum upload rate of all the peers in the torrent is limited to 50 kB/s. We observe that an RTT of 400 ms, between the seed in AS ′ 5 and the leechers in AS 5 , does not have a significant impact on the download completion time of the file. These observations are in line with the observations in Section III.

Figure 4, Figure 8, and Figure 9, confirm that the topology of homogeneous latency provides an upper bound on the download completion time of a file when the maximum latency between any two peers in a torrent is known. Further, the observations made in Section III can be used in experiments where the latency between a pair of peers is heterogeneous.

V. CONCLUSION

The network latency between the peers has a marginal impact on the download completion time when the peers have their upload rates limited to high values such as 50 kB/s and 100 kB/s; our experiments show that the ramp-up period which is required to attain the throughput equal to these upload rates has a marginal impact on the download completion time. When the peers are slow (upload rates limited to values less than or equal to 20 kB/s) we observe that the download completion time is affected by the network latency; however, the increase in the average download completion time is not more than 15% of the average download completion time when there is no network latency between the peers. As the network latency has a marginal impact on the time required to download a file, BitTorrent experiments can be performed on testbeds without explicitly emulating latency between the peers in a torrent.

VI. ACKNOWLEDGMENT

Experiments presented in this paper were carried out using the Grid'5000 experimental testbed, being developed under the INRIA ALADDIN development action with support from CNRS, RENATER and several Universities as well as other funding bodies (see https://www.grid5000.fr).

Fig. 2 :

 2 Fig.2: Impact of the number of leechers running on a machine on the RTT estimate of TCP. Error bars indicating the 5 th and 95 th percentile of RTT estimated by TCP for all the peers in five iterations. Increasing the number of leechers running on a given machine has a marginal impact on the RTT estimated by TCP.

Fig. 4 :

 4 Fig. 4: Impact of latency on the Download Completion Time. RTT values mentioned are twice the latency added on a link using Netem. The latency increases the download completion time by at most 15%.

Fig. 5 :

 5 Fig.5: Distribution of the inter-arrival time of data messages when the RTT between the peers is 0 ms.

Fig. 6 :Fig. 7 :

 67 Fig.6: Distribution of the inter-arrival time of data messages when the RTT between the peers is 400 ms.

AS

Fig. 8 : 5 Fig. 9 :

 859 Fig. 8: Download completion time of a 50 MB file by leechers in a given AS when the maximum upload rate of all the peers is 20 kB/s.

TABLE I :

 I Latency values for the ingress and egress of the loopback and ethernet device while emulating an AS on a machine.

TABLE II :

 II RTT between a pair of leechers. RTT between a leecher in AS 1 and a leecher in AS 5 is 210 ms.

		AS 1	AS 2	AS 3	AS 4	AS 5
	AS ′ 1	20 ms	40 ms	60 ms	210 ms 210 ms
	AS ′ 2	40 ms	60 ms	80 ms	230 ms 230 ms
	AS ′ 3	60 ms	80 ms	100 ms 250 ms 250 ms
	AS ′ 4	210 ms 230 ms 250 ms 400 ms 400 ms
	AS ′ 5	210 ms 230 ms 250 ms 400 ms 400 ms

TABLE

  Download completion time for leechers present in AS 1 , AS 3 , and AS 4 . Having peers in an AS with large latency (AS 4 ) and the initial seed in an AS with large latency (AS ′ 4 ) affects the download completion time. Download completion time for leechers present in AS 1 , AS 3 , and AS 5 . Having peers in an AS with large latency (AS 5 ) and the initial seed in an AS with large latency (AS ′ 5 ) affects the download completion time.

	Download Completion	Time of 50 MB file	by a leecher (seconds)	2500 3000 3500 4000	AS 1 AS 2 AS 3 Seed in AS' 1 Seed in AS' 2 Seed in AS' 3
					AS of a leecher present in the torrent
	(a) Download completion time for leechers present in AS 1 ,
	AS 2 , and AS 3 .
	Download Completion	Time of 50 MB file	by a leecher (seconds)	2500 3000 3500 4000	AS 1	AS 3 AS 4 Seed in AS' 1 Seed in AS' 4
					AS of a leecher present in the torrent
	(b) 2500 3000 3500 4000 Download Completion Time of 50 MB file by a leecher (seconds)	AS 1	AS 3 Seed in AS' 1 Seed in AS' 5	AS 5
					AS of a leecher present in the torrent
	(c)