Frans Schalekamp

René Sitters
email: r.a.sitters@vu.nl

Suzanne Van Der Ster
email: suzanne.vander.ster@vu.nl

Leen Stougie
email: l.stougie@vu.nl

Víctor Verdugo
email: vverdugo@dim.uchile.cl

Anke Van Zuylen

Split scheduling with uniform setup times

Keywords: Scheduling, Job splitting, Setup times, Complexity theory, Approximation algorithms

come

Introduction

We consider a scheduling problem with setup times and job splitting. Given a set of identical parallel machines and a set of jobs with processing times, the goal of the scheduling problem is to schedule the jobs on the machines such that a given objective, for example the makespan or the sum of completion times, is minimised. With ordinary preemption, feasible schedules do not allow multiple machines to work on the same job simultaneously. In job splitting, this constraint is dropped. Without setup times, allowing job splitting makes many scheduling problems trivial: both for minimising makespan and for minimising total (weighted) completion time, an optimal schedule is obtained by splitting the processing time of each job equally over all machines, and processing the jobs in arbitrary order on each machine in case of makespan, and in (weighted) shortest processing time first ((W)SPT) order in case of total (weighted) completion time. See Xing and Zhang [START_REF] Xing | Parallel machine scheduling with splitting jobs[END_REF] for an overview of several classical scheduling problems which become polynomially solvable if job splitting is allowed.

In the presence of release times, minimising total completion time with ordinary preemption is NP-hard [START_REF] Du | Minimizing mean flow time with release time constraint[END_REF], whereas it is easy to see that if we allow job splitting, then splitting all jobs equally over all machines and applying the shortest remaining processing time first (SRPT) rule gives an optimal schedule.

Triviality disappears when setup times are present, i.e., when each machine requires a setup time before it can start processing the next job (part). During setup, a machine cannot process any job nor can it set up the processing of any other job (part). Problems for which the setup times are allowed to be sequence-dependent are usually NP-hard, as such problems tend to exhibit routinglike features. For example, the Hamiltonian path problem in a graph can be reduced to the problem of minimising the makespan on a single machine, where each job corresponds to a node in the graph, the processing times are 1, and the setup time between job i and j is 0 if the graph contains an edge between i and j, and 1 otherwise. However, as we will see, adding setup times leads to challenging algorithmic problems, already if the setup times are assumed to be job-, machine-, and sequence-independent.

We encountered such problems in studying disaster relief operations [START_REF] Van Der Ster | The allocation of scarce resources in disaster relief[END_REF]. For example in modelling flood relief operations, the machines are pumps and the jobs are locations to be drained. Or in the case of earthquake relief operations, the machines are teams of relief workers and the jobs are locations to be cleared. The setup is the time required to install the team on the new location. Although, in principle, these setup times consist partly of travel time, which is sequence-dependent, the travel time is negligible compared to the time required to equip the teams with instructions and tools for the new location. Hence, considering the setup times as being location-and sequence-independent was in this case an acceptable approximation of reality.

In this paper we concentrate on a basic scheduling problem and consider the variation where we allow job splitting with setup times that are job-, machine-, and sequence-independent, to which we will refer here as uniform setup times; i.e., we assume a uniform setup time s. There exists little literature on this type of scheduling problem. The problem of minimising makespan on parallel identical machines is in the standard scheduling notation of Graham et al. [START_REF] Graham | Optimization and approximation in deterministic sequencing and scheduling: a survey[END_REF] denoted as P ||C max (see Section 8 for an instruction on this notation). This problem P ||C max , but then with job splitting and setup times that are job-dependent, but sequenceand machine-independent, is considered by Xing and Zhang [START_REF] Xing | Parallel machine scheduling with splitting jobs[END_REF], and Chen, Ye and Zhang [START_REF] Chen | Lot-sizing scheduling with batch setup times[END_REF]. Chen et al. [START_REF] Chen | Lot-sizing scheduling with batch setup times[END_REF] mention that this problem is NP-hard in the strong sense, and only weakly NP-hard if the number of machines is assumed constant. Straightforward reductions from the 3-Partition and Subset Sum problem show that these hardness results continue to hold if setup times are uniform. Chen et al. provide a 5/3-approxi-mation algorithm for this problem and an FPTAS for the case of a fixed number of machines. A PTAS for the version of P ||C max with preemption and job-dependent, but sequence-and machine-independent setup times was given by Schuurman and Woeginger [START_REF] Schuurman | Preemptive scheduling with job-dependent setup times[END_REF]. It remains open whether a PTAS exists with job splitting rather than preemption, even if the setup times are uniform. See [START_REF] Liu | Minimizing total completion time subject to job release dates and preemption penalties[END_REF] and [START_REF] Potts | Integrating scheduling with batching and lot-sizing: A review of algorithms and complexity[END_REF] for a more extensive literature on problems with preemption and setup times.

Our problem is related to scheduling problems with malleable tasks. A malleable task may be scheduled on multiple machines, and a function f j (k) is given that denotes the processing speed if j is processed on k machines. If k machines process task j for L time, then f j (k)L units of task j are completed. What we call job splitting is referred to as malleable tasks with linear speedups, i.e., the processing time required on k machines is 1/k times the processing time required on a single machine. We remark that job splitting with setup times is not a special case of scheduling malleable tasks, because of the discontinuity caused by the setup times. We refer the reader to Drozdowski [START_REF] Drozdowski | Scheduling for Parallel Processing[END_REF] for an extensive overview of the literature on scheduling malleable tasks.

The main algorithmic result of our paper considers the job splitting variant of the problem of minimising the sum of completion times on identical machines, with uniform setup times: given a set of m identical machines, n jobs with processing times p 1 , . . . , p n , and a setup time s, the objective is to schedule the jobs on the machines to minimise total completion time ( C j ) (where the chosen objective is inspired by the disaster relief application). The version of this problem with ordinary preemption and fixed setup time s is solved by the Shortest Processing Time first rule (SPT); the option of preemption is not used by the optimum. However, the situation is much less straightforward for job splitting. If s is very large, then an optimal schedule minimises the contribution of the setup times to the objective, and a job will only be split over several machines if no other job is scheduled after the job on these machines. It is not hard to see that the jobs that are not split are scheduled in SPT order. If s is very small (say 0), then each job is split over all machines and the jobs are scheduled in SPT order. However, for other values of s, it appears to be a non-trivial problem to decide how to schedule the jobs, as splitting a job over multiple machines decreases the completion time of the job itself, but it increases the total load on the machines, and hence the completion times of later jobs.

Consider the following instance as an example. There are 3 machines and 6 jobs, numbered 1, 2, . . . , 6, with processing times 1, 2, 3, 5, 11, 12, respectively, and setting up a machine takes 1 time unit. One could con- sider filling up a schedule in round-robin style, assigning the jobs to machine 1, 2, 3, 1, 2, 3, respectively. This schedule is given in the Gantt chart in Figure 1(a). The schedule has objective value 49. By splitting job 6 over machines 1 and 3, instead of processing it on machine 3 only, we can lower the completion time of job 6, and this improves the objective value since there are no jobs scheduled after job 6. In fact, to get the best improvement in objective value, we make sure that both job parts of job 6 finish at the same time, see Figure 1(b). The objective value of the schedule is 45.

Splitting jobs early in the schedule, may increase the objective value, as (many) later jobs may experience delays. For example, if we choose to split job 2 over machines 2 and 3, we will cause delays for jobs 3 and 6, while improving the completion times of jobs 2 and 5. If we require that job parts of the same job end at the same time, we get the schedule pictured in Figure 1(c) with objective value 46. Finally, Figure 1(d) depicts the optimal schedule with objective value 40.

This example illustrates the inherent trade-off in this problem mentioned earlier: splitting jobs will decrease the completion times of some jobs, but it also may increase the completion times of other jobs.

In Section 3 we present a polynomial-time algorithm for the case in which there are two machines. The al-gorithm is based on a careful analysis of the structure of optimal solutions to this problem. Properties of optimal solutions that hold under any number of machines are presented in a preliminary section. Though a first guess might be that the problem would remain easy on any fixed number of machines, we will show by some examples in Section 4 that nice properties, which make the algorithm work for the 2-machine case, fail to hold for three machines already. The authors are split between thinking that we have encountered another instance of Lawler's "mystical power of twoness" [START_REF] Karel | The mystical power of twoness: in memoriam Eugene L. Lawler[END_REF], a phrase signifying the surprisingly common occurance that problems are easy when a problem parameter (here the number of machines) is two, but NP-hard when it is three, or that we just lacked the necessary flash of insight to find a polynomial-time algorithm. We present a constant-factor approximation algorithm for the general case with any number of machines in Section 5, and in Section 6 we give a polynomial-time approximation scheme for the case of a fixed number of machines. We leave the complexity of the problem (even for only three machines) as a tantalising open problem for the scheduling research community. We show in Section 7 that introducing weights for the jobs makes the problem NP-hard, already on 2 machines. We finish the paper by giving a table with the state of the art for other split scheduling problems with uniform setup times. We summarize whether they are known to be NP-hard or in P, and present the best known approximation ratios.

Preliminaries

An instance is given by m parallel identical machines and n jobs. Job j has processing time p j , for j = 1, . . . , n. Each job may be split into parts and multiple parts of the same job may be processed simultaneously. Before a machine can start processing a part of a job, a fixed setup time s is required. During setup of a job (part) the machine cannot simultaneously process or setup another job (part). The objective is to minimise the sum of the completion times of the jobs (total completion time), which is equivalent to minimising the average completion time.

Here we derive some properties of an optimal schedule, which are valid for any number of machines. Some additional properties for the special case of two machines, presented in Section 3, will lead us to a polynomial-time algorithm for this special case. We show in Section 4 that the additional properties that make the 2-machine case tractable do not hold for the case of three machines.

Claim Let σ be a feasible schedule with job completion times

C 1 ≤ C 2 ≤ • • • ≤ C n .
Let σ be obtained from σ by rescheduling the job parts on each machine in order 1, 2, . . . , n. Then C j ≤ C j for j = 1, . . . , n.

Proof Let q ij be the time that j is processed on machine i in σ and let C ij be the time that j finishes on machine i. Let y ij = s + q ij if q ij > 0 an let y ij = 0 otherwise. Fix some job j and machine i.

Let k = arg max{C ik | 1 ≤ k ≤ j}. Then C j ≥ C k ≥ C ik ≥ j h=1 y ih = C ij ,
where the first inequality is by assumption and the last one by the fact that all work on jobs smaller than or equal to j has been done on machine i at time C ik . Since C j ≥ C ij for any machine i on which j is scheduled the proof follows.

The claim above has several nice corollaries. First, note that if in an optimal schedule

