Mohamed Rabie Oularbi

Abdeldjalil Aissa-El-Bey

Sebastien Houcke

Physical Layer IEEE 802.11 Channel Occupancy Rate Estimation

Keywords: Vertical handover, channel occupancy rate, IEEE 802.11

With the rapid growth of wireless communication networks, we are facing the challenge of integration of diverse wireless networks such as WLAN and WWAN. Therefore, it becomes important to think at vertical handoff solution where the user can move seamlessly among various type of networks, to provide the best QoS to the higher layers applications. In this paper, we propose a method to estimate the channel occupancy rate metric of an IEEE 802.11 network from the physical layer. This method can be applied for any network that use the CSMA/CA protocol. Our theorical results are validated using experimental measurements captured by RAMMUS RF Plateform.

I. INTRODUCTION

Wireless communication devices become more popular and pervasiveness. These devices are used for a wide range of multimedia applications, which continually require a high QoS (quality of service) to achieve the expected performance by the user.

In homogeneous cellular systems, the link quality is keeped continuously by searching the best neighboring cell to associate to. The signal to noise ratio is then used as a common metric. When the base station giving the best link quality is found, a horizontal handoff occurs and the connection is transfered between the two base stations.

When the number of users connected on a given networks increases, the effective available bandwidth is affected. To maintain the QoS required by the higher layer application, and benefiting from the coexistence of various wireless networks, the users should roam freely from one interface to an other without any disturbing. This important process is referred as vertical handover.

The vertical handover is a very important capability in the future wireless communication era, where an integrated network including multiple technologies will try to offer a global broadband access to mobile users. However, compared to the horizontal handover, the signal strength metric is sometimes not suited and often not sufficient to appropriately trigger the vertical handover : as heterogeneous networks have different system characteristics, their performance cannot be simply compared by using the signal strength of two cells [START_REF] Dai | Vertical handover criteria and algorithm in IEEE 802.11 and 820.16 hybrid networks[END_REF].

Within this framework, this paper aims at defining an efficient estimation method for a metric that can informs us on the channel occupancy rate of an IEEE 802.11 network. The originality of this method relies on the fact that this metric is estimated from the physical layer instead of the MAC layer as proposed in [START_REF] Dai | Vertical handover criteria and algorithm in IEEE 802.11 and 820.16 hybrid networks[END_REF], [START_REF] Guo | A seamless and proactive end-to-end mobility solution for roaming across heterogeneous wireless networks[END_REF], which represent a large gain of time and complexity.

In [START_REF] Dai | Vertical handover criteria and algorithm in IEEE 802.11 and 820.16 hybrid networks[END_REF], [START_REF] Guo | A seamless and proactive end-to-end mobility solution for roaming across heterogeneous wireless networks[END_REF] it has been highlighted that the usage of the channel bandwidth in a WiFi system can be approximated as the ratio between the time in which the channel status is busy according to the NAV (Network Allocation Vector) settings and the considered time interval. Indeed, prior to transmitting a frame, a station calculates the amount of time necessary to send the frame based on the frame's length and data rate. The station places a value representing this time in the duration field in the header of the frame. From the above description we can see that the NAV busy state can well reflect the traffic load. The higher the traffic, the larger the NAV busy occupation, and vice versa. Therefore, if we observe a NAV value during a certain time window, the available bandwidth and access delay can be estimated given a certain packet length [START_REF] Zhang | Efficient Mobility Management for Vertical Handoff between WWAN and WLAN[END_REF].

The matter with this method is that it requires to be connected to the access point in order to have access to the NAV duration from the header, this may increases the decision time if many standards or Access Point (AP) are detected. In this paper, we propose a method that requires no connection to the AP, and no NAV duration reading. This technique is based on a physical layer sensing : Considering that the medium is free when only noise is observed and occupied when signal plus noise samples are observed (data frame), we use a likelihood function that can distinguish the signal plus noise samples from the one corresponding to noise only. Once we get the number of signal plus noise samples, a simple ratio processing can inform us on the network occupancy rate.

II. MODEL STRUCTURE

In the rest of the paper, we assume that IEEE 802.11 access points are detected. An IEEE 802.11 communication is based on a collision avoidance medium access protocol. Between two consecutive frames we have different inter frame spacing (IFS) intervals which guarantee different types of priority. At the receiver side, the observed signal is a succession of frames of noise samples corresponding to the IFS intervals or idle periods and of data frames (signal plus noise).

For clarity reason, we assume in this section that we have only one data frame in the observation duration (Ns samples) and explain in section III the proposed algorithm to locate it.

Consider that our receiver is doted of N antennas and let y i = [yi(1), . . . , yi(Ns)] be a set of Ns observations on the i th antenna such that

8 < : yi(n) = wi(n) 1 < n < n1 -1 yi(n) = xi(n) n1 < n < n2 yi(n) = wi(n) n2 + 1 < n < Ns (1)
where xi(n) is the based band sample being received on the i th antenna at the instant n, expressed as

xi(n) = L-1 X k=0 hi(k)s(n -n1 -k) + wi(n) (2)
where s(n) denotes the n th transmitted symbol, with

E[|s(n)| 2 ] = σ 2 s . hi(k)
is the channel response between the source signal and the i th antenna. L is the order of the channel, and

P L-1 0 σ 2 h i (k) = 1.
wi(n) is a complex additive white gaussian noise with zero mean and variance σ 2 w .

III. FRAME LOCALISATION

As presented in the previous section, the vector y i can be divided into three parts : noise , signal plus noise and noise. Starting from the set of observation y i we want to find which samples correspond to noise and which ones correspond to signal plus noise. Since the samples are supposed to be independent in the noise areas, and correlated in the signal plus noise area we propose to use a likelihood function that informs us on the independance of the processed sample.

Let now Yi(u) denotes the following set of observations

Yi(u) = [yi(u), . . . , yi(Ns)] 1 ≤ u < Ns (3)
And let us define fY the joint probability density function of Yi(u). If Yi(u) is composed of only noise samples : fY (Yi(u)) = Q Ns m=u fw(yi(m)), where fw is the probability density function of a complex normal law centered of variance σ 2 w . The variance σ 2 w is assumed to be known or at least estimated by a subspace-based algorithm [START_REF] Xu | Subspace-Based Noise Variance and SNR Estimation for OFDM Systems[END_REF].

The log-likelihood that the vector Yi(u) is formed of (Nsu) noise independent samples is expressed as

Li(u) = log " Ns Y m=u fw(yi(m)) # (4)
Computing the mean of the N log-likelihood functions expressed on each sensor, we get a criterion J (u) that informs us on the independence of the processed samples

J (u) = 1 N N X i=1 Li(u) = -(Ns -u) log(πσ 2 w ) - 1 N σ 2 w N X i=1 Ns X m=u |yi(m)| 2 (5)
As u varies in the interval [1, n1), the number of noise samples composing Yi(u) decreases and so does J (u) until it reaches a minimum bound at n1. However, for u varying from n1 to n2 the number of signal plus noise samples decreases, therefore the ratio noise samples over signal plus noise samples increases and by the way J (u) increases. It reaches its maximum value while Yi(u) contains only noise samples, i.e when u = n2.

Finally for n2 < u < Ns, J (u) decreases again for the sames reasons than the one explained for 1 < u < n1.

IV. ESTIMATION OF THE CHANNEL OCCUPANCY RATE

We propose to get the channel occupancy rate by a physical layer sensing. Indeed, while observing a set of Ns samples, if we can estimate the number of samples corresponding to signal plus noise (i.e the length of the data frame), we can easily estimate the channel occupancy rate.

When we have only one data frame in the observed window the occupancy rate can easily be estimated thanks to the previous criterion by n2 -n 1

Ns . However, the assumption to have only one frame in the duration window is too restrictive. In practice we may get a signal as shown in figure 1 or with more frames. Based on the behavior of J (u), we can clearly see (fig 1) that the slope of J (u) is positive when u corresponds to the index of a signal plus noise sample and negative when u corresponds to the index of a noise sample. Therefore, we can take advantage of the gradient of J (u) to distinguish the nature of our observed samples. Introducing the function Φ(u) such that Here we denote by ∇ the gradient of J (u) processed using the central difference method, such that the derivative for any point of index u / ∈ {1, Ns} is processed as ∇(J (u)) = 1 2 (J (u + 1) -J (u -1)). For the first point, we use the forward finite difference, and at the left end element the backward difference is used. sign{.} denotes the sign operator. Therefore, Φ(u) equals to 1 when signal plus noise samples are present and zero when it is only noise, and the channel occupency rate is estimated by

Φ(u) = 1 2 [sign{∇(J (u))} + 1] (6)
d Cor = 1 Ns Ns X u=1 Φ(u) (7)
The difficulty is to estimate the channel occupency rate accurately for low signal to noise ratio. In fact, there are fluctuations that can mislead the decision for a given sample. To fix this problem, we propose to use a smoothing technique.

V. SMOOTHING METHODS

We can either smooth the criterion J (u) or the function Φ(u). The two smoothing methods are presented in the next subsections.

A. Criterion J (u) smoothing

As said previously, in the signal plus noise frame there are fluctuations that mislead the decision for some samples. We propose here to smooth J (u) thanks to a polynomial interpolation of first order. Given a certain length of the smoothing window we search for the extrema of J (u) in this window and join them by a line, once it is done we slide the window with no overlapping until all the samples are treated. The procedure is illustrated in algorithm 1. The choice of the length of the smoothing window W is very important. Indeed, the algorithm do not behave correctly if we have multiple maximum or minimum of J (u) in the same window. To avoid such situation, we choose W equal to the length of a SIFS (for Short IFS) which is the smallest inter frame gap. Thus, theoretically we can't get a set of successive noise samples of a length less than a SIFS, and by the way we avoid multiple extrema in the same window.

B. Function Φ(u) smoothing

Here we propose to smooth Φ(u) directly. As presented below, theoretically we can't get a set of successive noise samples of a length less than a SIFS. Then, if it happens, it means that the algorithm took the wrong decision and Φ(u) will be forced to 1 for those samples. Practically, to avoid confusion it is judicious to choose a smoothing window less than a SIFS.

The matter of fluctuations, drives us to search the conditions which makes the criterion an increasing function for signal plus noise samples and decreasing for noise samples. In the next section, we show that is possible only for a given range of the noise variance.

VI. CRITERION VALIDATION LIMITS

In this section, we propose to investigate the limits of the proposed criterion J (u). The aim is to find the dynamic where J (u) do well behave, i.e where its slope is positive for signal plus noise samples and negative for noise samples.

For

1 ≤ u ≤ n1 : J (u) decreases only if ∂E[J (u)]
∂u < 0, and therfore if

E[J (u)] = -(Ns -u) log(πσ 2 w ) - 1 σ 2 w [(n1 -u)σ 2 w + (n2 -n1)(σ 2 w + σ 2 s ) + (Ns -n2)σ 2 w ]
the derivative costs : ∂E[J (u)] ∂u = log(πσ 2 w ) + 1, and we get σ 2 w < 1 πe . For n1 ≤ u ≤ n2 : J (u) is an increasing function only if

∂E[J (u)] ∂u > 0, then if E[J (u)] = -(Ns -u) log(πσ 2 w ) - 1 σ 2 w [(n2 -u)(σ 2 w + σ 2 s ) + (Ns -n2)σ 2 w ]
the partial derivative is :

∂E[J (u)] ∂u = log(πσ 2 w ) + 1 σ 2 w (σ 2 w + σ 2 s ) (8)
and

J (u) increases only if σ 2 w > 1 πe (1+γ) , where γ = σ 2 s σ 2
w . For n2 ≤ u ≤ Ns : we get the same result as for 1 ≤ u ≤ n1. As a conclusion for an optimal behavior of J (u) , the noise variance must satisfy

1 πe (1+γ) < σ 2 w < 1 πe (9)
This inequality represents the limits of the proposed criterion. It means that, the performance of the proposed method depends on the noise variance value, and also on the signal to noise ratio. Therefore, if the noise variance do not satisfy (9), we can think to adjust it applying a certain gain on the received signal. This alternative requires to know the signal to noise ratio which is not always possible. Another approach is to introduce a new criterion that get around this matter, this criterion is the distance between J (u), and a Parzen estimator based criterion introduced in the next section.

VII. PARZEN ESTIMATOR BASED CRITERION

The proposed solution, consists in processing a new criterion that aims to minimize the distance between the true p.d.f of the noise and a Parzen estimator estimated p.d.f of the observed samples. Starting from the set of observations y i for i = 1, . . . , Ns and deviding the samples into their real part pi(n) and imaginary part qi(n), we get 2N Ns samples available for estimating the Parzen window density distribution. Given a sample yi(n) = pi(n) + j.qi(n) its Parzen window distribution will be f (yi(n)) = f (pi(n)). f (qi(n)) where

f (z) = 1 2N NsF 2NNs -1 X k=0 K " z -z k F " (10)
K is the Parzen window kernel and F is a smoothing parameter called the bandwidth. This kernel has to be a suitable p.d.f function, we use Gaussian kernels with standard deviation one. The new processed criterion is :

JK (u) = 1 N N X i=1 log " Ns Y m=u f (yi(m)) # (11)
Once we get JK (u), we measure the distance between J (u) and JK (u) to obtain a new criterion

R(u) = |J (u) -JK (u)| (12)
Then, R(u) is smoothed thanks to the proposed method in V-A, a new function ΦK (u) is processed thanks to the equation 6, substituting J (u) by R(u). The obtained ΦK (u) is also smoothed using the technique detailed in V-B, then the channel occupancy rate is processed thanks to the equation [START_REF]Ettus Research LLC website[END_REF]. As treated previously, the Channel occupancy rate is function of the behavior of J (u). In figure 2, we show the NMSE (Normalized Mean Square Error) of the estimation of the channel occupancy rate versus the SNR, the SNR is defined as SNR =

σ 2 s σ 2 w P L-1 l=0 σ 2 h i (l)
and is assumed to be constant on each sensor. The results are averaged over M = 1000 Monte Carlo runs, and the NMSE is here defined as

1 M P M 1 " d Cor k -Cor " 2 /Cor 2 ,

where d

Cor k is the channel occupancy rate estimated at the k th realization and Cor is the real channel occupency rate. We can clearly observe that for a high SNR the error tends to zero, and thus we achieve a good estimation. The proposed method is compared to the CFAR (Constant False Alarm Rate) method with a probability of false alarm P f a = 10 -3 and to the energy detector proposed by Urkowitz [START_REF] Urkowitz | Energy Detection of Unknown Dterministic Signals[END_REF], with a P f a = 10 -4 . The cognitive terminal is supposed to be doted of N = 2 antennas. The smoothed Φ(u) has the best performace followed by the Parzen estimator for low SNR and finally the smoothed J (u). Concerning R(u), it has the advantage of being independent of the noise variance value, its performance are the same for any value of σ 2 w even if its does not satisfy (9). The proposed blind detection approach is evaluated using RAM-MUS RF Plateform. Experiments were invistigated on the Channel 3 (2.422 Ghz) using the IEEE 802.11g norm. We tested differents schemes with diffent number of users varying the maximum bit-rate allocated to each one. The physical layer signal was captured thanks to an USRP2 device (Universal Software Radio Peripheral [START_REF]Ettus Research LLC website[END_REF]), the sampling rate was set to 1 Mega-samples/sec 1 . We varied the observation window from 1ms to 10 ms. The tested scenarios are presented in table I, and the results are shown in figure 3. The presented results were averaged over 1000 non-correlated experiments. As explained previously the aim of the algorithm is to triger a vertical handoff toward the access point where the traffic is lower. According to the figure we can clearly see that the channel occupancy rate is lower in the configurations where a lower bit-rate is required by users, and increases as the required bit-rate and number of users increase. We also observe that a 2 ms observation window length is sufficient to get

Fig. 1 .

 1 Fig. 1. (a) Absolute value of a wifi signal, (b) corresponding behavior of the criterion J (u)

  VIII. SIMULATIONS IEEE 802.11n signals are simulated. We recall that the IEEE 802.11n are 64 subcarriers OFDM signals with a cyclic prefix of length 16 [5]. The propagation channel {h(l)} l=0,••• ,L has an exponential decay profile for its non-null component (i.e., E[|h(l)| 2 ] = Ge -l/µ for l = 0, • • • , L) with L = 2/3.D, G is chosen such that P L l=0 E[|h k (l)| 2 ] = 1, and the RMS delay spread is set to 25 percent of the cyclic prefix. The channel is assumed to be time variant with a Doppler frequency equal to 10 Hz .

Fig. 2 .

 2 Fig. 2. NMSE of the channel occupency rate versus SNR


Yizhi Yin

Ramesh Pyndiah

Karine Amis

Performance of Random Linear Network Codes Concatenated with Reed-Solomon Codes Using Turbo Decoding

Keywords: Broadcast, network coding, random linear code, ordered statistics, block turbo code, error correction

published or not. The documents may come

Performance of random linear network codes concatenated with reed-solomon codes using turbo decoding

I. INTRODUCTION

In broadcast applications, the network coding [START_REF] Ahlswede | Network information flow[END_REF] is studied for its advantages over traditional store-and-forward strategy. Linear network coding [START_REF] Ahlswede | Network information flow[END_REF] [START_REF] Koetter | An algebraic approach to network coding[END_REF] resides in that the intermediate node can linearly combine its received packets (on a finite field) into one packet of the same length and send it to the following node. Each packet generated by the intermediate nodes will record in its header the information about the combination. By retrieving header information at the destination, we can establish the reverse operation (network decoding) to regenerate the original packets. This new strategy benefits from the correlation of different data flux heading to different sinks during the transmission. Different data flux can be fused and refined in an optimized way so that it will ease the passage of data through network bottlenecks. So the throughput of the network is improved and both energy and bandwidth are economized. According to [START_REF] Koetter | An algebraic approach to network coding[END_REF], it is sufficient to apply linear combinations over a Galois field (GF) to achieve the capacity of the error-free broadcast networks. Later, the random linear network (RLN) coding was proposed in [START_REF] Ho | The benefits of coding over routing in a randomized setting[END_REF] where the linear combination coefficients are chosen randomly, uniformly and independently over a GF with large enough cardinal. This This work has been financially supported by the Région Bretagne and Pracom.

randomness approach guarantees us a capacity-achieving configuration of network coding with a high probability that its reverse operation exists. More detailed discussion on network coding theory can be referred to [START_REF] Ho | Network Coding, An Introduction[END_REF].

In broadcast networks, the packets sent to a particular destination will travel through different routes and different delays will be introduced. Suppose that RLN coding is applied, then every packet received by the destination is a random linear combination of original packets. The destination needs to wait for at least K linearly independent packets before the network decoding. In this paper, we investigate a compromise between transmission delay and error rate performance. Fig. 1 illustrates the network where the source node S sends K packets to Access Points (AP) by broadcasting in an errorfree fixed network. The RLN coding is applied between S and APs to exploit the network redundancy. The end user T is attached to AP via a wireless channel, in which, the FEC code is used for improving transmission reliability.

The GF cardinal is supposed to be sufficiently large, so AP can accomplish network decoding with high probability as long as it has received at least K independent packets. A straightforward solution is to let AP do the network decoding. After getting the original packets, the AP encodes them with a FEC code and sends them to T. The system error rate performance is then determined by the FEC code. Its drawback is that the packet delay in the broadcast network can differ a lot [START_REF] Walsh | Optimal rate delay tradeoffs for multipath routed and network coded networks[END_REF] and each AP has to wait for K independent packets before network decoding and wireless transmission. Another problem happens when T moves from AP1 to AP2 during the wireless transmission. It will only receive M packets (M < K) from AP1. In order to recover all K original packets, T has to resort to AP2 either by indicating the (K-M ) missing packet numbers through signaling channel or with some central control unit over the APs to let AP2 know the missing packets. This solution suffers from additional cost.