C 1 ≤ C 2 ≤ • • • ≤ C n ,
then we maintain an optimal schedule with the same completion time for each job by scheduling the job parts on each machine in the order 1, 2, . . . , n. This allows to characterize an optimal schedule by a permutation of the jobs and the times that job j is processed on each machine i. The optimal schedule is then obtained by adding a setup time s for each non-zero job part and processing them in the order of the permutation on each machine. Consequently, in the optimal schedule obtained, each machine contains at most one part of each job.

We thus have the following lemma, which we will use throughout this work.

Lemma 1 There exists an optimal schedule such that each machine contains at most one part of each job.

In the sequel, given a schedule, we use M j to denote the set of machines on which parts of job j are processed. We will sometimes say that a machine processes job j, if it processes a part of job j.

Lemma 2 There exists an optimal schedule that satisfies the property of Lemma 1 such that on each machine the job parts are processed (started and completed) in SPT order of the corresponding jobs. Proof Among the optimal schedules that satisfy Lemma 1, we choose the schedule that minimises h p h C h . By the observations preceding Lemma 1, we may assume the jobs are numbered 1, . . . , n so that C 1 ≤ C 2 ≤ . . . ≤ C n , and each machine processes the job parts in the order given by the numbering of the jobs. Suppose by contradiction that there exist jobs j, k such that p j < p k and there exists some machine that processes job k before j, i.e., C k ≤ C j . Choose among such pairs of jobs j, k a pair that minimizes j -k. Note that any machine that processes both j and k must process k immediately before j, since if there is some job that is processed between them, then C k ≤ C ≤ C j , and either p > p j or p ≤ p j < p k , so either j, or , k should have been chosen instead of the pair j, k. We now show how to define a new optimal schedule for which h p h C h is strictly less than for the original schedule, thus contradicting the choice of our schedule.

Note that M j ∩ M k = ∅. We define a new schedule by rescheduling both jobs within the time slots these jobs occupy in the current schedule (including the slots for the setup times). First remove both jobs. Then consider the machines in M k one by one, starting with the machines in M k \M j and fill up the slots previously used by job k, until we have completely scheduled job j including the setup times. This is possible since p j < p k . We consider the remaining slots, which we note are single time intervals for each machine by our choice of j and k. We will show that they provide sufficient time for the processing and set up of job k, by showing that the combined number of setups for j and k does not increase.

Let M j and M k denote the sets of machines occupied by j and k, respectively, in the new schedule. We distinguish two cases. If job j cannot be rescheduled completely in the slots used by k in M k \ M j then we have

M k ⊆ M j . Together with M j ⊆ M k it follows that (M j ∩ M k ) ⊆ (M j ∩ M k ).
Hence, any machine containing both j and k in the new schedule did also contain both jobs in the old schedule, and therefore there are no extra setups on any machine needed. Now consider the case that job j is rescheduled completely in the slots used by k in M k \ M j . Then, after adding job k, the total number of setups needed for j and k does not increase since there is at most one machine of M k \ M j containing both jobs in the schedule, but none of the machines in M j ∩ M k is used by j in the new schedule.

We conclude that the remaining slots after scheduling job j provide sufficient room to feasibly schedule both the processing of job k and the required setups. Note that, if there is some machine on which the time allotted to k is at most s, then we can simply leave the machine idle for that time interval.

Let C denote the new completion times. We have

C j ≤ C k and C k ≤ max{C j , C k },
since in the new schedule j is processed only where job k was processed in the old schedule, and job k is processed in the new schedule only where either job j or job k was processed in the old schedule. For all other jobs, the completion time remains the same. Now, by assumption, we have that C j > C k , and hence C k ≤ C j . Therefore, the sum of completion times did not increase, and h p h C h < h p h C h , which contradicts the choice of the original schedule.

From now on, assume that jobs are numbered in SPT order, i.e., p 1 ≤ • • • ≤ p n . Given a schedule, we call a job balanced if it completes at the same time on all machines on which it is processed.

Lemma 3 There exists an optimal schedule that satisfies the properties of Lemma 1 and Lemma 2 in which all jobs are balanced.

Proof Consider an optimal schedule of the form of Lemma 1 and Lemma 2 with a minimum number of job parts. Let C j be the completion time of j in this schedule and define M j for this schedule as before. Consider the following linear program in which there is a variable x ij for all pairs i, j with i ∈ M j , indicating the amount of processing time of job j assigned to machine i:

min j C j s.t. i∈Mj x ij = p j ∀j = 1, . . . , n, k≤j: M k i (s + x ik ) ≤ C j ∀j = 1, . . . , n, ∀i ∈ M j , x ij ≥ 0, C j ≥ 0 ∀j = 1, . . . , n, ∀i ∈ M j .
Note that a schedule that satisfies Lemmas 1 and 2 gives a feasible solution to the LP, and on the other hand that any feasible solution to the LP gives a schedule with total completion time at most the objective value of the LP: if there exist some j and i ∈ M j such that x ij = 0, then the LP objective value is at least the total completion time of the corresponding schedule, as there is no need to set up for job j on machine i if x ij = 0. We know that a solution is a basic solution to this LP, only if the number of variables that are non-zero is at most the number of linearly independent tight constraints (not including the non-negativity constraints). By the minimality assumption on the optimal schedule, in any optimal solution to the LP all C j and x ij variables are non-zero, which gives a total of n + j |M j | variables. Since there are only n + j |M j | constraints, all constraints must be tight, which proves the lemma.

An O(n log n)-time algorithm for two machines

Given a feasible schedule, we call a job j a d-job, if |M j | = d. In this section we assume that the number of machines is two.

Lemma 4 Let σ be an optimal schedule for a 2-machine instance that satisfies the properties of Lemmas 1, 2 and 3. Let j < k be two consecutive 2-jobs. If there are 1jobs between j and k, then there is at least one 1-job on each machine. Also, the last 2-job is either not followed by any job or is followed by at least one 1-job on each machine.

Proof Let j and k be two consecutive 2-jobs and assume there is at least one in-between 1-job on machine 1 and none on machine 2. Let s 1 , s 2 be the start time of job j on respectively machine 1 and 2. We may assume without loss of generality that s 1 ≥ s 2 : otherwise we just swap the schedules of the two machines for the interval [0, C j ] and get the inequality. We change the schedule of j and k and the in-between 1-jobs as follows. Job j is completely processed on machine 2, starting from time s 2 , and the in-between 1-jobs are moved forward such that the first starts at time s 1 . Let ∆ be the amount of processing on job j that was previously assigned to machine 1, where we note that ∆ ≤ 1 2 p j . We increase the part of job k on machine 1 by ∆, and decrease the part of job k on machine 2 by ∆. This is possible, since the part of job k that was previously on machine 2 is at least 1 2 p k ≥ 1 2 p j ≥ ∆. The completion time of each of the in-between 1jobs decreases by ∆ + s, the completion time of job j increases by ∆ and the completion time of job k remains unchanged. The total completion time is thus reduced by at least s. If j is the last 2-job then we can make the same adjustment.

Lemma 5 In the case of two machines there are no 1jobs after a 2-job in an optimal schedule satisfying the properties of Lemmas 1, 2 and 3.

Proof Suppose the lemma is not true. Then there must be a 2-job j that is directly followed by a 1-job. By Lemma 4, there must be at least one such 1-job on each machine, say jobs h and k. Assume without loss of generality that p h ≤ p k . Let x 1j , x 2j be the processing time of j on machine 1 and 2, respectively. As argued before, without loss of generality we assume that x 1j ≥ x 2j . Let us define the starting time of j as zero, and let ∆ = x 1j -x 2j . Note that C j = 1 2 (∆ + p j + 2s). Then, the sum of the three completion times is

C j + C h + C k = C j + (C j + p h + s) + C k = ∆ + p j + 2s + p h + s + C k . (1)
We reschedule the jobs j, h and k as follows, while the remaining schedule stays the same. Place job j, the shortest among j, h and k, on machine 1 (unsplit), job h on machine 2 (unsplit), and behind these two, job k is split on machine 1 and 2, in such a way that it completes on one machine at time C h and at time C k on the other. The sum of the completion times of the three jobs becomes

(p j + s) + (∆ + p h + s) + C k ,
which is exactly s less than the sum of the three completion times in (1) from before the switch.

Given the previous lemmas, we see that the 2-jobs are scheduled in SPT order at the end. By Lemma 2, the first 2-job, say job k, is not shorter than the preceding 1jobs. But this implies that the 1-jobs can be scheduled in SPT order without increasing the completion time of job k and the following jobs. By considering each of the n jobs as the first 2-job, we immediately obtain a O(n 2 )-time algorithm to solve the problem. Carefully updating consecutive solutions leads to a faster method.

Theorem 1 There exists an O(n log n) algorithm for minimising the total completion time of jobs on two identical parallel machines with job splitting and uniform setup times.

Proof Suppose we schedule the first k jobs (for any 1 ≤ k ≤ n) in SPT order as 1-jobs and the other jobs in SPT order as 2-jobs. We would like to compute the change in objective value that results from changing job k from a 1-job to a 2-job. However, this happens to give a rather complicated formula. It is much easier to consider the change for job k -1 and k simultaneously.

The schedule for the 1-jobs j < k -1 does not change. To facilitate the exposition, suppose that job k-1 starts at time zero and job k starts at time a. Then C k-1 +C k = p k-1 +s+a+p k +s. After turning the jobs into 2-jobs, the new completion times become C k-1 = (a + p k-1 + 2s)/2 and C k = (a + p k-1 + p k + 4s)/2. Hence,

C k-1 + C k -C k-1 -C k = s -p k /2.
In addition, each job j > k completes s time units later. Hence, the total increase in objective value due to turning both job k -1 and k from a 1-job into a 2-job is

f (k) := (n -k + 1)s -p k /2.
Notice that f (k) is decreasing in k, since s > 0 and p k is nondecreasing in k. Hence, either there exists some k ∈ {2, . . . , n} such that f (k) < 0 and f (k -1) ≥ 0, or either f (n) ≥ 0, or f (2) < 0.

Suppose there exists some k ∈ {2, . . . , n} such that f (k) < 0 and f (k -1) ≥ 0. The optimal schedule is to have either k -1 or k -2 unsplit jobs, since the first inequality and monotonicity implies that a schedule with k -2 unsplit jobs has a better objective value than a schedule with k or more unsplit jobs, and the second inequality and monotonicity implies that a schedule with k -1 unsplit jobs has a better objective value than a schedule with k -3 or fewer unsplit jobs.

If f (n) ≥ 0 then the optimal solution is either to have only 1-jobs or have only job n as a 2-job. If f (2) < 0 then the optimal solution is either to have only 2-jobs or have only job 1 as a 1-job.

Straightforward implementation of the above gives the desired algorithm, the running time of which is dominated by sorting the jobs in SPT order.

Troubles on more machines

The properties exposed in Section 2 have been proven to hold for any number of machines. The properties presented in Section 3 were shown specifically for two machines only. In this section we investigate their analogues for three and more machines. We will present some examples of instances that show that the extension is far from trivial. It keeps the complexity of the problem on three and more machines as an intriguing open problem.

Lemma 5 shows that for two machines, there always exists an optimal schedule in which |M j | is monotonically non-decreasing in j. The following lemma shows that this does not hold for an arbitrary number of machines.

Lemma 6 There exist instances for which there is no optimal schedule in which |M j | is monotonically nondecreasing in j. Proof Consider the instance on three machines having 10 jobs with their vector of processing times p = (3, 10, 10, 10, 10, 50, 50, 50, 50, 50) (1 small job, 4 medium-sized jobs and 5 large jobs) and s = 0.7. We slightly perturb the processing times if necessary, obtaining p j < p j+1 for all j = 1, 2, . . . , n -1. We found all optimal solutions for this instance by exhaustive search. An optimal solution is depicted in Figure 2(a). As we see, job 2 is split over machines 2 and 3, but job 3 starting later than job 2 is not split. Jobs 4 and 5 are again what we call 2-jobs and are split over machines 2 and 3. The large jobs are all split over all three machines.

Below, we will describe all other optimal solutions to this instance. We will consider two solutions to be the same, if one solution can be obtained from the other by a relabelling of machines, and/or (repeatedly) swapping the schedule of two machines from some time t till the end of the schedule, if these two machines both complete processing of some job at time t.

The second optimal schedule, in Figure 2(b), is obtained by scheduling job 1 on machine 1, job 2 split on machines 2 and 3, job 3 on machine 2 (or 3), and jobs 4 and 5 as split jobs on the machines not used by job 3. The remaining jobs are again all split on all three machines. It is easily verified that the objective of this schedule is the same as the objective of the schedule in Figure 2(a): the completion time of job 3 increases by 2, and the completion times of jobs 4 and 5 each decrease by 1, and all other completion times remain the same. The remaining two optimal schedules, in Figures 2(c) and 2(d) are obtained by switching jobs 3 and 4 in the first two optimal schedules. We note that these schedules continue to be optimal if the processing times are slightly perturbed, as mentioned earlier.

All optimal solutions for this instance share the property that job 2 is a 2-job, and either job 3 or job 4 is a 1-job, which proves the lemma.

If we slightly change the instance from the proof of Lemma 6 by deleting one of the large jobs, then there is a unique optimal solution, which splits job 3 over machines 1 and 2 and continues with splitting job 4 over machines 1 and 3. Job 5 and the four large jobs are split over all three machines, see Figure 3. Lemma 6 and the fact that a subtle change in the problem instance causes such a substantial change in the optimal schedule bodes ill for an algorithmic approach like the one in Section 3. Fig. 3: Gantt chart depicting the unique optimal solutions to the 3-machine instance with processing times p = (3, 10, 10, 10, 10, 50, 50, 50, 50) (1 small job, 4 medium-sized jobs and 4 large jobs) and s = 0.7.

Approximation algorithm

We will now show a constant-factor approximation algorithm for our problem, for an arbitrary number of machines. We remark that we do not know whether this problem is NP-hard, but the examples in the previous section do show that the way a job is scheduled in an optimal schedule may depend on jobs that occur later in the schedule. Our approximation algorithm, on the other hand, is remarkably simple, and only uses a job's processing time and the setup time to determine how to schedule the job. We schedule the jobs in order of non-decreasing processing time. Let s > 0 and let α be some constant that will be determined later. Job j will be scheduled such that it completes as early as possible under the restriction that it uses at most j := min{ αp j /s , m} machines. Thus, the job will be scheduled on the at most j machines that have minimum load in the schedule so far. It is easy to see that a job is always balanced this way.

Theorem 2 The algorithm described above is a (2+α)approximation algorithm for minimising the total completion time with job splitting and uniform setup times, provided that α ≥ 1 4 ( √ 17 -1).

Proof Let σ be the schedule produced by the described algorithm. Note that the total load (processing times plus setup times) of all jobs in σ up to, but not including, job j is upper bounded by L j = k<j (p k + k s), since job k introduced at most k setups. Therefore, the average load on the j least loaded machines is upper bounded by L j /m. Since job j is balanced, we can thus upper bound the completion time Cj of job j in the schedule by L j /m + p j / j + s. Note that this is an upper bound on the completion time of job j when we try to schedule it on at most j machines. Noting that p j / j = p j / min{ αp j /s , m} ≤ p j / αp j /s + p j /m ≤ (1/α)s + p j /m, and

k s = min{ αp k /s , m}s < αp k + s, we obtain Cj ≤ L j /m + p j / j + s ≤ 1 m k<j (p k + k s) + p j / j + s < 1 m k<j (1 + α)p k + s + p j /m + (1 + 1/α)s ≤ 1 + α m k≤j p k + j -1 m + 1 + 1 α s.
We can lower bound the sum of completion times in an optimal schedule by j (s+ 1 m k≤j p k ): suppose we only needed a setup time for the first job to be processed on a machine, for any machine. Clearly, the optimal sum of completion times for this problem gives a lower bound on the optimum for the original problem. Now, the optimal schedule when we only need a setup time for the first job on a machine processes the jobs in SPT order and splits each job over all machines, which gives a sum of completion times of j (s + 1 m k≤j p k ). Also, in any schedule, at most m jobs are preceded by only one setup, at most another m by two setups, etc., giving a lower bound of j j/m s on the sum of completion times: this is exactly the optimal value when all processing times are 0. We will show below that j j m s ≥ j j-1 m s + 1 2 ns. Hence, by using 1 + α times the first bound, and 1 time the second bound, we get

(2 + α) j C j ≥ (1 + α) j s + 1 m k≤j p k + j j-1 m s + 1 2 ns = j 1+α m k≤j p k + (1 + α)s + j-1 m s + 1 2 s ,
which is at least as large as j Cj provided α > 0 and

3 2 + α ≥ 1 + 1 α , which is equivalent to α ≥ 1 4 ( √ 17 -1).
Next we show that j j m s ≥ j j-1 m s + 1 2 ns. Let j = qm + a for some q ≥ 0 and a ∈ {1, . . . , m}. Then

= r m a=1 1 - a -1 m + b a=1 1 - a -1 m = rm + b -r m a=1 a -1 m - b a=1 a -1 m = n -r(m -1)/2 - 1 2 (b -1)b/m ≥ n -r(m -1)/2 - 1 2 (b -1) = n -(rm + b)/2 + r/2 + 1/2 = n/2 + r/2 + 1/2 ≥ n/2.
Hence multiplying both sides with s yields

n j=1 j m s ≥ n j=1 j -1 m s+ 1 2 ns.
Corollary 1 There exists a 2 + 1 4 ( √ 17 -1) < 2.781approximation algorithm for minimising total completion time with job splitting and uniform setup times.

A polynomial-time approximation scheme

We give an approximation scheme which runs in polynomial time if the number of machines is assumed constant. The idea is simple: by splitting a job j, at most p j on its completion time can be saved. Denote by Opt the sum of completion times of an optimal schedule. It is easy to show that the value of a non-preemptive SPT schedule is no more than j p j larger than Opt. In particular, if we schedule the first K = n -m/ jobs by non-preemptive SPT then the extra cost is at most K j=1 p j . But, as we will see, this is only an -fraction of the total completion time of the last m/ jobs. These last jobs we schedule optimally given the schedule of the first K jobs. Now, we define the algorithm and its running time in more detail. Let, as before,

p 1 ≤ • • • ≤ p n . Let K = n -m/ . (If K ≤ 0 then n ≤ m/
and the optimal solution can be found in constant time.) Let ρ be an optimal schedule and let ρ(K) be the schedule ρ restricted to the jobs 1, 2, . . . , K. By Lemma 2 we may assume that ρ(K) has no idle time. Let t i (ρ) be the completion time of machine i in ρ(K). The algorithm makes an approximate guess about the values t i (ρ). That means, it finds values t i such that

t i (ρ) ≤ t i ≤ t i (ρ) + s + p K . (2)
Note that for any i, we have t i (ρ) ≤ K(s + p K ). Hence, we need to try only K m guesses for (t 1 , . . . , t m ). Assume from now that we guessed (t 1 , . . . , t m ) correctly, i.e., (2) is satisfied.

We apply SPT to the jobs 1, 2, . . . , K such that no machine i is loaded more than t i +s+p K . This can easily be done as follows: apply list scheduling in SPT order and close a machine once its load becomes t i or more. Let T i be the completion time of machine i in the resulting schedule. Then T i ≤ t i + s + p K ≤ t i (ρ) + 2(s + p K ). Next, we find a near-optimal completion of the schedule by guessing for each job j > K a set M j and apply linear programming. There are 2 m(n-K) possibilities for choosing such sets, which is a constant. The linear program works as follows. Note that the LP of Section 2 can be extended to do the following. Given a set M j for each job j and a time T i for each machine, we can find the optimal schedule among all schedules for which: (i) job parts are in SPT order on each machine, (ii) machine i does not start before T i , (iii) job j can only be scheduled on machines in M j , and (iv) job j has a setup time s for each machine in M j even when its processing time x ij is zero. Note that it is not clear if the LP gives us the real optimal completion since we have not proved that the SPT properties hold also for optimal schedules if an initial part is fixed, as we do here. However, we can show that the solution given by the LP is close to optimal.

Approximation ratio Let σ be the final schedule and let Cj be the completion time of job j. Here we use Opt to denote the objective value of optimal schedule ρ. For any h ∈ {1, . . . , n} define µ h = h k=1 (s + p k )/m. Then for any schedule, the h-th completion time is at least µ h . Hence,

Opt ≥ n h=1 µ h ≥ n h=K+1 µ h ≥ n h=K+1 µ K = m/ µ K ≥ 1 K k=1 (s + p k ).
For the rest the proof we will denote by C h the completion time of job h in the optimal schedule and by Ch the completion time of job h in the schedule produced by our algorithm. Further, we will use the notation C (h) for the h-th completion time of ρ, (h = 1, . . . , n). Notice that C (h) is not necessarily equal to C h . For h ≤ K it is easy to see that Ch ≤ C (h) + s + p h . This implies

K h=1 Ch ≤ K h=1 (C (h) + s + p h ) = K h=1 C (h) + K h=1 (s + p h ) (3)
≤ K h=1 C (h) + Opt.
So for the first K jobs we are doing fine. Next, we give a bound on the total completion time of the other jobs. Let M * j be the set of machines used by job j in the optimal schedule ρ. One of the guesses of the algorithm will be M j = M * j for j > K. We show that the corresponding LP-solution gives a near-optimal completion of the schedule.

A feasible LP solution is to take for x ij , j > K, the values that correspond to ρ and choose values C LP j = C j + 2(s + p K ), where we we remind that C j is the completion time of job j in the optimal schedule ρ. The latter is feasible since T i ≤ t i (ρ) + 2(s + p K ). Hence, we can bound the total completion times of jobs K + 1, . . . , n by Adding (3) we can bound the total completion time by (1 + 5 )Opt.