In this paper, we consider another solution with high flexibility when the end user switches AP during the transmission. The AP starts the wireless transmission of the network coded packets (after FEC encoding) as soon as they are received at AP. If T changes its connected AP, it simply continues to receive further network coded packets from a new AP without additional signaling nor central control. It is at T where we decode both the channel and network codes. This solution is not tightly restricted by the packet delay in broadcast network so the transmission efficiency will be improved.

The remainder of the paper is organized as follows. In Section II, we describe the system setup and the coding scheme. In Section III, we study the error correction performance of RLN code and we introduce the ordered statistics decoding (OSD) [START_REF] Fossorier | Soft-decision decoding of linear block codes based on ordered statistics[END_REF]. In Section IV, we describe how to use the OSD in a turbo decoding algorithm. In Section V, the system frame error rate (FER) under AWGN and non frequency-selective block Rayleigh fading channels is analyzed. In Section VI, we conclude this paper.

II. SYSTEM DESCRIPTION

The system setup is depicted in Fig. 1. An error-free network connects the source S and the Access Points (AP) and RLN coding is applied in this fixed network. Without loss of generality, we suppose there are at least N disjoint paths from S to each AP. The end user T is connected to an AP via a wireless channel and can switch AP during the transmission. For each transmission, S maps data into K packets of length K ′ denoted by x 1 , x 2 ,. . ., x K , K < N on GF(q) and sends them out. The network nodes perform RLN coding with coefficients chosen randomly and uniformly over the data symbol alphabet GF(q). So, the output packet of an intermediate node is a random linear GF(q)-combination of its received packets. Each AP receives N packets of length K ′ . For a particular AP, we denote y 1 , y 2 , . . ., y N its received packets and we define matrices X and Y whose rows are respectively transmitted and received packets. The linear combinations from S to AP can be modeled by an N × K transfer matrix F as described in [START_REF] Koetter | An algebraic approach to network coding[END_REF]. The transfer equation is:

Y = F • X. (1)
According to [START_REF] Ho | The benefits of coding over routing in a randomized setting[END_REF], if the cardinal basis of GF(q) tends to infinity, the probability that the RLN coding yields a solution tends to 1. Then the network decoding is:

X = (F T F) -1 F T Y, (2)
where F T denotes the transpose of F and F -1 its inverse.

After detecting the original data, the AP transmits them to T using traditional wireless protocols (e.g. Wi-Fi). This straightforward solution mentioned in Section I has the drawback that the transmission delay is determined by the longest route in the broadcast network. Since transmission delay is a major system concern, we consider an alternative solution.

Suppose T is initially connected to AP1. After AP1 has received several network coded packets, instead of waiting for the others, it encodes them using a FEC code. We consider a linear FEC code defined on GF(q) (binary dimension K ′ m, binary length N ′ m) where m = log 2 (q). The resulting codewords are sent to T immediately. When T moves to AP2, it just continues to receive FEC coded packets from AP2. We assume a BPSK transmission and z i is the BPSK modulated FEC codeword. Supposing an AWGN channel, let b i denote a N ′ m-length vector with independent identically distributed (i.i.d) zero mean Gaussian components, and let r i be the i-th packet of length N ′ m received by the end user T:

r i = z i + b i . ( 3
)
By stacking successively all the received packets r i (i = 1 . . . N ), we build the received matrix R of dimension N × N ′ m. In the same way, we build the transmitted GF(q) matrix Z of dimension N × N ′ whose rows are the FEC codewords z i (i = 1 . . . N ) converted to GF(q). The row vectors of matrix Z are the FEC codewords. We denote the i-th column vector of X by x c i and its i-th row vector by x r i . We consider a systematic FEC code whose redundancy comes first so that the last K ′ columns of Z are the result of network coding in (1) that remains unchanged during FEC encoding. For simplicity, we continue to use Y to represent these columns of Z. From (1):

(y c i ) T = (x c i ) T • F T . ( 4
)
So the column code of Y can be modeled as a FEC code whose generator matrix is F T . We name this code as the RLN code. Providing the RLN code can be made systematic, a twodimensional product code can be built by serial concatenation of the systematic RLN code and the systematic FEC code.

To transform F T into its systematic form, we only need to apply the row elementary operation and column permutation on F T . It leads to a permutation of the column codeword elements (a row permutation on Z). So at the end user side, after arranging the rows of R according to the permutation applied to Z, we will get a product codeword, whose row code is the systematic FEC code used at APs and its column code is a random linear block code whose generator matrix is the systematic form of F T . This permutation is invertible after decoding. In the following, for the sake of simplicity, we keep the notations R, Z and F T to denote the resulting matrices after the systematic transformation.

III. ORDERED STATISTICS DECODING (OSD) AND

RANDOM LINEAR NETWORK CODE Some researchers have studied the error correction capability of network coding in erroneous networks. Refs. [START_REF] Yeung | Network error correction, part I: Basic concepts and upper bound[END_REF]- [START_REF] Koetter | Coding for errors and erasures in random network coding[END_REF] introduce different theories about the error-control in Network

Coding. These emerging theories analyze the problem within the vector space domain and they handle the network code through hard decoding approach and the reliability information is not exploited.

In this paper, we investigate the soft decoding of RLN code. For soft decoding, we use the log-likelihood ratio (LLR) as a measure of the reliability. Under AWGN channel and BPSK modulation assumption, the decoder input is:

LLR ij in = ln Pr{z ij = +1/r ij } Pr{z ij = -1/r ij } . ( 5
)
As described in Section II, the network coding transfer matrix F is updated for each transmission of a data matrix X. The corresponding RLN code doesn't have a fixed minimum Hamming distance. Decoding this kind of code with ordinary decoding methods is not efficient due to its variable generator matrix. So we propose to use the soft ordered statistics decoding (OSD) algorithm. It was proposed in [START_REF] Fossorier | Soft-decision decoding of linear block codes based on ordered statistics[END_REF] to decode a binary linear block code depending on the reliability measure.

Here is a brief review.

The linear codeword is symbol-to-binary converted (systematic or not) with binary code length and code dimension equal to N b and K b , respectively. We reorder the bits of received codeword c 1 according to their reliabilities in decreasing order and get c 2 . The reordering defines a permutation function λ 1 to be applied on the columns of generator matrix G, resulting in matrix G . We find the first K b independent columns with highest associated reliabilities in G ′ and put them as the first K b columns of a matrix G" of dimension K b × N b . The remaining (N b -K b ) columns of G ′ complete G". This selection defines a permutation function λ 2 to be used on c 2 getting a sequence c 3 . The first K b components of c 3 are called its most reliable independent (MRI) bits. With the MRI bits as message part and the corresponding systematic generator matrix G", we get a hard decision codeword â. Its corresponding original codeword can be obtained by applying permutation λ 1 -1 λ 2 -1 . The order-l reprocessing aims at making all possible changes of i bits (1 ≤ i ≤ l) among the K b MRI bits and re-encoding them into candidate codewords with the generator matrix G". Finally, we choose the most probable one as the decoding result.

The OSD is compatible with the variable generator matrix of RLN code without much additional complexity. Here is an example of order-1 OSD to decode a RLN code of dimension 5 and code length form 7, 8 to 9 on GF [START_REF] Yeung | Network error correction, part II: Lower bounds[END_REF]. In Table I, we build a statistics by generating a large number of such random codes and classify them according to their minimum Hamming distances (d min ). We shall now calculate the union bound of the FER as:

FER| dmin=k 1 2 n ∑ d=k Wd | dmin=k • erfc ( √ dR a E b N 0 ) . (6) FER = max(dmin) ∑ dmin=1 (FER| dmin=k ) • Pr (d min = k) , ( 7
)
where FER| dmin=k is the average FER in condition that all the codes tested satisfy d min = k and Wd is the average multiplicity of codewords with weight d. erfc(x) is the complementary error function and R a is the code rate. The bounds and simulated FER of the RLN codes are illustrated in Fig. 2. We observe that the order-1 OSD successfully decodes such random codes as the simulated FER reaches its union bounds (in the three cases, both curves are superposed for FER less than 4 × 10 -3 ). Moreover, simulated results match with Table I : with each unity increase of code length, the percent of d min = 1 has decreased by a factor of about 10, resulting in an improvement of the error rates. d min = 1 is to be avoided as it means no error correction capability.

IV. ITERATIVE DECODING AT DESTINATION

Since we can piece together the received packets into a product code, we will apply turbo decoding at the end user. The turbo decoder receives a noisy BPSK modulated N ×N ′ m matrix R corresponding to the product codeword whose row code is a linear block FEC code and its column code is a RLN code. The Soft-in-Soft-out (SISO) decoding is performed iteratively on R using the Chase-Pyndiah algorithm [START_REF] Pyndiah | Near optimum decoding of product codes[END_REF] reinforced by the OSD for the RLN code. For the i-th row of R: r r i = (r i1 , r i2 , . . . , r in ) where n = N ′ m, the Chase algorithm [START_REF] Chase | A class of algorithms for decoding block codes with channel measurement information[END_REF] generates a subset Ω 1 of the most probable row codewords. Within Ω 1 , the most probable codeword is noted

d = (d 1 , d 2 , . . . , d n ), d k ∈ {±1}. We define the competing codeword c (j) = (c 1 , c 2 , . . . , c n ), c k ∈ {±1}, for a specific position j as the most probable codeword in Ω 1 with c (j) j = d j .
Then the soft output LLR of the decision d j is calculated as:

LLR ij out = ⟨LLR i in , d⟩ -⟨LLR i in , c (j) ⟩ 2 • d j , ( 8
)
where ⟨a, b⟩ = ∑ k a k b k is the correlation between two vectors a and b. The extrinsic information is obtained by:

W ij = LLR ij out -LLR ij in . ( 9
)
If there is no competing codeword c (j) in Ω 1 , then we use:

W ij = β × d j , ( 10
)
where β is a constant to be optimized through simulations.

For the columns of R, we use the SISO OSD [START_REF] Fossorier | Soft-input soft-output decoding of linear block codes based on ordered statistics[END_REF]. The order-1 reprocessing is used to generate the candidate column codeword subset Ω 2 to find the competing codeword. The rest of the column decoding procedure is the same as the row one.

We update the input LLR every half iteration by:

LLR in (m + 1) = LLR in (0) + α(m + 1) • W(m), (11)
where the LLR in (0) is the channel output LLR and W(m) is the extrinsic information obtained at the previous iteration. At the first iteration, W(0) is set to zero. α(m) is a scaling factor increasing with the iteration step m.

V. SIMULATION RESULTS

For each transmission of a data matrix X, the matrix F is constructed randomly. T will test if F has full rank. An additional packet will be requested until having K independent ones. The more packets are received by T, the more probable F is invertible. Another interest of increasing the packets received by T is to improve the d min distribution of the RLN code as shown in Table I and Fig. 2.

We investigate the performance of our solution on AWGN and non frequency-selective block Rayleigh fading channels using BPSK modulation. At destination, we suppose perfect synchronization, no intersymbol interference (ISI), perfect channel information and knowledge of the transfer matrix F. Our simulation is based on GF(32). The FEC code is the Reed-Solomon (RS) code(31, 29). Different sizes of generator matrix F T with K = 5, N = 7, 8 or 9 are tested.

The Chase algorithm considers the 4 least reliable bits whose indices compose the subset Θ. For an 8-iteration turbo decoding, we use the following scaling factor α and reliability factor β for the different half-iterations: α = {0, 0.1, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0.55, 0.6, 0.65, 0.7, 0.9, 1.0, 1.0}.

β = min   ∑ j∈Θ |r in j |, 5   .
|r in j | is the reliability of the j-th least reliable input binary component. For the order-1 reprocessed OSD, we choose the competing codeword among the 6 most probable codewords in Ω 2 to eliminate possible pseudo-competing codewords.

Here, we consider the FER without the error propagation due to network decoding that will influence the BER. Fig. 3 plots the FER in AWGN channel. The curve "RS(31, 29) soft" is the reference representing the soft decoding of RS(31, 29) using Chase algorithm with 16 test patterns. It can be viewed as the straightforward solution where network decoding is applied at AP before transmitting all the RS encoded original packets to T. As pointed out in Section II, this solution will suffer from long delays. For the proposed (7, 5) network code based scheme, the SNR gain is equal to 2dB for a FER of 10 -2 and increases with the number of packets received by T: 2.5dB (resp. 2.75dB) for N = 8 (resp. N = 9). In Fig. 3, we observe an error floor phenomenon for our proposed solution especially for high SNR. We analyze it using Fig. 4 where we consider only the case of RLN code [START_REF] Yeung | Network error correction, part I: Basic concepts and upper bound[END_REF][START_REF] Walsh | Optimal rate delay tradeoffs for multipath routed and network coded networks[END_REF]. The curve "Turbo (7, 5) d min = 1" (resp. "Turbo (7, 5) d min > 1") is the FER of our proposed solution when the RLN code has d min = 1 (resp. d min > 1). The overall average FER curve for RLN code (7, 5) with all possible minimum distances is the one with diamond marks.

Based on statistics tests, we estimate p = Pr(d min = 1) ≈ 6.8 × 10 -3 for the RLN code (7, 5) on GF(32), and use the approximation of the FER given in [START_REF] Pyndiah | Near optimum decoding of product codes[END_REF]:

FER = max(dmin) ∑ dmin=1 (FER| dmin=k ) • Pr (d min = k) = FER| dmin=1 • p + FER| dmin>1 • (1 -p). ( 12
)
The condition d min refers to the minimum Hamming distance of the random component code in the product code. Since erfc( √ x) decreases exponentially with x, from ( 6) and ( 12), at low SNR where FER| dmin=1 ≈ 1 and FER| dmin>1 ≫ p:

FER ≈ 1 • p + FER| dmin>1 ≈ FER| dmin>1 . ( 13
)
It corresponds to the superposition of the overall FER curve and the curve "Turbo (7, 5) d min > 1" below 3.5dB in Fig. 4.

At high SNR, FER| dmin=1 cannot be approximated by 1 and FER| dmin>1 becomes negligible to p • FER| dmin=1 :

FER ≈ p•FER| dmin=1 +FER| dmin>1 ≈ p•FER| dmin=1 . (14)
The factor p in (14) corresponds to the vertical gap at high SNR in Fig. 4 between the overall average FER curve and the "Turbo (7, 5) d min = 1" curve. This factor causes the error floor effect. If it was possible to eliminate the configuration d min = 1, the error floor effect would disappear as shown by the curve with pentagrams. But since the network coding is randomly generated by the network nodes, there is no evident solution. An intermediate solution is to check if the network configuration leads to a RLN code with d min = 1 and request additional packets if necessary. We observe that the percent of RLN code [START_REF] Yeung | Network error correction, part I: Basic concepts and upper bound[END_REF][START_REF] Walsh | Optimal rate delay tradeoffs for multipath routed and network coded networks[END_REF] having d min = 1 decreases with the Galois field cardinal q (10% for q = 8 against 0.7% for q = 32 and 0.17% for q = 64). It may be interesting to deeply investigate the influence of q over d min .

Fig. 4. Proposed scheme performance on AWGN channel, with a (7,5) network code, BPSK and GF(32) Fig. 5 plots the FER for a flat Rayleigh block fading channel. Each row of the product code is attenuated with an independent Rayleigh coefficient. At 10 -2 , the coding gain compared to the AP network decoding solution is about 8dB with a RLN code (7, 5). There is an average improvement of 2.5dB per additional received packet.

VI. CONCLUSION

In this paper, we propose a flexible transmission scheme for the broadcast network where network coding is applied. As the random linear network code can be modeled as a random linear block FEC code, it can be softly decoded by the ordered statistics decoding algorithm. Simulations proved that, the more coded packets collected from the network, the better the error correction performance. Increasing the received packet number ameliorates the minimum Hamming distance distribution of the resulting random code. The error floor effect is analyzed using the union bounds and one possible solution is to avoid the network coding configurations leading to a minimum Hamming distance of 1.

Fig. 1 .

 1 Fig. 1. One-source multicast network structure using network coding

Fig. 2 .

 2 Fig. 2. Union bounds / FER simulated of OSD(1) on GF(8) in AWGN channel with BPSK modulation

Fig. 3 .

 3 Fig. 3. Proposed scheme performance on AWGN channel, BPSK / GF(32)

Fig. 5 .

 5 Fig. 5. Proposed scheme performance on Rayleigh block fading channel with BPSK and GF(32)

ACKNOWLEDGMENT

The authors thank Raphaël Le Bidan for the fruitful discussions and his interesting suggestions.


Pedro Casas
email: pcasas@fing.edu.uy

Federico Larroca
email: federico.larroca@telecom-paristech.fr

Sandrine Vaton

Robust Routing Mechanisms for Intradomain Traffic Engineering in Dynamic Networks

Keywords: Traffic Uncertainty, Proactive Traffic Management, Robust Optimization, Stable and Reactive Robust Routing, End-to-End Routing Performance Evaluation

Internet traffic is highly dynamic and difficult to predict in current network scenarios. This makes of Traffic Engineering (TE) a very challenging task for network management and resources optimization. We study the problem of Intradomain Routing Optimization under this traffic uncertainty. Recent works have proposed robust optimization techniques to tackle the problem, conceiving the Robust Routing (RR) approach. RR copes with traffic uncertainty in an off-line preemptive fashion, computing a single static routing configuration that is optimized for traffic variations within some predefined uncertainty set. Despite achieving routing reliability with relatively low performance loss, RR presents various drawbacks and conception problems as it is currently proposed. This paper brings insight into the different Robust Routing shortcomings, introducing new mechanisms that improve previous proposals and alleviate these problems. Among others, we propose and evaluate new optimization objectives to attain better global performance from an end-to-end quality of service perspective.

I. INTRODUCTION

As network services and Internet applications evolve, network traffic is becoming increasingly complex and dynamic. The convergence of data, telephony and television services on an all-IP network directly translates into a much higher variability and complexity of the traffic injected into the network. Recent Internet traffic studies from major network technology vendors like Cisco Systems forecast the advent of the Exabyte era [START_REF]Global IP Traffic Forecast and Methodology 2006-2011[END_REF], a massive increase in network traffic driven by high-definition video. Furthermore, current evolution and deployment-rate of broadband access technologies (e.g. Fiber To The Home technology) is such that the old assumption of infinitely provisioned core links will soon become obsolete; market research reports like [START_REF]Worldwide Markets for Fiberoptics in Broadband Access Networks: Fiberoptic Demand Forecast in Broadband Access Networks[END_REF] forecast a value of bandwidth demand per user as high as 50 Gb/sec in 2030. In this context, simply upgrading link capacities may no longer be an economically viable solution to cope with dynamic traffic. To make matters worse, the presence of unexpected events such as network equipment failures, largevolume network attacks, flash crowd occurrences and even external routing modifications induces large uncertainty in traffic patterns. In the light of this traffic scenario, we study the problem of intradomain routing optimization under traffic uncertainty. This uncertainty is assumed to be an exogenous traffic modification, meaning that traffic variations are not produced within the domain for which routing is optimized but are due to external and difficult to predict events.

Recent works [START_REF] Ben-Ameur | Routing of Uncertain Traffic Demands[END_REF]- [START_REF] Wang | COPE: Traffic Engineering in Dynamic Networks[END_REF] have proposed a plausible solution to the routing optimization under traffic uncertainty problem: the Robust Routing (RR) approach. In RR, traffic uncertainty is taken into account directly within the routing optimization, computing a single routing configuration for all traffic demands within some uncertainty set where traffic is assumed to vary. This uncertainty set can be defined in different ways, depending on the available information: largest values of links load previously seen, a set of previously observed traffic demands (previous day, same day of the previous week), etc. The criteria to search for this unique routing configuration in these works is to minimize the maximum link utilization (i.e., the utilization of the most loaded link in the network) for all traffic demands of the corresponding uncertainty set. While this routing configuration is not optimal for any single traffic demand within the set, it minimizes the worst case performance over the whole set.