Hardness for weighted completion times

We prove that introducing weights for the jobs in our problem makes it strongly NP-hard for any number of machines and weakly NP-hard for 2 machines.

Theorem 3 The problem of minimising total weighted completion time with job splitting and uniform setup times on parallel identical machines (P |s, split| w j C j ) is strongly NP-hard.

Proof We reduce from 3-Partition: given 3n positive numbers a 1 , . . . , a 3n and a number A such that a 1 + • • • + a 3n = nA, does there exist a partition A 1 , . . . , A n of {1, . . . , 3n} such that |A i | = 3 and j∈Ai a j = A for all i? Given an instance of 3-Partition, we construct the following instance of our scheduling problem: We have n machines and 3n jobs. We set p j = a j and w j = a j + s for all j = 1, . . . , 3n, where the setup time s is some large enough number, to be defined later.

The idea behind the reduction is the following: the large setup time will make sure that exactly three jobs are scheduled (unsplit) per machine. The weights are chosen such that a schedule where all machines complete at exactly the same time is optimal, if such a schedule is feasible.

Suppose we schedule the jobs unsplit where A i is the set of jobs processed on machine i. Then, the cost of the schedule is: where l i is the total load on machine i. Note that the second term is independent of the schedule. This cost is minimised when l i = l h for all i and h and this can be realised if a perfect 3-partition exists. Let us denote this minimum by Opt 3P .

Fig. 1 :

 1 Fig.1: Gantt charts depicting the schedules for the instance described in Section 2. The grey blocks indicate the setup times, the numbered blocks are scheduled job parts. Each row of blocks gives the schedule for a machine.

Fig. 2 :

 2 Fig.2: Gantt charts depicting the optimal solutions to the 3-machine instance with processing times p = (3, 10, 10, 10, 10, 50, 50, 50, 50, 50) (1 small job, 4 medium-sized jobs and 5 large jobs) and s = 0.7. The grey blocks indicate the setup times, the numbered blocks are scheduled job parts. Each row of blocks gives the schedule for a machine.

j m - j - 1 m

 1 = (q + 1) -(qm + a -1)/m = 1 -(a -1)/m.Now assume that n = rm + b, for some integer r ≥ 0 and b ∈ {1, . . . , m}. Then

C

  h + 2(n -K)(s + p K ).To bound the second term in the right-hand side of (4) we derive another bound on Opt: K)(s + p K )/ .Combining this with (4) we getn h=K+1 Ch ≤ n h=K+1 C h +4 Opt ≤ n h=K+1 C (h) +4 Opt.

  a j ) 2 ,

Acknowledgements

The authors wish to acknowledge an anonymous reviewer for detailed and helpful comments on an earlier version of this manuscript.

Research supported by EU-IRSES grant EUSACOU

If no perfect 3-partition exists, then any schedule where no jobs are split has strictly higher cost than Opt 3P . It remains to prove that also any schedule with at least one split job has a strictly higher cost than Opt 3P .

First observe that

(s + a j ) 2 = 6ns 2 + O(ns).

Now assume that at least one job is split, then there are at least 3n + 1 setup times of s each. Consider the extreme case where all 3n values a j are zero. In this case it is easy to see that the weighted sum of the 3n completion times is at least (6n + 1)s 2 . Clearly, this bound holds as well for arbitrary value a j . For large enough s we have (6n + 1)s 2 > Opt 3P .

Theorem 4

The problem P 2|s, split| w j C j is weakly NP-hard.

Proof We now reduce from a restricted form of the Subset Sum problem: Given 2n positive integers a 1 , . . . , a 2n such that a 1 + • • • + a 2n = 2A, is there a set I ⊂ {1, . . . , 2n} such that |I| = n and i∈I a i = A? Given an instance of Subset Sum, we construct the following instance of our scheduling problem. We have 2 machines and 2n jobs. We set p j = a j and w j = a j + s for j = 1, . . . , 2n, where the setup time s is some large enough number, to be defined later. The proof follows the same reasoning as the previous proof: the large setup time will now make sure that exactly n jobs are scheduled (unsplit) per machine, and the weights will make sure that a schedule where the two machines complete at exactly the same time is optimal, if such a schedule is feasible. Suppose we schedule the jobs unsplit. Then, just as in the proof above for an arbitrary number of machines we have that the cost of the schedule is:

where l i is the total load on machine i. Note that the second term is independent of the schedule. This cost is minimised when l 1 = l 2 and this can be realised if a perfect subset I exists. Let us denote this minimum by Opt S .

If no perfect subset exists, any unsplit schedule has strictly higher cost. It remains to prove that also any schedule with at least one split has a strictly higher cost than Opt S .

First observe that

Now assume that at least one job is split, then there are at least 2n + 1 setup times of s each. Consider the extreme case where all 2n values a j are zero. In this case it is easy to see that the weighted sum of the 2n completion times is at least (n 2 + n + 1)s 2 . Clearly, this bound holds as well for arbitrary values a j . For large enough s we have (n 2 + n + 1)s 2 > Opt S .

Epilogue

In the following table we gather the state of the art on scheduling problems with job splitting and uniform setup times. For describing the problems in the first column of the table we use the standard three-field scheduling notation [START_REF] Graham | Optimization and approximation in deterministic sequencing and scheduling: a survey[END_REF]. In the first field, expressing the processor environment, we only consider parallel identical machines, denoted by P , possibly with the number of parallel machines mentioned additionally. In the second field, expressing job characteristics, the term 'pmtn' denotes ordinary preemption, 'split' denotes job splitting as we consider in this paper and s denotes the presence of uniform setup times. Though this paper is mainly concerned with problems with a total completion time objective, indicated by C j in the third field, expressing the objective, we will also show the state of the art on the total weighted completion time (indicated by w j C j ) and on the makespan (indicated by C max ).

In the second column, we summarize the complexity status of these problems. A question mark indicates that the complexity of the problem is unknown. In the third column we give the best approximation guarantee known, where a '-' indicates that no algorithm with a performance guarantee is known. If we consider it relevant, we also present, as a footnote, the knowledge on the comparable version with preemption instead of splitting.


Hassan S Bakouch
email: hnbakouch@yahoo.com

Christophe Chesneau
email: christophe.chesneau@unicaen.fr

Nauman Muhammad

Khan

On a weighted exponential distribution with a logarithmic weight: Theory and applications

Keywords: Statistical distributions, Survival function, Hazard rate function, Record values, Parametric estimation. 2000 MSC: 60E05, 62E15, 62F10

A new weighted exponential distribution is proposed using a logarithmic weight function. Mathematical properties are studied, including moments, residual life function, order statistics and record values. A parametric estimation is performed using the maximum likelihood method, the weighted least-square method and the maximum product of spacing estimators method. Applications are provided using four real data sets. It is shown that our new distribution works better than a number of classical and recent distributions.

Introduction

The general weighted exponential family of distributions has a probability density function (pdf) of the form:

f (x) = Cw(x)g(x),
where g(x) denotes a pdf of the exponential distribution with parameter λ > 0, w(x) ≥ 0 is a weight function and C represents the normalization constant given by C -1 = ∫ +∞ -∞ w(x)g(x)dx. Recent weighted exponential distributions with great ability in statistical applications include the exponentiated exponential (EE) distribution introduced by [START_REF] Gupta | Exponentiated Exponential Family; An Alternative to Gamma and Weibull[END_REF], the weighted exponential (WE) distribution introduced by [START_REF] Gupta | A new class of weighted exponential distribution[END_REF], the Gamma exponentiated exponential (GEE) distribution introduced by [START_REF] Ristić | The gamma-exponentiated exponential distribution[END_REF], the weighted generalized exponential (WGE) distribution introduced by [START_REF] Mahdavi | Two Weighted Distributions Generated by Exponential Distribution[END_REF] and the exponentiated weighted exponential (EWE) distribution introduced by [START_REF] Oguntunde | On the Exponentiated Weighted Exponential Distribution and Its Basic Statistical Properties[END_REF]. A brief review on weighted distributions (exponential or non-exponential) can be found in [START_REF] Saghir | Weighted Distributions: A Brief Review, Perspective and Characterizations[END_REF]. The statistical literature on proposing new weighted exponential distributions is vast and growing fast.

In this paper, we introduce a new weighted exponential distribution characterized by the pdf:

f (x) = λβ (β + 1) log(β + 1)
log(β + e λx )e -λx , x, β, λ > 0.

(

) 1
Thus the associated weight is w(x) = log(β + e λx ) and C = β (β + 1) log(β + 1)

. To the best of our knowledge, the consideration of this weight is new in the literature, opening the door of new perspectives of applications. We shall refer to the distribution given by (1) as the logarithmic weighted exponential (LWE). Also, note that w(x) is an increasing and concave up, hence one of the possible shapes of the hazard rate function (hrf) of this distribution is increasing ( [START_REF] Jain | Relations for reliability measures of weighted distributions[END_REF]), as it shall be seen later.

Moreover, some plots are given to show the graphical behavior of the pdf and hrf of the distribution, illustrating its flexibility in data modeling. Other motivations to the LWE distribution are as follows.

• The pdf (1) can be written as a two-component mixture of two pdfs:

f (x) = pf 1 (x) + (1 -p)f 2 (x), x > 0,
where p = β (β + 1) log(β + [START_REF] Ahsanullah | Record Statistics[END_REF] is the mixing proportion, f 1 (x) = λ 2 xe -λx is the pdf of the Gamma distribution with parameters (2, λ), f 2 (x) = λβ (β + 1) log(β + 1) -β log

( β + 1 -β(1 -e -λx
) ) e -λx is an extension of the pdf of gamma-G distribution introduced by [START_REF] Zografos | On families of beta-and generalized gamma-generated distributions and associated inference[END_REF]. To be more specific, it is associated to a cumulative distribution function (cdf) of the form H(G(x)), belonging to the family of composition of two cdfs: H(y) and G(x) (see [START_REF] Al-Hussaini | Composition of cumulative distribution functions[END_REF]). Here G(x) depends only on λ; it is the cdf associated to the exponential distribution of parameter λ and H(y) depends only on β:

H(y) = 1 - 1 (β + 1) log(β + 1) -β [(1 + β(1 -y)) log (1 + β(1 -y)) -β(1 -y)] , y ∈ [0, 1].
One can show that H(y) is a slowly decreasing according to β > 0 with H(y) ∼ 1 -(1 -y) 2 when β → 0 and H(y) ∼ y when β → +∞.

• Let X follow the LWE distribution and consider the transformation Y = β + e λX . Then distribution of Y is characterized by the pdf given by

f (y) = β (β + 1) log(β + 1) log(y) (y -β) 2 , y ≥ β + 1, β > 0,
which is a new extension of the Pareto distribution and may be called log-Pareto distribution. The rest of the paper is organized as follows. Section 2 presents the cdf and hrf related to the LWE distribution, with some plots. Section 3 is devoted to some of its statistical properties. Residual life function is studied in Section 4. Order statistics and record values are investigated in Section 5. Three different estimation procedures are presented in Section 6. Finally, four real data applications are presented in Section 7 to show the performance of the distribution.

Functions related to the LWE distribution

Survival, cumulative density and hazard rate functions

Proposition 1. Let X be a random variable with pdf f (x) [START_REF] Ahsanullah | Record Statistics[END_REF]. The associated survival function (sf ) is given by :

S(x) = P (X > x) = log(1 + βe -λx ) + βe -λx log(β + e λx ) (β + 1) log(β + 1) , x, β, λ > 0. ( 2
)
Proof. Let us set c = λβ (β+1) log(β+1) and remark that

f (x) = c log(β + e λx )e -λx = c log ( e λx (1 + βe -λx ) ) e -λx = c ( λxe -λx + log(1 + βe -λx )e -λx
) .

Using definition of S(x), we find that

S(x) = ∫ +∞ x f (t)dt = c (∫ +∞ x λte -λt dt + ∫ +∞ x log(1 + βe -λt )e -λt dt
) .

Using this decomposition, integration by parts, change of variables and after some algebra, we obtain

S(x) = log(1 + βe -λx ) + βe -λx log(β + e λx ) (β + 1) log(β + 1) .
This ends the proof.

The cumulative distribution function follows immediately:

F (x) = 1 -S(x) = 1 - log(1 + βe -λx ) + βe -λx log(β + e λx ) (β + 1) log(β + 1) , x, β, λ > 0. (3)
The hazard rate function is given by

h(x) = f (x) 1 -F (x) = λβ log(β + e λx )e -λx log(1 + βe -λx ) + βe -λx log(β + e λx ) , x, β, λ > 0. ( 4
)

Behaviour of the functions

First of all, let us investigate some features of the functions inherent to the LWE distribution. For the pdf (1), we have

f (x) ∼ λβ β + 1 e -λx → λβ β + 1 , x → 0, f (x) ∼ λ 2 β (β + 1) log(β + 1) xe -λx → 0, x → +∞.
For the sf (2), we have

S(x) ∼ 1 - λ log(β + 1) x → 1, x → 0, S(x) ∼ β (β + 1) log(β + 1) (1 + λx)e -λx → 0, x → +∞.
And for the cdf (3), we have

F (x) ∼ λ log(β + 1) x → 0, x → 0, F (x) ∼ 1 - β (β + 1) log(β + 1) (1 + λx)e -λx → 1, x → +∞.
For the hrf (4), we have

h(x) ∼ λβ β + 1 e -λx → λβ β + 1 , x → 0, h(x) ∼ λ 2 x 1 + λx → λ, x → +∞.
Let us now study the change points of h(x). By denoting y = log(β + e λx ), we have

h ′ (x) = βλ 2 e -λx [ -λxe λx -(β + e λx y(-1 -λx + (1 + 2βe λx )y)) ] (β + e λx )(xλ -(1 + βe λx )y) 2 .
The change points of the hrf function are given by b * by solving the two following equations: b * λ -(1 + βe λb * )y 1 = 0 and -λb * e λb * -(β + e λb * y 1 (-1 -λb * + (1 + 2βe λb * )y 1 )) = 0. Figures 1 and2 show the graphical features of the pdf and hrf of the LWE distribution for several choices of parameters (λ, β). In particular, we see that the pdf is very flexible and highly right skewed. Interesting comportment of the hrf, with possible change-point(s), can be observed.

Statistical properties

This section is devoted to some statistical properties of the LWE distribution.

Expansion of the pdf

Proposition 2. Let f (x) be (1). Then we have the following expansion:

f (x) = λ 2 β (β + 1) log(β + 1) xe -λx + ∞ ∑ k=0 k ∑ ℓ=0 A k,ℓ e -λ(ℓ+1)x , ( 5
)
with

A 0,0 = λβ (β + 1) log(β + 1) log(β + 1), A k,ℓ = λβ (β + 1) log(β + 1) ( k ℓ ) 1 k ( β β + 1 ) k (-1) ℓ+1 ,
for k > 0 or ℓ > 0.

β β + 1 (1 -e -λx )
)) e -λx = c ( log(β + 1) + λx + log

( 1 - β β + 1 (1 -e -λx )
)) e -λx .