The RR approach can be used as a preemptive TE technique to deal with dynamic and uncertain traffic demands. It can handle unexpected traffic variations with relatively low performance loss, depending on the size of the uncertainty set. However, RR presents some conception problems and serious shortcomings in its current state which we highlight and try to ease in this work. The first drawback of current RR is related to the objective function it intends to minimize. Optimization under uncertainty is generally more complex than classical optimization, which forces the use of simpler optimization criteria such as maximum link utilization (MLU). The MLU is not the most suitable network-wide optimization criterion; setting the focus too strictly on MLU often leads to worse distributions of traffic, adversely affecting the mean network load and thus the total network end-to-end delay, an important QoS indicator. It is easy to see that the minimization of the MLU in a network topology with heterogeneous link capacities may lead to poor results as regards global network performance. The second drawback of RR we identify is its inherent dependence on the definition of the uncertainty set of traffic demands: larger sets allow to handle a broader group of traffic demands, but at the cost of routing inefficiency; conversely, tighter sets produce more efficient routing schemes, but subject to poor performance guarantees. Thus, considering a unique RR configuration to address both traffic in normal operation and unexpected traffic variations is an inefficient strategy: a single routing can not be suitable for both situations.

A. Related Work

There is a large literature on routing optimization with uncertain traffic demands. Traditional algorithms rely on a single or a small group of expected traffic demands to compute optimal and reliable routing configurations. An extreme case is presented in [START_REF] Roughan | Traffic Engineering with Estimated Traffic Matrices[END_REF], where routing is optimized for a single estimated traffic demand and is then applied for daily routing. Traffic uncertainty is characterized by multiple traffic demands in [START_REF] Zhang | On Optimal Routing with Multiple Traffic Matrices[END_REF] (set of traffic demands from previous day, same day of previous week, etc.), where different mechanisms to find optimal routes for the set are presented. In the previously described scenario this perspective is no longer suitable (see [START_REF] Casas | An Adaptive Multi Temporal Approach for Robust Routing[END_REF] for further arguments in this sense).

A different approach has emerged in the recent years to cope with the traffic increasing dynamism and the need for costeffective solutions, Dynamic Load-Balancing (DLB) [START_REF] Elwalid | MATE: MPLS Adaptive Traffic Engineering[END_REF]- [START_REF] Larroca | Minimum-Delay Load-Balancing Through Non-Parametric Regression[END_REF]. In DLB, traffic is split among a priori established paths in order to avoid network congestion. The two most well-known proposals in this area are MATE and TeXCP. In MATE [START_REF] Elwalid | MATE: MPLS Adaptive Traffic Engineering[END_REF], a convex link cost function is defined, which depends on the link capacity and the link load. The objective is to minimize the total network cost, for which a simple gradient descent method is proposed. TeXCP [START_REF] Kandula | Walking the Tightrope: Responsive yet Stable Traffic Engineering[END_REF] proposes a somewhat simpler objective: minimize the biggest utilization each traffic demand obtains in its paths. Another DLB scheme which has the same objective but a relatively different mechanism is REPLEX [START_REF] Fischer | REPLEX: dynamic traffic engineering based on wardrop routing policies[END_REF]. In [START_REF] Larroca | Minimum-Delay Load-Balancing Through Non-Parametric Regression[END_REF], we use a link cost function based on measurements of the queueing delay, which results in better global performance from a QoS perspective. DLB presents a desirable property, that of keeping routing adapted to dynamic traffic. However, DLB algorithms present a trade-off between adaptability and stability which might be particularly difficulty to address under significant and abrupt traffic changes. Besides, network operators are reluctant to use dynamic mechanisms and prefer stable routing configurations, as they claim they get a better feeling of what is going on in the network.

The last category of algorithms consists of Robust Routing techniques [START_REF] Ben-Ameur | Routing of Uncertain Traffic Demands[END_REF]- [START_REF] Juva | Robust Load Balancing[END_REF]. The objective in RR is to find a unique static routing configuration that fulfills a certain criterion for a broad set of traffic demands, generally the one that minimizes the maximum link utilization over the whole set of demands. In [START_REF] Ben-Ameur | Routing of Uncertain Traffic Demands[END_REF], authors capture traffic variations by introducing a polyhedral set of demands, which allows for easier and faster linear optimization. [START_REF] Johansson | Data-driven Traffic Engineering: techniques, experiences and challenges[END_REF] applies this robust technique to compute a robust MPLS routing configuration without depending on traffic demand estimation, and discusses corresponding methods for robust OSPF optimization. Oblivious Routing [START_REF] Applegate | Making Intra-Domain Routing Robust to Changing and Uncertain Traffic Demands: Understanding Fundamental Tradeoffs[END_REF] also defines linear algorithms to optimize worst-case MLU for different sizes of traffic uncertainty sets. [START_REF] Juva | Robust Load Balancing[END_REF] analyses the use of robust routing through a combination of traffic estimation techniques and its corresponding estimation error bounds, in order to shrink the set of traffic demands. In [START_REF] Wang | COPE: Traffic Engineering in Dynamic Networks[END_REF] authors introduce COPE, a RR mechanisms that optimizes routing for predicted demands and bounds worst-case MLU to ensure acceptable efficiency under unexpected traffic events. As we mentioned before, RR presents two important shortcomings, the former related to the objective function it intends to minimize (i.e., MLU minimization) and the latter as an inherent consequence of its stability property (i.e., using a single routing configuration for all traffic events).

B. Contributions of the Paper

In this paper we propose and evaluate new variants of the RR approach to alleviate the two problems identified in current proposals. As regards the objective function to minimize, we firstly propose to minimize the mean link utilization instead of the MLU. The mean link utilization provides a better image of network-wide performance, as it does not depend on the particular load or capacity of each single link in the network but on the average value. A direct minimization of the mean link utilization does not assure a bounded MLU, which is not practical from an operational point of view. Thus, we minimize the mean link utilization while bounding the MLU by a certain utilization threshold a priori defined. This adds a new difficult to set constraint to the problem, namely how to define this utilization threshold. We further improve our proposal by providing a multiple objective optimization criterion, where both the MLU and the mean link utilization are minimized simultaneously. We evaluate the improvements of our proposals from a QoS perspective, using the mean path end-to-end queuing delay as a measure of global performance.

Regarding the trade-off between routing performance and routing reliability, the only previous work that has evidenced the problem and proposed some solution is COPE [START_REF] Wang | COPE: Traffic Engineering in Dynamic Networks[END_REF]. Nevertheless, COPE proposes a single routing configuration to handle expected as well as large and abrupt traffic variations, which is clearly not the best solution. In [START_REF] Casas | Robust and Reactive Traffic Engineering for Dynamic Traffic Demands[END_REF] we have recently proposed a solution to manage this trade-off, known as Reactive Robust Routing (RRR). Basically, RRR consists of constructing a RR configuration for expected traffic in nominal operation, adapting this nominal routing configuration after the detection and localization of a large and long-lived traffic modification. RRR provides good performance for both nominal operation and unexpected traffic, but it is difficult to deploy in a real implementation, because of the routing reconfiguration step. Reconfiguring the routing of an entire Autonomous System is a nontrivial task. In this paper we modify the RRR approach, using a preemptive Load Balancing algorithm to balance traffic among prestablished paths after the localization of a large volume traffic modification.

The remainder of this paper is organized as follows. In section II we recall the traditional Robust Routing approach, exposing the above mentioned problems with some real evaluations. Section III presents the proposed variants to the former RR mechanism, designed to improve the detected shortcomings. The evaluation of the proposed algorithms under different traffic scenarios and the discussion of results are provided in section IV. Finally, section V concludes this work and presents some future perspectives.

II. STABLE ROBUST ROUTING

Let us begin by introducing the notation used in this paper. The network topology is defined by n nodes and a set L = {l 1 , . . . , l q } of q links, each with a corresponding capacity c i , i = 1, . . . , q. The Traffic Matrix (TM) X = {x i,j } denotes the traffic demand between every origin node i and every destination node j (i = j) of the network; we shall note each of these origin-destination pairs as OD pairs, and each origindestination traffic demand x i,j as OD flows. Let X = {x k } be the vector representation of the TM, where we have reordered OD flows by index k = 1, . . . , m (m = n.(n -1)). Let N = {OD 1 , . . . , OD m } be the set of m OD pairs. We consider a multi-path network topology, where each OD flow x k can be arbitrarily split among a set of p k origin-destinations paths P k . In this sense, we shall call r k p the portion of traffic flow x k sent through path p ∈ P k , where 0 r k p 1 and p∈P k r k p = 1. Let λ p l be an indicator variable that takes value 1 if path p traverses link l and 0 otherwise, and Y = {ρ 1 , . . . , ρ q } a vector representation of links traffic load. Then X and Y are related through the routing matrix R, a q ×m matrix R = {r k l } where r k l = p∈P k λ p l . r k p . The variable r k l indicates the fraction of OD flow x k routed through link l; this results in the following relation:

Y = R . X (1)
Given X, the multi-path routing optimization problem consists in choosing the set of paths P k for each OD pair k and computing the routing matrix R, in order to optimize a certain objective function f (X, R). A simplified version of this problem is the optimal load-balancing problem which, given a set of paths, calculates R. The most popular TE objective function f (X, R) has traditionally been the maximum link utilization u max , defined as:

u max (X, R) = max l∈L {u l } (2)
where u l = ρ l /c l stands for the link utilization; a value of u l close to 1 indicates that the link is operating near its capacity. Network operators usually prefer to keep links utilization relatively low in order to support sudden traffic increases and link/node failures.

Finding a multi-path routing configuration minimizing u max is an instance of the classical multi-commodity flow problem which can be formulated as a simple linear program [START_REF] Ahuja | Network Flows: Theory, Algorithms, and Applications[END_REF]. For a single known traffic matrix X, the problem can be easily solved by linear programming techniques [START_REF] Mitra | A Case Study of Multiservice, Multipriority Traffic Engineering Design for Data Netorks[END_REF]. However, as we have previously discussed, traffic demands are uncertain and difficult to predict in current scenario, and all we can expect is to find them within some bounded uncertainty set.

In a robust perspective of the multi-path routing optimization problem, demand uncertainty is taken into account within the routing optimization, computing a single routing configuration for all demands within some uncertainty set. In this work we consider a polyhedral uncertainty set X, more precisely a polytope as in [START_REF] Ben-Ameur | Routing of Uncertain Traffic Demands[END_REF], based on the intersection of several half-spaces that result from linear constraints imposed to traffic demand. minimize umax [START_REF] Ahuja | Network Flows: Theory, Algorithms, and Applications[END_REF] subject to:

k∈N p∈P k λ p l . r k p . x(k) umax.c l ∀ l ∈ L, ∀ X ∈ X p∈P (k) r k p = 1 ∀ k ∈ N r k p 0 ∀ p ∈ P k , ∀ k ∈ N umax 1
As an example, let us define an uncertainty set X based on a given routing matrix R o and the peak-hour links traffic load Y peak obtained with this routing matrix:

X = X ∈ R m , R o .X Y peak , X 0
Observe that this definition of the uncertainty set has a major advantage: routing optimization can be performed from easily available links traffic load Y without even knowing the actual value of the traffic demand X. Figure 1 depicts the obtained uncertainty set, based on the convex intersection of q halfspaces of the form

r i o • X ρ peak i
, ∀i ∈ L, where r i o stands for the i-th row of the routing matrix R o . The traditional Robust Routing Optimization Problem (RROP) defined in (3) consists of minimizing the maximum link utilization u max , considering all demands within X. The solution to the problem is twofold: on the one hand, a routing configuration R robust , and on the other hand, a worst-case performance threshold u robust max :

r 1 o r 2 o r 3 o r i o r q-1 o r q o X X
R robust = argmin R max X∈X u max (X, R) u robust max = max X∈X u max (X, R)
Given a suitable definition of the uncertainty set, the obtained robust routing configuration R robust is applied during long periods of time; in this sense, we refer to Robust Routing as Stable Robust Routing (SRR). The authors of [START_REF] Ben-Ameur | Routing of Uncertain Traffic Demands[END_REF] have shown that the RROP can be efficiently solved by linear programming techniques, applying a combined columns and constraints generation method. This method iteratively solves the problem, progressively adding new constraints and new columns to the problem. The new constraints are the extreme points of the uncertainty set X, and the new columns represent new paths added to reduce the objective function value. Only extreme points of X are added as new constraints, as it is easy to see that every traffic demand X ∈ X can be expressed as a linear combination of these extreme demands. Regarding new added paths, the algorithm in [START_REF] Ben-Ameur | Routing of Uncertain Traffic Demands[END_REF] may not be the best choice from a practical point of view since the number of paths for each OD pair is not a priori restricted and the characteristics of added paths are not controlled. For example, it would be interesting to have disjoint paths to route traffic from each single OD pair, improving resilience. For this reason we modify the algorithm to select new paths, both limiting the maximum number of paths in P k and taking as new candidates the shortest paths with respect to link weights w i l :

w i l = 1 + 1 -r k l i (4)
where r k l i corresponds to the fraction of traffic flow x k that traverses link l after iteration i and is a small constant that avoids numerical problems. If OD pair k uses a single path p at iteration i, r k l i = 1 for every link l ∈ p, and so this path is removed from the graph where new shortest paths are computed (w l → ∞, ∀l ∈ p). While this may result in a sub-optimal performance, it allows a real and practical implementation. In case there are no disjoint paths for OD pair k, we use the column constraint generation method used in [START_REF] Ben-Ameur | Routing of Uncertain Traffic Demands[END_REF] to add new paths for OD pair k.

A. Shortcomings of Stable Robust Routing: An Example

In this section we shall present some simulations that will help us gain insight into the Stable Robust Routing mechanism and highlight the previously discussed shortcomings. We use the Abilene Network as the environment for simulations. Abilene is a high-speed Internet2 backbone network, connecting 12 router-level nodes through 30 optical links (we only consider intra-domain links). The used router-level network topology and traffic demands are available at [START_REF] Zhang | Abilene Dataset 04[END_REF]. Traffic data consists of 6-month traffic matrices collected every 5' via Netflow from the Abilene Observatory [START_REF]The Abilene Observatory[END_REF]. As measured traffic demands do not significantly load the network, we rescaled them by multiplying all their entries by a constant. The dataset in [START_REF] Zhang | Abilene Dataset 04[END_REF] also provides the static routing configuration R o deployed in Abilene during 6-month TMs measurement campaign.

Let us discuss the issue related to the objective function used in the traditional SRR algorithm. As we stated in the Introduction, the maximum link utilization is a local performance indicator, and a routing configuration minimizing u max may often lead to a worse distribution of traffic, adversely affecting the global performance of the network. Besides, while it is true that overloaded links tend to cause QoS degradation (e.g., larger delays and packet losses, throughput reduction, etc.), u max does not represent a direct QoS indicator, a desirable property in the context of QoS provisioning.

In order to evaluate SRR from a network-wide QoS perspective, let us consider a performance indicator directly related to QoS: the path end-to-end (e2e) delay. The e2e delay on a path is the sum of the delays on each link of the path. The delay on each link consists of two components, namely the queuing

f M/M/1 l (ρ l ) ρ l
Fig. 2. Mean queue size, measurements and approximations delay (i.e., buffer and service delay) and the link propagation delay. The former depends on the link load, while the latter is constant. In this sense and as a simplification to the problem, we shall consider the e2e path queuing delay as a measure of performance. Assume that queuing delay on link l is given by the function d l (ρ l ). Given this function, we can compute the e2e queuing delay of path p as d p = l∈p d l (ρ l ). In order to evaluate the network-wide performance of SRR, we define the expected e2e path queuing delay d mean as follows:

d mean (X, R) = k∈N p∈P k r k p . x k d p = l∈L ρ l . d l (ρ l ) (5)
That is to say, a weighted mean e2e queuing delay, where the weight for each path is how much traffic is sent through it (r k p . x k ), or in terms of links, the weight for each link is how much traffic is traversing it (ρ l ). A large mean e2e queuing delay translates into bad performance for all the traffic and not only for the traffic that traverses a particular loaded link. We prefer a weighted mean queuing delay to a simple total delay because it reflects more precisely performance as perceived by traffic. Two situations where the total delay is the same, but in one of them most of the traffic is traversing heavily delayed links should not be considered as equivalent. Note that, by Little's law, the value f l (ρ l ) = ρ l . d l (ρ l ) is proportional to the volume of data in the queue of link l. We will then use this last value as the addend in the last sum in [START_REF] Roughan | Traffic Engineering with Estimated Traffic Matrices[END_REF], since it is easier to measure than the queuing delay.

The function f l (ρ l ) is unknown and in the literature it is generally estimated using a classical M/M/1 model, where [START_REF] Kleinrock | Queueing Systems[END_REF]. However, in [START_REF] Larroca | Minimum-Delay Load-Balancing Through Non-Parametric Regression[END_REF] we show that a simple M/M/1 model has little to do with reality, and so we propose to use a non-parametric regression technique to estimate f l (ρ l ) from measurements without assuming any given model. Figure 2 depicts the real mean queue size of an operational network link at Tokyo obtained from [START_REF] Cho | WIDE-TRANSIT 150 Megabit Ethernet Trace 2008-03-18[END_REF], together with the M/M/1 estimation f M/M/1 l (ρ l ) and the non-parametric regression fl (ρ l ). It is clear that f M/M/1 l (ρ l ) consistently underestimates the real queue size value, while fl (ρ l ) provides quite accurate results.

f M/M/1 l (ρ l ) = ρ l /(c l -ρ l )
Let us evaluate the performance of SRR as regards both u max and d mean . From now on we shall use RROP as a reference to SRR, recalling that the robust routing optimization problem is the one described in [START_REF] Ahuja | Network Flows: Theory, Algorithms, and Applications[END_REF]. In this evaluation we consider a traffic scenario that presents and abrupt and large volume increase due to an external routing modification. This corresponds to the TMs with indexes between 1050 and 1200 from dataset X23 in [START_REF] Zhang | Abilene Dataset 04[END_REF]. Figure 3 depicts the described traffic scenario.

The evaluation starts with a normal low traffic load situation, but after the 100th minute one of the OD flows abruptly increases its traffic volume, loading the links it traverses until the end of the evaluation. Based on the static routing matrix of Abilene R o we define two different polytopes, the former adapted to the Low Traffic Load period (LTL period, before the 100th minute) and the latter adapted to the High Traffic Load period (HTL period, after the 100th minute):

X LTL = {X ∈ R m , R o .X Y LTL , X 0} X HTL = {X ∈ R m , R o .X Y HTL , X 0}
We assume that traffic is known in advance in both definitions, and take Y LTL and Y HTL as the maximum link load values observed during the LTL and HTL periods respectively. We compute two different robust routing configurations for both polytopes; RROP LTL corresponds to the SRR configuration for polytope X LTL , and RROP HTL for polytope X HTL . In this evaluation, both RROP LTL and RROP HTL use the same set of paths, namely the paths obtained from (3) for polytope X LTL . Solving (3) for a given set of paths consists of only adding new extreme points of polytope X (i.e., only new constraints are added).

Figure 3 depicts (a) the maximum link utilization u max and (b) the mean end-to-end queuing delay d mean during the evaluation period. Both u max and d mean are updated every 5 minutes, when a new TM is measured. As a reference for comparison, we also compute the minimum values u opt max (X) and d opt mean (X) for every single TM X of the evaluation period, using once again the same set of paths. u opt max (X) is computed with a simplified version of (3), where there is only one TM X instead of a set X. d opt mean (X) is computed using the algorithms in [START_REF] Larroca | Minimum-Delay Load-Balancing Through Non-Parametric Regression[END_REF]. Figure 4 optimal performance as regards u max (a difference smaller than 4%), RROP HTL obtains a queuing delay that constantly exceeds the optimum by almost 40% under a fairly network load. Such a difference may not be even acceptable from a QoS perspective, where end-to-end delays are even more important than network congestion. As we will show latter, this loss in performance is a direct consequence of the local criterion used in RROP.

The second interesting observation comes from the difference between RROP LTH and RROP LTL performances before and after the abrupt traffic volume increase; figure 3(a) shows that, despite an almost negligible network load, RROP LTL outperforms RROP HTL by almost 50% of relative utilization during the LTL period, while the opposite happens during the HTL period. The difference is not that big as regards delay before the 100th minute, but it becomes really important after the volume increase, where RROP LTL obtains a really bad performance. These results are somehow expected given the polytopes definition, and evidences both the dependence of RROP on the uncertainty set definition and the inherent consequence of using a single static routing configuration under large traffic variations. A final remark about this simple evaluation and the definition of the uncertainty set; we have considered that traffic was known in advance for the definition of both polytopes X LTL and X HTL . While traffic during the LTL period is easy to predict, the definition of X HTL in a real traffic scenario is a challenging task. We will come back to this issue in the following section.

III. NEW MECHANISMS FOR ROBUST ROUTING

In this section we shall present three enhanced mechanisms to overcome the problems of SRR evidenced in the previous evaluation. We will first introduce and evaluate two similar mechanisms to attain better global performance as regards endto-end delay. Then we will present a mechanism to manage the problem of defining the uncertainty set under unexpected and large traffic variations, previously introduced in [START_REF] Casas | Robust and Reactive Traffic Engineering for Dynamic Traffic Demands[END_REF]. In this work we provide some slight modifications to the former algorithm which allow a real deployment of the proposal. minimize umean [START_REF] Casas | An Adaptive Multi Temporal Approach for Robust Routing[END_REF] subject to:

l∈L k∈N p∈P k 1 c l λ p l . r k p . x(k) umean.q ∀ X ∈ X k∈N p∈P k λ p l . r k p . x(k) u thres max .c l ∀ l ∈ L, ∀ X ∈ X p∈P (k) r k p = 1 ∀ k ∈ N r k p 0 ∀ p ∈ P k , ∀ k ∈ N

A. Improving Network-Wide Performance

As we showed in figure 3(b), the minimization of u max leads to a distribution of traffic that results in an excessive end-to-end delay. Using the mean delay d mean (X, R) as the objective function in (3) would be an interesting approach to ease the problem; however, f l (ρ l ) is a non-linear function and the optimization problem becomes too difficult to solve. As we previously said, optimization under uncertainty is more complex than classical optimization and simple optimization criteria should be used. Let us consider a very simple networkwide linear objective function, namely the mean link utilization u mean (X, R), defined as:

u mean (X, R) = 1 q l∈L u l (6)
The mean link utilization considers at the same time the load of every link in the network and not only the utilization of the most loaded link; as we will show in the results, such an objective function provides a better global performance as regards end-to-end delay. However, a direct minimization of u mean does not assure a bounded maximum link utilization, which is not practical from an operational point of view. In this sense, we propose to change the objective function in (3) by u mean , while bounding the maximum link utilization by a certain threshold u thres max a priori defined. This results in the Robust Routing Mean Utilization Optimization Problem (RRMP) defined in [START_REF] Casas | An Adaptive Multi Temporal Approach for Robust Routing[END_REF].

Problem [START_REF] Casas | An Adaptive Multi Temporal Approach for Robust Routing[END_REF] is solved in the same way as (3), using the same recursive algorithm proposed in [START_REF] Ben-Ameur | Routing of Uncertain Traffic Demands[END_REF]. Note that (7) adds only a new constraint per each new traffic demand in X (in fact, for each extreme point of X). The drawback of ( 7) is on its dependence on the value of u thres max , which directly influences the routing performance as we will shortly see. An interesting choice for u thres max would be to use the output of (3), namely u robust max . To some extent this would result in a similar routing solution but with better traffic balancing.

A alternative approach would be to minimize both the value of u max and u mean at the same time, what constitutes a problem of multi-objective optimization. The problem with multi-objective optimization is that traditional single-objective optimization techniques can not be directly applied. An intuitive and easy approach to solve a multi-objective optimization problem using single-objective optimization techniques is to construct a single aggregated objective function (AOF) that combines both objective functions. We define a weighted linear minimize uaof = β . umax + (1 -β) . umean [START_REF] Ben-Ameur | Routing of Uncertain Traffic Demands[END_REF] subject to:

l∈L k∈N p∈P k 1 c l λ p l . r k p . x(k) umean.q ∀ X ∈ X k∈N p∈P k λ p l . r k p . x(k) umax.c l ∀ l ∈ L, ∀ X ∈ X p∈P (k) r k p = 1 ∀ k ∈ N r k p 0 ∀ p ∈ P k , ∀ k ∈ N
combination of u max and u mean as the new objective function u aof = β . u max + (1 -β) . u mean , where 0 β 1 is the combination fraction. Despite its simple form, this new objective is very effective and provides accurate results for both performance indicators. The new obtained optimization problem is the Robust Routing AOF Optimization Problem (RRAP), defined in [START_REF] Ben-Ameur | Routing of Uncertain Traffic Demands[END_REF]. Once again, problem ( 8) is solved with the same algorithms used in [START_REF] Ahuja | Network Flows: Theory, Algorithms, and Applications[END_REF].

B. Comparison between RRMP and RRAP

We will now evaluate both the RRMP and RRAP versions of SRR in the same traffic scenario previously used in section II-A. In order to appreciate the dependence of RRMP on the maximum link utilization threshold u thres max , two different thresholds are used in the evaluation: u thres max1 = 1, which corresponds to the constraint u max 1 in (3), and u thres max2 = u robust max , where u robust max is the output of RROP HTL in section II-A. In the case of RRAP, the weight β is set to 0.5, namely an even balance between u max and u mean . Figures 5 and6 depicts the comparison as regards (a) maximum link utilization and (b) mean end-to-end queuing delay. Let us focus the attention on the operation after the 100th minute, as all robust routing configurations use X HTL as uncertainty set. To be as fair as possible, both RRMP and RRAP use the same set of paths as those used by RROP in figure 3. RRMP performance clearly depends on the threshold u thres max ; in the first case, the attained maximum link utilization is well beyond the optimal values, reaching almost a 70% of relative performance degradation. This overload directly translates into huge mean end-to-end queuing delays. Results are quite impressive when considering the second threshold, both as regards u max and d mean . RRMP u thres max2 provides a highly efficient robust routing configuration, showing that it is possible to improve current implementations of SRR with a slight modification of the objective function. However, this dependence on the threshold u thres max introduces a new tunable parameter, something undesirable when looking for solutions that simplify network management.

As regards RRAP, obtained results are slightly worse than those obtained by RRMP u thres max2 , but still very close to the optimal performance, with a relative performance degradation of about 10% as regards u max and d mean w.r.t. an optimal routing configuration. Nevertheless, RRAP has no tunable parameter apart from the combination factor β, which in fact is set to a half independently of the traffic situation.

C. The Reactive Robust Routing

As we showed in section II-A, the definition of the uncertainty set has a major impact on the performance of SRR. In particular, we saw that using a single definition of uncertainty set under highly variable traffic can not provide routing efficiency for both normal operation traffic and unexpected traffic events. Despite being one of its most important features, using a single SRR configuration is not the best strategy.

In [START_REF] Casas | Robust and Reactive Traffic Engineering for Dynamic Traffic Demands[END_REF] we proposed an adaptive version of SRR, known as the Reactive Robust Routing (RRR). The basic idea in RRR consists of computing a primal robust routing configuration R o robust for expected traffic variations in normal operation within a primal polytope X o . This polytope is defined as in section II, based on a certain fixed routing configuration R o and the expected links traffic load we shall call

Y o = {ρ o i }.
Additionally, a set of m anomaly polytopes X j are defined, and a preemptive robust routing configuration R j robust is computed for each of these anomaly polytopes.

Let us explain the concept of an anomaly polytope. In figure 3, the abrupt increase in traffic volume is caused by a single anomalous OD flow x k that unexpectedly carries a many times bigger traffic load θ due to an external routing modification. After this exogenous unexpected event, the normal operation traffic demand X takes the value X ' = X + θ.δ k , where

δ k = (δ 1,k , . . . , δ k,k , . . . , δ m,k ) T , δ i,k = 0 if i = k and δ k,k = 1.
We shall designate this unexpected traffic increase in OD flow x k as anomalous traffic event A k . The anomaly polytope X k results from expanding the primal polytope X o in the directions of the links that traverses the anomalous OD flow x k , with respect to R o . The reader should bear in mind that the kind of unexpected traffic events we deal with are independent of the intradomain routing; these events originate outside the network and propagate between origin-destination nodes. This justifies the relevance of the polytope expansion with respect to R o . The obtained polytope X k is the smallest polytope that contains the unexpected traffic demand X ' and thus, the corresponding robust routing configuration R k robust provides accurate performance under its occurrence. Figure 7 explains the idea of the multiple anomaly polytope expansion. As before, r i o stands for the i-th row of the routing matrix R o . Note that in a real scenario it is not possible to predict the size of the anomalous traffic θ. As a consequence, the primal polytope X o is expanded to the limits of link capacities, obtaining the following anomaly polytope for each anomalous traffic event A k :

Xo X1 X2 Xj Xm r i o r i o r i o r 2 o r 2 o r q-1 o r 1 o r q o r q o r 3 o c i c i c i c 2 c 2 c q-1 c 1 cq cq c 3
X k = X ∈ R m , R o .X Y A k , X 0 , ∀k ∈ N (9)
In ( 9), the i-th component of

Y A k takes the value ρ oi if r i,k o = 0, or the value c i if r i,k o > 0, being r i,k o the element (i, k) of R o .
Given the primal and the m preemptive robust routing configurations R o robust and R j robust , RRR uses different anomaly detection/localization sequential algorithms to detect the occurrence of an anomalous event A k , switching routing from R o robust to R k robust . In particular, RRR defines a recursive anomaly detection function we shall call g(t), and m recursive anomaly localization functions s j (t), j ∈ N . These functions are updated with every new traffic measurement at time t, and when the difference between g(t) and g(t-1) exceeds a certain anomaly detection threshold ∆ detection 1 an anomaly is declared. Functions s j (t) are used to locate the anomalous OD flow; under the presence of an anomaly in OD flow x k at time t, the function s k (t) dramatically increases its value w.r.t. the rest of the s j (t) j =k functions, highlighting the anomalous OD flow. Additionally, RRR uses a similar detection algorithm to detect the end of the anomaly at OD flow x k , switching routing from R j robust to R o robust when normal operation is regained. We refer the reader to [START_REF] Casas | Robust and Reactive Traffic Engineering for Dynamic Traffic Demands[END_REF] for additional details on the implementation of RRR.