Using the expansion: log(1 [START_REF] Ahsanullah | Record Statistics[END_REF], noticing that 0 ≤ β β+1 (1 -e -λx ) < 1, and the binomial series expansion, we have log

-y) = - ∞ ∑ k=1 y k k , y ∈ [-1,
( 1 - β β + 1 (1 -e -λx ) ) = - ∞ ∑ k=1 1 k ( β β + 1 ) k (1 -e -λx ) k = ∞ ∑ k=1 k ∑ ℓ=0 ( k ℓ ) 1 k ( β β + 1 ) k (-1) ℓ+1 e -λℓx .
We end the proof by putting these equalities together.

Quantile function

We can obtain the quantile function using the cdf (3); it is the function Q(x) satisfying the nonlinear equation:

F (Q(x)) = x ⇔ log(1 + βe -λQ(x) ) + βe -λQ(x) log(β + e λQ(x) ) = (1 -x)(β + 1) log(β + 1). (6)
Variates of the LWE distribution can be simulated by using X = Q(U ) where U is a random variable following the uniform distribution on [0, 1].

Moments Proposition 3. Let X be a random variable with pdf f (x) given by (1). Consider the Gamma function: Γ(ν) =

∫ +∞ 0

x ν-1 e -x dx and the coefficient A k,ℓ of expansion [START_REF] Barlow | A Bayesian analysis of stress rupture life of Kevlar 49/epoxy spherical pressure vessels[END_REF]. Then the r-th moment of X is given by

E(X r ) = βΓ(r + 2) (β + 1) log(β + 1)λ r + ∞ ∑ k=0 k ∑ ℓ=0 A k,ℓ Γ(r + 1) λ r+1 (ℓ + 1) r+1 . ( 7
)
Proof. Let us set c = λβ (β+1) log(β+1) and consider the expansion [START_REF] Barlow | A Bayesian analysis of stress rupture life of Kevlar 49/epoxy spherical pressure vessels[END_REF]. Hence

E(X r ) = ∫ +∞ -∞ x r f (x)dx = ∫ +∞ 0 x r ( cλxe -λx + ∞ ∑ k=0 k ∑ ℓ=0 A k,ℓ e -λ(ℓ+1)x ) dx = cλ ∫ +∞ 0 x r+1 e -λx dx + ∞ ∑ k=0 k ∑ ℓ=0 A k,ℓ ∫ +∞ 0 x r e -λ(ℓ+1)x dx = c Γ(r + 2) λ r+1 + ∞ ∑ k=0 k ∑ ℓ=0 A k,ℓ Γ(r + 1) λ r+1 (ℓ + 1) r+1 .
The proof is completed.

Moment generating function Proposition 4.

Let X be a random variable with pdf f (x) given by [START_REF] Ahsanullah | Record Statistics[END_REF]. Consider the coefficient A k,ℓ of expansion [START_REF] Barlow | A Bayesian analysis of stress rupture life of Kevlar 49/epoxy spherical pressure vessels[END_REF]. Then the moment generating function of X is given by

M (t) = λ 2 β (β + 1) log(β + 1)(λ -t) + ∞ ∑ k=0 k ∑ ℓ=0 A k,ℓ 1 λ(ℓ + 1) -t , t < λ.
Proof. Let us set c = λβ (β+1) log(β+1) and consider the expansion [START_REF] Barlow | A Bayesian analysis of stress rupture life of Kevlar 49/epoxy spherical pressure vessels[END_REF]. Hence, for t < λ,

M (t) = E(e tX ) = ∫ +∞ -∞ e tx f (x)dx = ∫ +∞ 0 e tx ( cλxe -λx + ∞ ∑ k=0 k ∑ ℓ=0 A k,ℓ e -λ(ℓ+1)x ) dx = cλ ∫ +∞ 0 e -(λ-t)x dx + ∞ ∑ k=0 k ∑ ℓ=0 A k,ℓ ∫ +∞ 0 e -(λ(ℓ+1)-t)x dx = c λ λ -t + ∞ ∑ k=0 k ∑ ℓ=0 A k,ℓ 1 λ(ℓ + 1) -t .
This ends the proof.

Conditional moments and mean deviations

Proposition 5. Let X be a random variable with pdf f (x) given by [START_REF] Ahsanullah | Record Statistics[END_REF]. Consider the lower incomplete Gamma function Γ(t, ν) = ∫ t 0 x ν-1 e -x dx and the coefficient A k,ℓ of expansion [START_REF] Barlow | A Bayesian analysis of stress rupture life of Kevlar 49/epoxy spherical pressure vessels[END_REF]. Then we have

∫ t 0 x r f (x)dx = βΓ(λt, r + 2) (β + 1) log(β + 1)λ r + ∞ ∑ k=0 k ∑ ℓ=0 A k,ℓ Γ(λ(ℓ + 1)t, r + 1) λ r+1 (ℓ + 1) r+1 , t > 0. ( 8
)
Proof. Let us set c = λβ (β+1) log(β+1) and consider the expansion [START_REF] Barlow | A Bayesian analysis of stress rupture life of Kevlar 49/epoxy spherical pressure vessels[END_REF]. Therefore

∫ t 0 x r f (x)dx = ∫ t 0 x r ( cλxe -λx + ∞ ∑ k=0 k ∑ ℓ=0 A k,ℓ e -λ(ℓ+1)x ) dx = cλ ∫ t 0 x r+1 e -λx dx + ∞ ∑ k=0 k ∑ ℓ=0 A k,ℓ ∫ t 0 x r e -λ(ℓ+1)x dx = c Γ(λt, r + 2) λ r+1 + ∞ ∑ k=0 k ∑ ℓ=0 A k,ℓ Γ(λ(ℓ + 1)t, r + 1) λ r+1 (ℓ + 1) r+1 .
The proof is completed.

The r-th conditional moments of the LWE distribution is given by

E(X r | X > t) = 1 1 -F (t) ∫ +∞ t x r f (x)dx = 1 1 -F (t) ( E(X r ) - ∫ t 0 x r f (x)dx
) .

The r-th reversed moments of the LWE distribution is given by

E(X r | X ≤ t) = 1 F (t) ∫ t 0 x r f (x)dx.
Both moments can be expressed using [START_REF] Jain | Relations for reliability measures of weighted distributions[END_REF]. The mean deviations of X about the mean µ = E(X) can be expressed as δ = 2µF (µ) -2 ∫ µ 0 xf (x)dx and the mean deviations of X about the median M has the form

η = µ -2 ∫ M 0 xf (x)dx.

Residuals life function

Let X be a random variable having the pdf f (x) given by ( 1). The residual life is described by the conditional random variable R (t) = X -t | X > t, t ≥ 0. It naturally appears in life length studies (see [START_REF] Bryson | Some criteria for aging[END_REF] and [START_REF] Muth | Reliability models with positive memory derived from the mean residual life function[END_REF]). Using the sf S(x) (2), the sf of the residual lifetime R (t) is given by

S R (t) (x) = S(x + t) S(t) = log(1 + βe -λ(x+t) ) + βe -λ(x+t) log(β + e λ(x+t) ) log(1 + βe -λt ) + βe -λt log(β + e λt ) , x > 0.
The associated cdf is given by

F R (t) (x) = log ( 1+βe -λt 1+βe -λ(x+t) ) + βe -λt ( log(β + e λt ) -e -λx log(β + e λ(x+t) ) ) log(1 + βe -λt ) + βe -λt log(β + e λt ) , x > 0.
Then, the corresponding pdf is given by

f R (t) (x) = λβ log(β + e λ(x+t) )e -λ(x+t) log(1 + βe -λt ) + βe -λt log(β + e λt ) , x > 0.
The associated hrf is given by

h R (t) (x) = λβ log(β + e λ(x+t) )e -λ(x+t) log(1 + βe -λ(x+t) ) + βe -λ(x+t) log(β + e λ(x+t) ) , x > 0.
The mean residual life is defined as

K(t) = E(R (t) ) = E(X -t | X > t) = 1 S(t) ∫ +∞ t xf (x)dx -t = 1 S(t) ( E(X) - ∫ t 0 xf (x)dx ) -t,
where f (x) is given by (1), S(t) is given by ( 2), E(X) is given by ( 7) with r = 1 and

∫ t 0 xf (x)
dx is given by ( 8) with r = 1.

On the other hand, the variance residual life is given by

V (t) = V ar(R (t) ) = V ar(X -t | X > t) = 2 S(t) ∫ +∞ t xS(x)dx -2tK(t) -[K(t)] 2 = 1 S(t) ( E(X 2 ) - ∫ t 0 x 2 f (x)dx ) -t 2 -2tK(t) -[K(t)] 2 ,
where E(X 2 ) is given by ( 7) with r = 2 and ∫ t 0 x 2 f (x)dx is given by ( 8) with r = 2. The reverse residual life is described by the conditional random variable R (t) = t -X | X ≤ t, t ≥ 0. Using the cdf (3), the sf of the reversed residual lifetime R (t) is given by

S R (t) (x) = F (t -x) F (t) = (β + 1) log(β + 1) -log(1 + βe -λ(t-x) ) + βe -λ(t-x) log(β + e λ(t-x) ) (β + 1) log(β + 1) -log(1 + βe -λt ) + βe -λt log(β + e λt ) , 0 < x ≤ t.
The associated cdf is given by

F R (t) (x) = log(1 + βe -λ(t-x) ) + βe -λ(t-x) log(β + e λ(t-x) ) -log(1 + βe -λt ) + βe -λt log(β + e λt ) (β + 1) log(β + 1) -log(1 + βe -λt ) + βe -λt log(β + e λt ) , 0 < x ≤ t.
Therefore, the corresponding pdf is given by

f R (t) (x) = λβ log(β + e λ(t-x) )e -λ(t-x) (β + 1) log(β + 1) -log(1 + βe -λt ) + βe -λt log(β + e λt ) , 0 < x ≤ t.
The associated hrf is given by

h R (t) (x) = λβ log(β + e λ(t-x) )e -λ(t-x) (β + 1) log(β + 1) -log(1 + βe -λ(t-x) ) + βe -λ(t-x) log(β + e λ(t-x) ) , 0 < x ≤ t.
The mean reversed residual life is defined as

L(t) = E(R (t) ) = E(t -X | X ≤ t) = t - 1 F (t) ∫ t 0 xf (x)dx,
where f (x) is given by (1), F (t) is given by (3) and ∫ t 0 xf (x)dx is given by ( 8) with r = 1. The variance reversed residual life is given by

W (t) = V ar(R (t) ) = V ar(t -X | X ≤ t) = 2tL(t) -[L(t)] 2 - 2 F (t) ∫ t 0 xF (x)dx = 2tL(t) -[L(t)] 2 -t 2 + 1 F (t) ∫ t 0 x 2 f (x)dx,
where ∫ t 0 x 2 f (x)dx is given by ( 8) with r = 2.

Order statistics and record values

Order statistics naturally appear in many areas of statistics, as reliability and quality control testing. Let X 1 , X 2 , . . . , X n be n i.i.d. random variables having the pdf (1). Let us consider its order statistics is X 1:n , X 2:n , . . . , X n:n . Using the expressions of f (x) given by (1) and F (x) (3), the pdf of the i-th order statistic X i:n is given by

f i:n (x) = n! (i -1)! (n -i)! [F (x)] i-1 [1 -F (x)] n-i f (x) = n! (i -1)! (n -i)! [ 1 - log(1 + βe -λx ) + βe -λx log(β + e λx ) (β + 1) log(β + 1) ] i-1 × [ log(1 + βe -λx ) + βe -λx log(β + e λx ) (β + 1) log(β + 1)
] n-i λβ (β + 1) log(β + 1) log(β + e λx )e -λx , x > 0.

In particular, let us mention that the pdf of X 1:n = inf(X 1 , . . . , X n ) is given by

f 1:n (x) = n[1 -F (x)] n-1 f (x) = λβ (β + 1) log(β + 1) n [ log(1 + βe -λx ) + βe -λx log(β + e λx ) (β + 1) log(β + 1)
] n-1 log(β + e λx )e -λx and the pdf of X n:n = sup(X 1 , . . . , X n ) is given by

f n:n (x) = n[F (x)] n-1 f (x) = λβ (β + 1) log(β + 1) n [ 1 - log(1 + βe -λx ) + βe -λx log(β + e λx ) (β + 1) log(β + 1)
] n-1 log(β + e λx )e -λx .

The cdf of the i-th order statistic X i:n is given by

F i:n (x) = n! (i -1)! (n -i)! n-i ∑ k=0 ( n -i k ) (-1) k i + k [F (x)] i+k = n! (i -1)! (n -i)! n-i ∑ k=0 ( n -i k ) (-1) k i + k [ 1 - log(1 + βe -λx ) + βe -λx log(β + e λx ) (β + 1) log(β + 1)
] i+k , x > 0.

For i < j, the joint pdf of (X i:n , X j:n ) is given by

f (i:n,j:n) (x i , x j ) = n! (i -1)! (n -j)! (j -i -1) [F (x i )] i-1 [F (x j ) -F (x i )] j-i-1 [1 -F (x j )] n-j f (x i )f (x j ) = n! (i -1)! (n -j)! (j -i -1) [ 1 - log(1 + βe -λxi ) + βe -λxi log(β + e λxi ) (β + 1) log(β + 1) ] i-1 ×   log ( 1+βe -λx i 1+βe -λx j ) + β(e -λxi log(β + e λxi ) -e -λxj log(β + e λxj )) (β + 1) log(β + 1)   j-i-1 × [ log(1 + βe -λxj ) + βe -λxj log(β + e λxj ) (β + 1) log(β + 1)
] n-j ( λβ (β + 1) log(β + 1)

) 2
× log(β + e λxi ) log(β + e λxj )e -λ(xi+xj ) , 0 < x i < x j .

As order statistics, record values arise in many situations in applied statistics. Theoretical and practical aspects of record values can be found in [START_REF] Ahsanullah | Record Statistics[END_REF] and [3]. Let X 1 , X 2 , . . . , be a sequence of i.i.d. random variables having the pdf (1). We define a sequence of record times U (n) as follows: U (1) = 1, U (n) = min{j; j > U (n -1), X j > X U (n-1) } for n ≥ 2. We define the i-th upper record value by R i = X U (i) , with R 1 = X 1 . Using the pdf f (x) (1) and the cdf F (x) (3), the pdf of R i is given by

f Ri (x) = [-log(1 -F (x))] i-1 (i -1)! f (x) = 1 (i -1)! [ -log ( log(1 + βe -λx ) + βe -λx log(β + e λx ) (β + 1) log(β + 1
)

)] i-1 × λβ (β + 1) log(β + 1)
log(β + e λx )e -λx , x > 0.

Using pdf (1) and the hrf (4), the joint pdf of (R 1 , . . . , R n ) is given by

f (R1,...,Rn) (x 1 , . . . , x n ) = f (x n ) n-1 ∏ k=1 h(x k ) = λβ (β + 1) log(β + 1) log(β + e λxn )e -λxn × n-1 ∏ k=1 λβ log(β + e λx k )e -λx k log(1 + βe -λx k ) + βe -λx k log(β + e λx k ) , 0 < x 1 < . . . < x n .
For i < j, the joint pdf of (R i , R j ) is given by

f (Ri,Rj ) (x i , x j ) = [-log(1 -F (x i ))] i-1 (i -1)! [ log ( 1-F (xi) 1-F (xj ) )] j-i-1 (j -i -1)! h(x i )f (x j ) = 1 (i -1)! [ -log ( log(1 + βe -λxi ) + βe -λxi log(β + e λxi ) (β + 1) log(β + 1) )] i-1 × 1 (j -i -1)! [ log ( log(1 + βe -λxi ) + βe -λxi log(β + e λxi ) log(1 + βe -λxj ) + βe -λxj log(β + e λxj ) )] j-i-1
× λβ log(β + e λxi )e -λxi log(1 + βe -λxi ) + βe -λxi log(β + e λxi ) × λβ (β + 1) log(β + 1) log(β + e λxj )e -λxj , 0 < x i < x j .

Estimation procedures

This section is devoted to three estimation procedures for the LWE model parameters: the maximum likelihood estimation, weighted least square estimation and maximum product of spacings.

Maximum likelihood estimators

Let X 1 , X 2 , . . . , X n be a random sample from the LWE distribution with parameter vector Θ = (β, λ) and x 1 , x 2 , . . . , x n are the corresponding observed values. By considering f (x) (1), the likelihood function is given by

L(Θ) = n ∏ i=1 f (x i ) = λ n β n (β + 1) n (log(β + 1)) n [ n ∏ i=1 log(β + e λxi ) ] [ n ∏ i=1 e -λxi
] .

The log-likelihood function can be expressed as

ℓ(Θ) = n log(λ) + n log(β) -n log(β + 1) -n log(log(β + 1)) + n ∑ i=1 log(log(β + e λxi )) -λ n ∑ i=1 x i .
The nonlinear log-likelihood equations given by ∂ℓ(Θ) ∂Θ = 0 are listed below

∂ℓ(Θ) ∂β = n β - n β + 1 - n (β + 1) log(β + 1) + n ∑ i=1 1 (β + e λxi ) log(β + e λxi ) = 0 and ∂ℓ(Θ) ∂λ = n λ + n ∑ i=1 x i e λxi (β + e λxi ) log(β + e λxi ) - n ∑ i=1 x i = 0.
Solving the equations above simultaneously, we obtain the maximum likelihood estimators (MLEs) of the model parameters. Some numerical iterative methods may be used to estimate the model parameters and the global maxima of the log-likelihood can be justified by considering different starting values for the parameters. To show the likelihood equations have a unique solution in the parameters, we plot the profile log-likelihood function of the model parameters in Figure 3, for one of the data sets (D1) used in the application section.

Using the asymptotic distribution of the MLEs, the information matrix is used to establish the confidence intervals for the model parameters. The elements of the 2 × 2 observed information matrix J(Θ) = Jrs(Θ) for r, s ∈ {β, λ} can be obtained from the authors upon request.

Weighted least-square estimators

Let x 1 , x 2 , . . . , x n be an ordered sample of the random sample of size n from the LWE distribution. Let F (x) be the cdf (3). Then the weighted least square estimators (WLSE) can be obtained by minimizing

S(Θ) = n ∑ i=1 w i {F (x i ) -E(F (X i:n ))} 2 = n ∑ i=1 w i { 1 - log(1 + βe -λxi ) + βe -λxi log(β + e λxi ) (β + 1) log(β + 1) - i n + 1 } 2 ,
with respect to the unknown parameters of the LWE distribution, where

w i = 1 V ar(F (X i:n )) = (n + 1) 2 (n + 2) i(n -i + 1
) .

The associated nonlinear equations ∂S(Θ) ∂Θ = 0 are given by

∂S(Θ) ∂β = 2 n ∑ i=1 w i { 1 - log(1 + βe -λxi ) + βe -λxi log(β + e λxi ) (β + 1) log(β + 1) - i n + 1 } η 1 (x i , Θ) = 0 and ∂S(Θ) ∂λ = 2 n ∑ i=1 w i { 1 - log(1 + βe -λxi ) + βe -λxi log(β + e λxi ) (β + 1) log(β + 1) - i n + 1 } η 2 (x i , Θ) = 0,
where

η 1 (x i , Θ) = - 1 (β + 1) 2 (log(β + 1)) 2 ( log(β + 1) + β log(β + 1) -β log(β + e λxi ) + log(β + 1) log(β + e λxi ) -e λxi (1 + log(β + 1)) log(1 + βe -λxi ) ) e -λxi (9)
and

η 2 (x i , Θ) = β (β + 1) log(β + 1)
x i e -λxi log(β + e λxi ). [START_REF] Mahdavi | Two Weighted Distributions Generated by Exponential Distribution[END_REF] Solving this system of nonlinear equations simultaneously will yield the WLSE's of the distribution parameters.

Maximum product of spacings estimators

The maximum product of spacings (MPS) method is a powerful alternative to MLE for the estimation of the unknown parameters of continuous univariate distributions. Let x 1 , x 2 , . . . , x n be an ordered sample of the random sample of size n from the LWE distribution. Let F (x) be the cdf (3) and D i (Θ) = F (x i ) -F (x i-1 ) for i = 1, 2, . . . , n + 1, be the uniform spacings of a random sample from the LWE distribution, where F (x 0 ) = 0 and F (x n+1 ) = 1. The maximum product spacings estimators ΘMP S are obtained by maximizing the geometric mean of the spacings given by

G (Θ) = [ n+1 ∏ i=1 D i (Θ) ] 1 n+1
, with respect to Θ. It is equivalent to maximize the logarithm of the geometric mean of sample spacings given by

ℓ (Θ) = log (G (Θ)) = 1 n + 1 n+1 ∑ i=1 log (D i (Θ)) .
Using η 1 (x i , Θ) (9) and η 2 (x i , Θ) [START_REF] Mahdavi | Two Weighted Distributions Generated by Exponential Distribution[END_REF], the estimates ΘMP S can be obtained by solving the nonlinear equations

∂ℓ(Θ) ∂β = 1 n + 1 n+1 ∑ i=1 1 D i (Θ) [η 1 (x i , Θ) -η 1 (x i-1 , Θ)] = 0 and ∂ℓ(Θ) ∂λ = 1 n + 1 n+1 ∑ i=1 1 D i (Θ) [η 2 (x i , Θ) -η 2 (x i-1 , Θ)] = 0.

Applications to real data

In this section, the potentiality of the LWE distribution is highlighted. We show that our distribution has a better fit than the compared distributions under all estimation procedures: MLE, WLSE and MPSE, for four real data sets.

We fit the two-parameter gamma (gamma), exponential logarithmic (EL) [START_REF] Tahmasbi | A two-parameter lifetime distribution with decreasing failure rate[END_REF], exponentiated exponential (EE) [START_REF] Gupta | Exponentiated Exponential Family; An Alternative to Gamma and Weibull[END_REF], weighted exponential (WE) [START_REF] Gupta | A new class of weighted exponential distribution[END_REF] and LWE distributions to four real data sets. The pdf of these models are as follows.

• The pdf of the gamma distribution is given by

f (x) = x k-1 e -x θ θ k Γ(k) , x, k, θ > 0 .
• The pdf of the exponential logarithmic distribution is given by

f (x) = - β (1 -p) e -βx log(p) (1 -(1 -p) e -βx ) , x, β > 0, p ∈ (0, 1) .
• The pdf of the exponentiated exponential distribution is given by

f (x) = α λ e -λ x ( 1 -e -λ x ) α-1 , x, α, λ > 0 .
• The pdf of the weighted exponential distribution is given by

f (x) = (α + 1) λ e -λx ( 1 -e -α λx ) α ,
x, α, λ > 0 . We check the adequacy of the fitted models via the statistics of Anderson-Darling and the Cramérvon Mises. They allow to determine how closely a specific distribution fits the associated empirical distribution for a given data set. The smaller statistics give the better fit.

The description of the considered real data sets is as follows. Time to failure in hours for non-repairable Item data Set (D1). The first data set we consider represents the time to failure in hours for non-repairable items. The data set consists of 10 sample points and is reported at [START_REF] Murthy | Weibull models[END_REF]Page 279].

Failure times data Set (D2). The second data set we consider represents the failure times of 20 components. The data set can be found in [START_REF] Murthy | Weibull models[END_REF]Page 154].

Total monthly rainfall data Set (D3). The third data set is obtained from the department of water resources and power agency manager of water resources of the State of Sao Paulo. The data represents total monthly rainfall during April at Sao Carlos from 1960 to 2014. The data set can be found at [START_REF] Bakouch | Binomial-exponential 2 distribution: Different estimation methods with weather applications[END_REF].

Fatigue fracture data Set (D4). The fourth data represents the life of fatigue fracture of Kevlar 373/epoxy subjected to constant pressure at 90% stress level until all had failed. For previous studies with the data see [START_REF] Barlow | A Bayesian analysis of stress rupture life of Kevlar 49/epoxy spherical pressure vessels[END_REF].

Table 1 provides descriptive statistics of the data sets. From this table, it can be seen that all the used data sets are right skewed. Some moments of the LWE distribution are given in Table 2. We computed the results using the definitions, where the parameters of the LWE are replaced by MLEs for D1 and WLSEs for the other three data sets respectively from Tables 3456. From Tables 1 and2, it can be noted that the considered moments of the LWE distribution are approximately close to the sample moments for all the four data sets.

Tables 3-6 list the parameter estimates for the three estimation procedures discussed in Section 6, for the four data sets with goodness-of-fit statistics. From these tables, it can be noted that LWE distribution is a powerful competitor to all the compared distributions. It provides a rather flexible mechanism for fitting a wide spectrum of positive real data sets with shape property of being right skewed. It is also observed that for small sample size MLE performs better than the WLSE and MPSE but for larger sample size WLSE performs better than MLE and MPSE.

The conclusions drawn above are confirmed by Figure 4; it can be seen that the estimated LWE density superimposed on the histogram of the data sets based on MLEs for D1 and WLSEs for D2-D4 provides a good fit to these data.

Figure 1 :

 1 Figure 1: Plots of the LWE pdf.

Figure 2 :

 2 Figure 2: Plots of the LWE hrf.

Figure 3 :

 3 Figure 3: The profile of log-likelihood function for (λ, β).

Figure 4 :

 4 Figure 4: The fitted LWE pdfs superimposed on the histogram of data.

Table 1 :

 1 Descriptive statistics of the data sets

					D1
	n	Mean	Median	SD	Skewness	Kurtosis	M1	M2
	10	100.1	71.6	90.5638	0.943922	2.86235	70.4	47.8
					D2
	20	8.42945	8.662	5.32206	0.176969	2.43092	4.2029	3.663
					D3
	53	80.883	76	52.2441	0.494003	2.51167	43.2915	41.1
					D4
	76	1.95924	1.73615	1.57398	1.97956	8.16079	1.07026	0.7202
	SD = Standard Deviation , M1 = Mean deviation about the mean,
	M2 = Mean deviation about the median

Table 2 :

 2 Some moments of the LWE

	D1					D2
	Mean	Median	SD	M1	M2	Mean	Median	SD	M1	M2
	99.9137	74.4296	91.026	69.3594	66.3962	9.13596	7.51929	7.02205	5.38507	5.22879
	D3					D4
	Mean	Median	SD	M1	M2	Mean	Median	SD	M1	M2
	84.2705	69.9799	62.7895	48.1234	46.7521	1.89683	0.593826	1.35356	1.03643	1.00689

Table 3 :

 3 Comparison of fit of LWE using different methods of estimation for Time to Failure in Hours for Nonrepairable Item data

			MLEs
	Distributions	Estimates		A *	W *
	Gamma(k, θ) 1.139028	87.881928	0.144832	0.016461
	EL(p, β)	0.999999	0.009990	0.146562	0.017179
	WE(λ, α)	0.010233	39.953472	0.149407	0.017899
	EE(λ, α)	0.010826	1.136636	0.143582	0.016449
	LWE(λ, β)	0.015096	1.451571	0.128042 0.014967
			WLSEs
	Distributions	Estimates		A *	W *
	Gamma(k, θ) 0.796078	147.040746	0.187528	0.021962
	EL(p, β)	0.346263	0.006459	0.176709	0.020845
	WE(λ, α)	0.009125	1549.999999 0.129782	0.0148178
	EE(λ, α)	0.007323	0.789169	0.187466	0.021972
	LWE(λ, β)	0.010721	7.699068	0.134756 0.015240
			MPSEs
	Distributions	Estimates		A *	W *
	Gamma(k, θ) 0.788220	146.918918	0.189559	0.022436
	EL(p, β)	0.353425	0.006614	0.172185	0.020647
	WE(λ, α)	0.008841	499.999279	0.146562	0.017179
	EE(λ, α)	0.007347	0.781233	0.189547	0.022428
	LWE(λ, β)	0.010536	7.490030	0.143201 0.016290