RRR can handle large and unexpected traffic variations in single OD flows quite effectively (the case of multiple simultaneous anomalies is beyond the scope of RRR). However, given the difficulty involved in modifying the routing configuration of a large scale network in an on-line fashion, the contributions of RRR are mainly theoretical. This problem can be solved by using a load balancing technique instead of a complete routing reconfiguration. In load balancing, we keep the same set of paths P k for each OD pair k, and only modify the fractions of traffic sent through each path. Load balancing can be easily performed on-line and does not require any additional modifications in current path-based networks such as MPLS. We shall refer to the load balancing variant of RRR as Reactive Robust Load Balancing (RRLB), stressing the difference between routing reconfiguration and load balancing.

RRLB uses the same set of anomaly polytopes X j defined in RRR, but the computation of the m preemptive robust routing configurations R j robust is slightly modified. Parameter φ is used as an anomaly flag variable that takes value 0 if traffic is in normal operation and value k under the occurrence of an anomalous traffic event A k . The detection function h φ (t) is used to detect the end of the anomalous traffic event in OD flow φ. Function h φ (t) takes negative values during the presence of an anomaly in OD flow φ at time t, increasing its value above a certain detection threshold ∆ detection 2 when the anomaly has ended. Both detection thresh- olds ∆ detection1 and ∆ detection2 are tuned in order to maximize the detection probability while bounding the false alarm rate, see [START_REF] Casas | Robust and Reactive Traffic Engineering for Dynamic Traffic Demands[END_REF] for additional details.

IV. EVALUATION AND DISCUSSION

In this section we evaluate the performance of the Robust Routing enhanced mechanisms presented in this work, considering both normal operation and anomalous traffic situations. This allows for performance comparison at different levels of traffic variability. As both RRAP and RRMP provide similar results (when u thres max is correctly defined for RRMP), we will only consider the RRAP mechanism in the evaluation.

For the normal operation scenario, we shall compare the efficiency of RRAP against the traditional RROP mechanism. As regards the anomalous traffic scenario, we consider the same situation depicted in figure 3 and compare the execution of RRLB for both RRAP and RROP. In the evaluation, we shall use RRLB-OP and RRLB-AP to designate the Reactive Routing Load Balancing variants of RROP and RRAP respectively.

A. Normal Operation Traffic Scenario

The first case-scenario corresponds to traffic in normal operation. The only variability is due to typical daily fluctuations. Figure 8 presents the evolution of u max and d mean for RROP and RRAP, using a set of 260 TMs from dataset X01 in [START_REF] Zhang | Abilene Dataset 04[END_REF]. Both mechanisms perform similarly as regards maximum link utilization, depicted in figure 8(a). This may be further appreciated in the boxplot summaries presented in figure 9(a), where values are relative to those obtained with an optimal routing configuration. Note that the relative performance degradation is around 10% in both cases.

Figures 8(b) and 9(b) show that results are quite different as regards mean queuing delay. While RRAP has a relative degradation smaller than 10% w.r.t. the optimal delay, RROP systematically obtains an important difference, attaining a performance degradation close to 40%. These results further highlight the limitations of RROP as previously discussed: using u max as a performance objective results in a relatively low maximum utilization, but neglects the rest of the links, impacting the network-wide performance.

B. Anomalous Traffic Scenario

The second case-scenario is the one considered in section II-A, where there is a sudden and abrupt increase of the traffic volume carried by one OD flow. Both RRLB-OP and RRLB-AP use the RRLB mechanism previously described to adapt traffic balancing after the detection of the anomalous traffic variation. As a difference with respect to the evaluation in figure 3, where traffic was assumed known in advance, this case-scenario corresponds to a real situation where traffic anomalies can not be forecast.

Figure 10(a) depicts the attained maximum link utilization before (LTL period) and after the anomalous event (HTL period). Let us first discuss the execution of RRLB, comparing the performance obtained with RRLB-OP against the one obtained by RROP LTL and RROP HTL in figures 3(a) and 4(a). Figure 12 provides a boxplot summary of the relative maximum link utilization (w.r.t. the optimal values) obtained with the 3 algorithms during the (a) LTL period and (b) HTL period. The reader should remember that RROP LTL represents the optimal robust routing configuration for the LTL period as regards RROP, and similarly RROP HTL during the HTL period. It is interesting to see that RRLB-OP obtains a slightly worse performance than the optimum robust routing configurations during the whole evaluation period, with a difference close to 5%. Nevertheless, RRLB-OP represents a real situation where traffic is not assumed known in advance.

Regarding the comparison between RRLB-OP and RRLB-AP, figure 11 shows that both algorithms obtain similar results as regards maximum link utilization, with a relative performance degradation smaller than 15% during the whole evaluation period. Note that while important, this performance degradation is surprisingly small if we consider that traffic increases more than 500% in less than 10 minutes. As regards end-toend queuing delay, RRLB-AP clearly outperforms RRLB-OP, achieving a relative mean queuing delay almost 30% smaller. These results reinforces once again our observations about the difficulty in RROP to attain global performance, and the advantages of using a simple network-wide objective function in a robust routing algorithm.

V. CONCLUSIONS AND FUTURE WORK

In this work we have studied the problem of intradomain routing optimization under highly variable and difficult to predict traffic demands. We have presented a comprehensive analysis of a plausible solution to the problem, namely the Robust Routing approach, evaluating its performance under different traffic scenarios. From this analysis we have obtained tangible evidence to highlight two important shortcomings of current Robust Routing implementations. On the one hand, we saw that using a local performance criterion such as the maximum link utilization (MLU) does not provide a suitable RROP LTL RROP HTL RRLB-OP objective function as regards network-wide performance and QoS provisioning. In particular, we showed that an almost optimal robust routing configuration w.r.t. MLU can experience rather high mean end-to-end queuing delays, thus impacting on delay sensitive traffic. The maximum link utilization is widely used in current network optimization problems, particularly in most Robust Routing proposals, thus we believe that this simple evidence can help to enhanced future implementations.

On the other hand we have shown that using a single routing configuration is not a cost-effective solution when traffic is relatively dynamic. Stable Robust Routing obtains quite poor performance either when faced with non considered traffic demands (tight uncertainty sets) or when designed to manage as many traffic demands as possible (big uncertainty sets). It is clear from our study that some form of dynamism is necessary.

We have proposed solutions to both detected problems of current Robust Routing implementations. To begin with, we have shown that objective optimization functions can be kept simple yet better network-wide performance can be attained. By using a simple combination of performance indicators such as the maximum and the mean link utilization, we obtained a robust routing configuration that definitely outperforms current implementations from a global end-to-end perspective, while achieving almost identical results as regards worst-case link utilization. Additionally, we have extended a previous proposal of our own that correctly manages normal and anomalous traffic situations, improving its execution towards a practical implementation. The Reactive Robust Load Balancing introduces a dynamic approach to deal with unexpected traffic events, balancing load between prestablished paths only when large traffic variations are detected. The evaluations showed than combining this reactive scheme with network-wide performance objectives can offer proper solutions to deal with current and future traffic scenarios.

The framework of Aggregated Objective Functions (AOF) provides interesting results as regards multi-objective optimization, particularly in the context of robust optimization. An AOF approach can be used to construct better objective functions from simple performance indicators and thus deserves further analysis. The estimated queuing delay function provides a difficult to optimize objective, but the problem can be solved provided certain realistic simplifications; this is currently part of our ongoing work.

Fig. 1 .

 1 Fig. 1. The uncertainty set X as a polytope.

Fig. 3 .

 3 Fig. 3. Maximum link utilization and mean end-to-end queuing delay. Traffic demand volume abruptly increases after the 100th minute.

Fig. 4 .

 4 Fig. 4. Maximum link utilization and mean end-to-end queuing delay, boxplot performance summary. Depicted results are relative to the optimal values.

Fig. 5 .Fig. 6 .

 56 Fig. 5. Maximum link utilization and mean end-to-end queuing delay for RRMP and RRAP.

Fig. 7 .

 7 Fig. 7. Different anomaly polytopes for preemptive robust routing computation.

Fig. 8 .

 8 Fig. 8. Maximum link utilization and mean end-to-end queuing delay under normal operation. RROP and RRAP are compared against the optimal values.

Fig. 9 .

 9 Fig. 9. Maximum link utilization and mean end-to-end queuing delay under normal operation, boxplot performance summary. Depicted results are relative to the optimal values.

Fig. 10 .Fig. 11 .

 1011 Fig. 10. Maximum link utilization and mean end-to-end queuing delay under anomalous traffic. RRLB-OP and RRLB-AP are compared against the optimal values.

Fig. 12 .

 12 Fig. 12. Maximum link utilization for RROP LTL, RROP HTL and RRLB-OP, boxplot performance summary. Depicted results are relative to the optimal values.

  The same set of paths P k obtained during the computation of R o robust is used in every R j robust . As it was done in sections II-A and III-B, routing configurations R j robust are obtained with a simplified version of the former optimization algorithm, where only new traffic demands are progressively added and no extra paths are created. The following schema gives a high-level description of RRLB: : compute R o robust and P k , ∀k ∈ N for primal polytope Xo 2: using P k , compute R j robust , ∀ Xj , ∀j ∈ N 3: set anomaly Flag φ ← 0 4: for t = 1:∞ do

	Algorithm 1 Reactive Robust Load Balancing (RRLB)
	5:	if φ = 0 then
	6:	update anomaly detection function g(t)
	7:	for i = 1:m do
	8:	update anomaly localization function si(t)
	9:	end for
	10:	if g(t) -g(t -1) > ∆detection 1 then
	11:	j = argmax	sj(t)
		j∈N
	12: 13:	balance traffic according to R j robust φ ← j
	14:	end if
	15:	else
	16:	update anomaly-end detection function h φ (t)
	17:	if h φ (t) > ∆detection 2 then
	18:	balance traffic according to R o robust
	19:	φ ← 0
	20:	end if
	21:	end if
	22: end for

1