Table 4 :

 4 Comparison of fit of LWE using different methods of estimation for failure times of 20 components data

			MLEs
	Distributions	Estimates	A *	W *
	Gamma(k, θ) 1.276065	6.605812 0.964807	0.173748
	EL(p, β)	0.999999	0.118626 1.20985	0.240738
	WE(λ, α)	0.237238	0.000210 1.06544	0.116016
	EE(λ, α)	0.134624	1.235388 0.997362	0.182967
	LWE(λ, β)	0.211521	0.411468 0.707688 0.126298
			WLSEs
	Distributions	Estimates	A *	W *
	Gamma(k, θ) 1.872033	4.914790 0.925161	0.076585
	EL(p, β)	0.999999	0.104051 0.948973	0.161832
	WE(λ, α)	0.217647	0.000038 1.002540	0.074213
	EE(λ, α)	0.158061	1.918960 0.963278	0.081900
	LWE(λ, β)	0.202654	0.238868 0.605758 0.076161
			MPSEs
	Distributions	Estimates	A *	W *
	Gamma(k, θ) 0.994138	8.962964 1.06508	0.200422
	EL(p, β)	0.999999	0.112255 1.057910	0.198708
	WE(λ, α)	0.228869	0.000494 1.004650	0.092105
	EE(λ, α)	0.109947	0.968710 1.096200	0.207921
	LWE(λ, β)	0.195684	0.587401 0.662525 0.111726

Table 5 :

 5 Comparison of fit of LWE using different methods of estimation for Rainfall data

			MLEs
	Distributions	Estimates		A *	W *
	Gamma(k, θ) 1.539744	52.530164 0.877195	0.130219
	EL(p, β)	0.999999	0.012363	2.36462	0.425471
	WE(λ, α)	0.024726	0.000042	0.820019	0.092466
	EE(λ, α)	0.015789	1.514078	0.973941	0.147181
	LWE(λ, β)	0.023096	0.213791	0.596077 0.090642
			WLSEs
	Distributions	Estimates		A *	W *
	Gamma(k, θ) 1.823373	46.318640 0.665771	0.060652
	EL(p, β)	0.999999	0.011138	1.83096	0.2781
	WE(λ, α)	0.023778	0.000020	0.758456	0.064999
	EE(λ, α)	0.017183	1.901809	0.727332	0.067376
	LWE(λ, β)	0.022642	0.137041	0.483161 0.053447
			MPSEs
	Distributions	Estimates		A *	W *
	Gamma(k, θ) 1.379133	59.777856 1.00204	0.149678
	EL(p, β)	0.999999	0.012082	2.19852	0.382795
	WE(λ, α)	0.024345	0.000001	0.777771	0.0783983
	EE(λ, α)	0.014531	1.353717	1.10728	0.168091
	LWE(λ, β)	0.022328	0.277365	0.5744	0.080558

Table 6 :

 6 Comparison of fit of LWE using different methods of estimation for fatigue data

			MLEs
	Distributions	Estimates	A *	W *
	Gamma(k, θ) 1.641076	1.193875 0.680647	0.113341
	EL(p, β)	0.999999	0.510401 3.01876	0.574565
	WE(λ, α)	0.711243	1.541313 0.580135	0.0836306
	EE(λ, α)	0.702793	1.709493 0.677809	0.113052
	LWE(λ, β)	0.970107	0.150408 0.59646	0.100525
			WLSEs
	Distributions	Estimates	A *	W *
	Gamma(k, θ) 2.003304	0.944853 0.583048	0.063036
	EL(p, β)	0.999999	0.456246 2.60599	0.456466
	WE(λ, α)	0.847727	0.613244 0.557909	0.067250
	EE(λ, α)	0.815028	2.150388 0.569838	0.062684
	LWE(λ, β)	1.046110	0.021336 0.542988 0.066131
			MPSEs
	Distributions	Estimates	A *	W *
	Gamma(k, θ) 1.476364	1.358301 0.906648	0.160206
	EL(p, β)	0.999998	0.496530 2.83668	0.527345
	WE(λ, α)	0.650868	2.265658 0.659052	0.101915
	EE(λ, α)	0.635134	1.658918 0.844588	0.15135
	LWE(λ, β)	0.930942	0.217076 0.727137 0.1306
			14


Petr Velan
email: petr.velan@cesnet.cz

Radek Krejčí

Flow Information Storage Assessment Using IPFIXcol

Keywords: collector, flow, IPFIX, NetFlow, network monitoring

Network monitoring has became a significant part of network management. Each environment and type of network have their specific, different needs. To allow network traffic monitoring in various environments, a necessity of flexible approach thus grows. The current generation of flow collectors provides only a limited flexibility, mainly due to limits of their data storage formats. Moreover, it is quite a challenging task to compare particular storage formats and their suitability for the specific environment and usage. In this paper we present IPFIXcol -a flow collector framework designed for easy data storage formats changing. This way, we plan to evaluate performance and suitability of various data storage formats for specific tasks. Results can be used to build the most appropriate data storage for the specific production environments.

Introduction

Flow traffic monitoring becomes widespread and demands for its additional features grow. Sometimes, we significantly miss some specific information in the stored flow records. The currently used flow records relate to IP networks only. However, flow information can be beneficial in non-IP networks too, like building automation and control system (BACS) networks or supervisory control and data acquisition (SCADA) networks. The additional information fields not included in standard flow records are needed for applications such as classification of network traffic based on time characteristics [START_REF] Piskač | Using of Time Characteristics in Data Flow for Traffic Classification[END_REF]. These needs can be satisfied by metering and collecting processes supporting the IP Flow Information Export (IPFIX) protocol [START_REF] Claise | Specification of the IP Flow Information Export (IPFIX) Protocol for the Exchange of IP Traffic Flow Information[END_REF].

In this paper, we are introducing IPFIXcol -IPFIX protocol collecting framework for receiving, processing and storing various flexible flow information. We are using this framework in our research to store the same flow information data in multiple formats. This way, the efficiency and suitability of each format can be compared. Based on the results, we plan to optimise current formats or to propose new ones to store network traffic flow information for various scenarios.

The paper is organized as follows. After the Introduction, some of the currently used IPFIX implementations and related works are summarized in Section 2. The IPFIXcol architecture is described in Section 3. Future use of the framework for comparing and designing flow information storage formats is discussed in Section 4.

Related Work

A list of exporters supporting IPFIX protocol includes, for example, the Flow-Mon probe1 , nProbe2 , Vermont (VERsatile MONitoring Toolkit) 3 or YAF (Yet Another Flowmeter) 4 . IPFIX collector implementations are provided by nTop5 toolkit, which comes from the same authors as nProbe, and the above-mentioned Vermont toolkit. There is also the python library called ripfix6 that provides support for building IPFIX collecting and exporting applications. Results of the IPFIX implementations interoperability testing can be found in [START_REF] Trammell | DEMONS IPFIX Interoperability Event --Final Report[END_REF].

Our research is inspired by the comparison of nfdump7 flat file format and MySQL database [START_REF] Hofstede | The Network Data Handling War: MySQL vs. NfDump[END_REF]. We will take the next step and provide a general comparison of more types of data storage formats along with complex benchmarks.

The research of some new possibilities in storing of streamed network data, such as network traffic flow information, is done by Fusco [START_REF] Fusco | Real-time creation of bitmap indexes on streaming network data[END_REF]. IPFIXcol provides framework for testing and results deployment of such a research while using the real network traffic data.

IPFIXcol System Architecture

The IPFIXcol is composed of a core program which is further extended by plugins (see Figure 1). Incoming flow data pass through an input plugin associated with a transport protocol. The output of the input plugin is parsed by the core to get some general information for further use by storage plugins. If specified, data can be modified by internal plugins. Storage plugins then store flow records in a specific data format. A more detailed description of each plugin type follows.

Input Plugins get the data from different flow information sources such as local files or network flow record streams from routers or standalone probes. The currently implemented input plugins provide support for an IPFIX file format [START_REF] Trammell | Specification of the IP Flow Information Export (IPFIX) File Format[END_REF] and IPFIX data transferred over SCTP, UDP, TCP or TLS over TCP protocols.

Internal Plugins extend functionality of the IPFIXcol core. With internal plugins, the IPFIXcol is able to serve as an IPFIX mediator [START_REF] Kobayashi | IP Flow Information Export (IP-FIX) Mediation: Framework[END_REF] that can modify incoming data. Internal plugins can be used i.e. for data anonymization or flowbased sampling. Storage Plugins are used to process received data. This usually includes storing the flow records into files or databases. Furthermore, storage plugins can also resend the data to other destinations or process them to get overall statistics. Configuration of the plugins allows specifying rather complex cross connections between the flow information sources and the data stores.

Currently, the IPFIXcol storage plugins allow storing data in the IPFIX flat file format, PostgreSQL database and FastBit column-based database [START_REF]FastBit[END_REF].

Conclusion and Future Work

This paper has introduced the IPFIXcol flow collecting framework. The IPFIXcol is already in a testing phase, including several basic input and output plugins.

The IPFIXcol is not just another IPFIX collector. It provides a framework for testing various procedures of flow information pre-processing or storage data formats. The IPFIXcol with an effective storage data format can be used for monitoring high-speed networks. By changing the storage data format or an internal processing plugin to support specific information fields, the IPFIXcol can be effectively used for monitoring research and non-standard environments such as building management or SCADA networks [START_REF] Krejčí | Traffic Measurement and Analysis of Building Automation and Control Networks[END_REF].

Future Work

To measure and compare characteristics of various flow data storage formats, we need to establish a set of test scenarios. We will take several queries typical of flow information data, such as data selection based on specified properties, data aggregation or combinations of both, so that the stored data are accessed sequentially or randomly. Then we will transform these queries into a specific format used by the data storage management tools.

New IPFIXcol storage plugins will be used to store the same data in different formats, so that our test queries return the same results. This way we can also measure the throughput of different database solutions, which can be interesting in high-speed networks.

With these tests, we will be able to determine which type of data storage is most suitable for a specific environment and an application. We can also compare widely used nfdump flat file format, which is limited to the NetFlow v9, with more flexible databases that support IPFIX.

The acquired results will allow us to select and eventually improve the primary IPFIXcol storage format and its query tool. Currently, the usage of the FastBit library for storing network traffic flow information is very promising. The library is being developed at Lawrence Berkeley National Laboratory and we can contribute to it with specific demands from our use cases. Based on the test results, we might need to develop our own specific purpose database dedicated to storing flexible flow data.

Fig. 1 .

 1 Fig. 1. System Architecture

http://www.invea-tech.com/products-and-services/flowmon/flowmon-probes

http://www.ntop.org/nprobe/nprobe-complies-with-ipfix-specification/

http://vermont.berlios.de/

http://tools.netsa.cert.org/yaf/

http://www.ntop.org/

http://ripfix.rubyforge.org/

http://nfdump.sourceforge.net/

Acknowledgement This work is supported by the "CESNET Large Infrastructure" project LM2010005 funded by the Ministry of Education, Youth and Sports of the Czech Republic.