Antoine Loiseau
email: antoine.loiseau@cea.fr

Jacques J A Fournier
email: jacques.fournier@cea.fr

Binary Edwards Curves for intrinsically secure ECC implementations for the IoT

Keywords: IoT, Elliptic Curves Cryptography, Binary Edwards Curves

Even if recent advances in public key cryptography tend to focus on algorithms able to survive the post quantum era. At present, there is a urgent need to propose fast, low power and securely implemented cryptography to address the immediate security challenges of the IoT. In this document, we present a new set of Binary Edwards Curves which have been defined to achieve the highest security levels (up to 284-bit security level) and whose parameters have been defined to fit IoT devices embedding 32-bit general purpose processors. We optimized the choice of the point generator with the w-coordinate to save a multiplication in the addition and doubling formulae. We manage to compute one step of the Montgomery Ladder in 4 multiplications and 4 squares. On top of the performance benefits, cryptography over such curves have some intrinsic security properties against physical attacks.

Introduction

The rising amount of connectivity and number of connected devices (with the associated services) has given rise to the concept of the "Internet of Things". The technical challenges pertaining to the IoT are daunting [START_REF] Stankovic | Research directions for the internet of things[END_REF]. Among those challenges, there is the management of a massively huge network (hundreds of billions by the next decade) and the even bigger amount of associated data to manage, store and analyse, the heterogeneity of today's different systems that have to co-exist, the resilience of such systems, the openness of such an infrastructure and, last but not least, security & privacy [START_REF] Fournier | Security Challenges facing the IoT[END_REF]. Security and privacy have already been identified as vital issues of the IoT and as key to its adoption and deployment [START_REF] Skarmeta | Internet of Things[END_REF]) [START_REF] Rubens | Internet of Things a Potential Security Disaster[END_REF]. The main difficulties with deploying security for the IoT come from the heterogeneity of this "system of systems", lack of standard for security and inter-device communications, lack of offthe-shelf trusted IoT devices and tools (respecting all power, size and security constraints) and lack of appropriate business model ('low end' devices will potentially be communicating with 'high end' ones. . . ).

One of the key aspects for such systems is that end-to-end security has to be enforced, which means that the end-nodes of an IoT infrastructure shall also embed security features and applications. Asymmetric cryptography is one of such security tools and given the size, power and performance constraints pertaining to IoT end-nodes (which are usually lowend devices with limited resources) Elliptic Curves Cryptography (ECC) is by far the best fit for such use cases and this despite the recently growing fear of quantum computers: since the latter are thought to be able one day to efficiently solve the Discrete Logarithm Problem (DLP), post-quantum cryptosystems, like lattice based cryptography or cryptography based on isogenies are currently being studied. However, these cryptosystems as they exist today, cannot be used in the constrained IoT nodes due to the large keys and long computing times. ECC is still, in our opinion, the best short and mid term answer to address the immediate security needs for the IoT.

Elliptic Curve Cryptography for the IoT

ECC was introduced by Koblitz and Miller in [START_REF] Koblitz | Elliptic curve cryptosystems[END_REF] [START_REF] Miller | Use of Elliptic-Curves in Cryptography[END_REF] where the Elliptic Curves Discrete Logarithm Problem (ECDLP) is presented and its exploitation for doing cryptography over elliptic curves is explained. The main standardization effort for ECC has been through the FIPS 186-4 (Information Technology Laboratory, 2013) by the NIST. Alternative approaches have been proposed like [START_REF] Brainpool | ECC Brainpool Standard Curves and Curves Generation v1[END_REF] which introduces Brainpool curves, and more recently the Edwards curve introduced in [START_REF] Bernstein | High-speed high-security signatures[END_REF] by Bernstein et al. based on works by Edwards [START_REF] Edwards | A normal form for elliptic curves[END_REF]. The latter pieces of work open a new vision in the landscape of applied elliptic curve cryptography.

The so-called "NIST curves" are largely deployed in several devices with different goals in today's IT systems ranging from low power devices to modern computers. However these curves were proposed in 1999 where "low power devices" were not as common as today and hence these curves were not designed for today's IoT constraints. Elliptic Curves usually used for cryptography are defined by the Weierstrass equation

y 2 + a 1 xy + a 3 y = x 3 + a 2 x 2 + a 4 x + a 6
and defined over a prime field with large prime characteristic or over binary extension fields of characteristic 2. The underlying mathematics in each kind of field have different arithmetics, each with its own advantages and drawbacks.

To start with, the arithmetics of the Weiesrtrass form are not adapted to the low resource constraints of IoT devices. Moreover those curves were neither designed nor chosen with the additional constraint of being resistant to side channel attacks and fault attacks: for example, due to the incomplete law group of these NIST curves, we have one formula for adding points and another formula for doubling points (except for the specific case Weierstrass curves over prime fields with p > 3 where we have a new complete formula for the addition [START_REF] Renes | Complete addition formulas for prime order elliptic curves[END_REF]). Such constructions introduce in an ECC implementation weaknesses exploited in Simple Power Analysis (SPA) or Timing attacks (Fan et al., 2010a). In such cases, adding countermeasures against side channel attacks severely downgrade the performance and efficiency of such ECC implementations, which can be a heavy drawback for the IoT end-nodes. This is why it would be more interesting to look for curves that could avoid such weaknesses. For example curves like Curve25519, [START_REF] Bernstein | High-speed high-security signatures[END_REF] [START_REF] Edwards | A normal form for elliptic curves[END_REF] have a complete group law providing one same formula for point addition and point doubling.

Another important aspect for an ECC implementation is the choice between curves over F p or F 2 m . Curves over large prime characteristics are usually preferred to those over binary fields due to claims that there may be sub-exponential algorithms for solving the discrete logarithm problem for elliptic curves on binary fields. However even if the best optimistic conjectures of the DLP claim to have a complexity of L 2 n ( 23 )1 or even L 2 n ( 1 2 ), [START_REF] Galbraith | Recent progress on the elliptic curve discrete logarithm problem[END_REF]) recalls that they are still "hypotheses" and that recent experiments "raise the doubt about the subexponentiality claims". One of the main advantages of arithmetics in binary fields is the carry-less operation because in practice managing carry propagations as needed in operations over prime fields downgrades computation times especially when it comes to implementing dedicated hardware accelerators. Moreover carry propagations can introduce weaknesses against side channel attacks [START_REF] Fouque | The Carry Leakage on the Randomized Exponent Countermeasure[END_REF].

We have a large set of models on large prime fields [START_REF] Doche | Efficient scalar multiplication by isogeny decompositions[END_REF] [START_REF] Edwards | A normal form for elliptic curves[END_REF] (Castryck et al., ) [START_REF] Bernstein | Twisted hessian curves[END_REF] [START_REF] Hisil | Jacobi quartic curves revisited[END_REF]) [START_REF] Feng | Twisted Jacobi intersections curves[END_REF] (Moloney et al., ) (Bernstein et al., 2008a).

Our approach and contribution

In [START_REF] Edwards | A normal form for elliptic curves[END_REF], Harold Edwards introduced a very efficient model with a complete law group:

x 2 + y 2 = 1 + dx 2 y 2
However, with the aforementioned drawbacks of working on curves in large prime fields, Bernstein et al. [START_REF] Bernstein | Binary edwards curves[END_REF] introduced the binary version of Edwards curves with also the advantage of a complete law and obviously without the carry propagation problems :

d 1 (x + y) + d 2 (x 2 + y 2 ) = xy + xy(x + y) + x 2 y 2
Such a Binary Edwards Curve is, given the constraints of IoT devices, one of the best choices for implementing ECC, resistant to physical attacks, in current and future embedded devices : they are fast, have a complete group law and the underlying carryless arithmetics would require low power and compact dedicated hardware accelerators.

So far, only one such BEC has been proposed in the literature, which does not even offer what is the admitted reference of "128-bit security level". Hence our approach has been to generate other BECs, with higher security levels and with parameters tailored to fit the hardware constraints of IoT devices. To build those new elliptic curves, we split the parameters of the curve into three parts, each of which has been optimized separately. The first part is the finite field representation. The second is the curve's parameters which depend on the finite field representation and elliptic curve model chosen. The third part is the choice of the point generator of the elliptic curve group.

In this paper we provide the newly-generated BECs after explaining how they have been generated, i.e. by detailing the constraints and optimisations used in each of the three steps of the BECs generation. This new set of curves targets IoT implementations. Each curve was chosen to achieve from 112-bit to 284-bit security levels and the associated parameters have been chosen to be efficiently suited for hardware architectures working on data paths of 32 bits or 64 bits like those traditionally found in off-the-shelf IoT devices.

The methodology used for generating the curves and their associated parameters is depicted in Figure 1. We split security requirements in two sets. The first refers to criteria on the number of points of the curve and its twist, which is critical in the selection of a new elliptic curves. The second set of criteria contains all remaining criteria in 3.3. These criteria are checked after the construction of a optimized point generator.

The rest of this paper is organised as follows. Section 2 presents the concept of friendly polynomials from [START_REF] Scott | Optimal Irreducible Polynomials for GF(2 m ) Arithmetic[END_REF] and the associated advantages. We select a set of new binary fields to target different security levels to allow a good security/efficiency tradeoff. Section 3 introduces the Binary Elliptic Curves as defined in [START_REF] Bernstein | Binary edwards curves[END_REF] with a focus on those curves having efficient arithmetics. We also review the different security requirements to select a BEC. Section 4 details how we built the generator of the elliptic curve subgroup used and how the choices made provide performance gains. Section 5 summarizes the expected intrinsic security properties of such BEC-based implementations against physical attacks. Some concluding remarks and perspectives are finally provided in the last section of this paper. Finally, the section 6 provides a first overview of performances of this kind of curves. The complete set of newly generated curves and parameters are given in the appendix.

Optimal finite fields for efficient binary field arithmetics

The first step for generating a new set of elliptic curves for cryptography is to build "friendly" (i.e. for which calculations can be done rapidly) finite fields with an optimal arithmetic (i.e. for which some calculation "tricks" can be used for fast and compact implementations) for the targeted IoT devices. Arithmetics over elliptic curves depend on the representation of the base finite field which in turn depends on the choice of the modulus, represented under a polynomial form, used to define the finite field. If we want Trinomials and pentanomials generation with Scott's method [START_REF] Scott | Optimal Irreducible Polynomials for GF(2 m ) Arithmetic[END_REF] BECs generation with Bernstein's method [START_REF] Bernstein | Binary edwards curves[END_REF] Finding efficient group generators Selection and validation of BECs with security criteria Architecture constraints (data width) First security criteria [START_REF] Bernstein | Batch binary Edwards[END_REF] Second security criteria [START_REF] Bernstein | Batch binary Edwards[END_REF]

Modulus

BECs

Figure 1: Overview of the generation scheme of BEC.

to consider the constraints due to the hardware platform used for the implementation (in our case, an IoT device embedding a 32-bit processor), choosing an irreducible polynomial at random, i.e. with a random number of monomials, seems the wrong approach if we want the best trade-off between performance and security. The main operation in such implementations is the reduction of a finite field element by the modulus and if we use a random polynomial, this operation could be expensive in terms of execution time and power consumption. A choice usually made for polynomials in binary fields F 2 n is a trinomial or a pentanomial of the type:

p(x) = x n + x a + 1, n > a > 0 p(x) = x n + x a + x b + x c + 1, n > a > b > c > 0
The NIST (Information Technology Laboratory, 2013) for example selected its polynomials with the ANSI X9.62 rules also provided by Certicom [START_REF] Johnson | The Elliptic Curve Digital Signature Algorithm (ECDSA)[END_REF]. A trinomial with the minimal coefficient a shall be selected if it exists. If no trinomial exists a pentanomial is chosen with minimal coefficients a, b and c. However, as suggested in [START_REF] Scott | Optimal Irreducible Polynomials for GF(2 m ) Arithmetic[END_REF], the construction of the binary field is important for the efficiency of the implementation and this depends on the architecture used. Let w be the length in terms of bits of the registers of the targeted architecture, Scott [START_REF] Scott | Optimal Irreducible Polynomials for GF(2 m ) Arithmetic[END_REF] then defines three kinds of polynomials:

• Lucky Trinomial: A trinomial with na ≡ 0 mod w.

• Lucky Pentanomial: A pentanomial with na ≡ 0 mod w, nb ≡ 0 mod w and nc ≡ 0 mod w.

• Fortunate Pentanomial: A pentanomial with at least two of these equations are verified : n-a ≡ 0 mod w, nb ≡ 0 mod w and nc ≡ 0 mod w.

A lucky trinomial is rare and a lucky pentanomial does not exit if 4|n or n ≡ ±3 mod 8. These kinds of polynomials can save several shift operations compared to a classical polynomial used in cryptography. With a simple script in Sage (The Sage Developers, 2017) we searched all lucky trinomials and pentanomials for fields of bit lengths between 163 and 571: we found only 6 lucky trinomials compliant with a 32-bit architecture, i.e. whose representation is a multiple of 32. As for pentanomials fitted for 32-bit architectures, we have a larger number of them over several finite fields.

In the case of the field F 2 257 , we found the lucky trinomial x 257 + x 65 + 1, but on a 32-bit architecture (or any other architecture actually), a field element shall be stored onto 9 registers with only one bit in the last (uppermost) register. We might "waste" many store and load operations for only one bit. Such choices could have severe drawbacks for implementation of ECC in IoT devices. So, we select our own fields and polynomials to maximize the number of bits in the last register so as to minimize the number of load and store instructions used (even if we later keep in the end F 2 257 as intermediate field between F 2 223 and F 2 313 ).

We first selected a set of lucky trinomials for 32 architectures with different sizes to target different security levels. We split our searches/findings into layers of minimum security levels of 112,128,192,224 and 256 bits.

We also realised that most of those polynomials would also be compatible with 64-bit platforms (in case we wanted to anticipate for the day when such platforms would be available for IoT devices), except for the trinomials of F 2 479 and F 2 521 . So we did the additional effort of looking for polynomials compatible with 64-bit data widths for the latter levels of security: hence we ended up with two additional pentanomials in F 2 487 and F 2 569 .

The list of polynomials hence generated is as follows:

F 2 223 : x 223 + x 159 + 1 F 2 257 : x 257 + x 65 + 1 F 2 313 : x 313 + x 121 + 1 F 2 431 : x 431 + x 303 + x 239 + x 111 + 1 F 2 479 : x 479 + x 255 + 1 F 2 487 : x 487 + x 295 + x 167 + x 39 + 1 F 2 521 : x 521 + x 489 + 1 F 2 569 : x 569 + x 441 + x 313 + x 121 + 1
Table 1 compares these polynomials with those from the NIST in terms of number of logical operations (xor, shift, and) for a general purpose 32-bit architecture. We targeted a generic 32-bit architecture (excluding ARM-like architectures like where preshift operations are available to replace the sequence of logical instructions XOR-AND).

At this stage, we have selected a set of finite fields to match different security levels from 112 to 256 bits. With these selected binary fields we can now build the corresponding Binary Edwards Curves.

Binary Edwards Curves

Among the recently proposed models of elliptic curves for cryptography, Binary Edwards Curves (BEC) exhibit several interesting properties, [START_REF] Bernstein | Binary edwards curves[END_REF] [START_REF] Bernstein | Batch binary Edwards[END_REF]. First, these curves favour implementations that will be intrinsically safe against a set of side channel and fault attacks (see Section 5). In fact, these curves have a complete group law and can be used with the Montgomery ladder algorithm. Moreover, the binary field arithmetics avoid leakages exploiting the carry propagation. Second, these curves allow an efficient arithmetic when using the correct coordinates' representations and fast arithmetics in binary fields. Thanks to the latter characteristics, BECs might, in our opinion, be interesting candidates for IoT applications.

Construction and definition

In [START_REF] Edwards | A normal form for elliptic curves[END_REF], Harold Edwards introduced a new model for elliptic curves allowing efficient computations for the elliptic curves' arithmetics. The main interesting feature for us is the complete group law on this curve giving one same formula for addition and doubling . The Edwards equation is:

x 2 + y 2 = 1 + dx 2 y 2
E d 1 ,d 2 : d 1 (x + y) + d 2 (x 2 + y 2 ) = xy + xy(x + y) + x 2 y 2
From this definition, we have the following interesting proprieties. This curve is symmetric and we have ∀P ∈ E(K), P = (x, y) ⇒ -P = (y, x). With this assumption, we define the neutral point of the group law by (0, 0). Definition 2. Projective closure: The projective closure of the BEC E d 1 ,d 2 is:

d 1 (X +Y )Z 3 +d 2 (X 2 +Y 2 )Z 2 = XY Z 2 +XY (X +Y )Z +X 2 Y 2
The last important aspect of this construction is that each BEC is birationally equivalent to a Weierstrass form:

v 2 + uv = u 3 + (d 2 1 + d 2 )u 2 + d 4 1 (d 4 1 + d 2 1 + d 2 2 )
We can hence define the mapping φ(x, y) → (u, v) for this equivalence by:

u = d 1 (d 2 1 + d 1 + d 2 ) x + y xy + d 1 (x + y) v = d 1 (d 2 1 + d 1 + d 2 ) x xy + d 1 (x + y) + d 1 + 1
And the inverse mapping φ -1 (u, v) → (x, y) by:

x = d 1 (u + d 2 1 + d 1 + d 2 ) u + v + (d 2 1 + d 1 )(d 2 1 + d 1 + d 2 ) y = d 1 (u + d 2 1 + d 1 + d 2 ) v + (d 2 1 + d 1 )(d 2 1 + d 1 + d 2 )

Arithmetics on BEC

There are special cases where the arithmetics associated with the BEC can be optimised. For the case where d 1 = d 2 we have the following equation for a BEC:

d(x + y + x 2 + y 2 ) = xy + xy(x + y) + x 2 y 2 (1)
with the condition t ∈ K ⇒ d = t 2 +t. In the following sections, we consider this particular case of BECs.

Addition and doubling: Efficient arithmetics on elliptic curves depend on the choice for the coordinates' representation. In the case of BECs, the projective wcoordinate with differential additions and doublings [START_REF] Bernstein | Binary edwards curves[END_REF] is the fastest way for doing points addition and point doubling. The idea of the differential w-coordinate is to represent a point P(x, y) of the curve by w(P) = x + y. In this way, we can compute w(2P) and w(Q + P) given w(P), w(Q) and w(Q -P). Let P and Q be two points of the curve, and let w(P) = W 2 /Z 2 , w(Q) = W 3 /Z 3 and w(Q -P) = w 1 , we compute w(2P) = W 4 /Z 4 and w(Q + P) = W 5 /Z 5 with Algorithm 1. This kind of coordinates is well suited for implementing the Montgomery ladder (see Algorithm 3) to compute the exponentiation of a point.

For an "Add and Double operation" we have only 5 field multiplications, 4 squares and 2 multiplications by d. The Montgomery ladder computes w(kP) and w(kP + P) with an addition and a doubling at each step. The use of w-coordinates appears to be the best way to compute scalar exponentiations. In the Algorithm 1 for the addition and doubling we use a scalar w 1 = w(Q -P). In the case of the Montgomery ladder (Joye and [START_REF] Yen | The Montgomery powering ladder[END_REF] as given in Algorithm Algorithm 1 w-coordinates Adding and Doubling

Require: W 2 , Z 2 ,W 3 , Z 3 , w 1 1: C ← W 2 .(Z 2 +W 2 ) 2: V ← C.W 3 .(Z 3 +W 3 ) 3: W 4 ← C 2 4: Z 4 ← d(Z 2 2 ) 2 +W 4 5: Z 5 ← V + d(Z 2 .Z 3 ) 2 6: W 5 ← V + Z 5 .w 1 7: return W 4 , Z 4 ,W 5 , Z 5
3, the difference R 1 -R 0 is invariant at each step of the ladder. In other words the scalar w 1 is a constant scalar.

Recently, new formulae for w-differential adding and doubling have been proposed to save a multiplication by d [START_REF] Kim | Binary edwards curves revisited[END_REF] [START_REF] Koziel | Low-resource and fast binary edwards curves cryptography[END_REF]: an additional squaring can be saved with the Co-Z trick which performs

W 2 Z + W 3 Z = W 4 Z and 2 × W 2 Z = W 5 Z .
Finally, the addition and doubling operations can be written as :

W 4 Z = (W 2 +W 3 ) 2 + 1 w 1 (W 2 +W 3 ) 2 Z 2 + 1 w 1 (W 2 +W 3 ) 2 W 5 Z = (W 2 (W 2 + Z)) 2 dZ 4 + (W 2 (W 2 + Z)) 2
Algorithm 2 gives the computation of these coordinates. In this case the Adding and Doubling operation takes 5 multiplications, 1 multiplication by d and 4 squarings.

Algorithm 2 w-coordinates Adding and Doubling revisited with the Co-Z trick.

Require: W 2 ,W 3 , Z, 1 w 1 1: C ← (W 2 +W 3 ) 2 2: D ← Z 2 3: E ← 1 w 1 C 4: U ← E +C 5: V ← E + D 6: S ← (W 2 (Z +W 2 )) 2 7: T ← S + dD 2 8: W 4 ← UT 9: W 5 ← V S 10: Z ← V T 11: return W 4 ,W 5 , Z
With this kind of coordinates (with the Co-Z trick), we loose the completeness of the addition group law in favour of an "almost" complete addition law. In fact, if w(P-Q) = 0 we cannot define 1 w 1 used in the Algorithm 2. Recently, (Rezaeian Farashahi and Hosseini, 2017) defined a new kind of coordinates also "almost" complete with the same efficiency

Algorithm 3 Montgomery Ladder

Require: w(P), k = (k t-1 , ..., k 0 ) 2 1: R 0 ← O 2: R 1 ← P 3: for j = t -1 to 0 do 4:

if k j = 0 then 5:

R 1 ← R 0 + R 1 6: R 0 ← 2R 0 7: else 8: R 0 ← R 0 + R 1 9:
R 1 ← 2R 1 10:

end if 11: end for 12: return R 0 = w(kP), R 1 = w(kP + P)

as the mixed differential w-coordinates with Co-Z trick. However, (Rezaeian Farashahi and Hosseini, 2017) uses an advanced definition of w by w(x, y) = x+y d(x+y+1) and does not use the Co-Z trick. In this way, to recover affine coordinates from the Montgomery Ladder's output we have to compute two inversions against only one for mixed differential w-coordinates with Co-Z trick.

Recovering the affine coordinates from the wcoordinate: The Montgomery Ladder (Algorithm 3) computes w 2 = w(kP) and w 3 = w(kP + P). We can calculate the x coordinate of kP with the output of the Montgomery Ladder and with P = (x 1 , y 1 ) using the following equation from [START_REF] Bernstein | Binary edwards curves[END_REF]:

x 2 + x = w 3 (d + w 1 w 2 (1 + w 1 + w 2 ) + w 2 1 w 2 2 ) + d(w 1 + w 2 ) + (y 2 1 + y 1 )(w 2 1 + w 2 ) w 2 1 + w 1
This formula requires 1 inversion, 4 multiplications and 4 squares. From there we have to solve a equation x 2 + x = A if the trace of A is zero. For a fixed generator point, we can save the inversion operation if w 2 1 + w 1 is precomputed. We can also recover the y coordinate with

y 2 + y = d(x + x 2 ) d + x + x 2

Generating secure Binary Edwards Curves

The first issue we wanted to resolve was to have curves that provide the required level of security (i.e. at minimum at 128-bit level of security as recommended by the NIST since December 2015). So far, to the best of our knowledge, the only BEC curve proposed by Bernstein in [START_REF] Bernstein | Batch binary Edwards[END_REF]) has a security level of 125 bits.

The first things we started implementing were methods and tools to generate and propose a set of BECs with different field sizes for different security levels in order to match a targeted device's requirements and the best security/efficiency (power and performance) trade-offs. In [START_REF] Bernstein | Batch binary Edwards[END_REF] the authors enumerate the following security requirements for a BEC (and ECC in general):

• Prime degree of the finite field extension: The choice of the binary field extension is essential for fast computing but an extension field with a non-prime degree can be attacked by the GHS algorithm [START_REF] Gaudry | Constructive and destructive facets of weil descent on elliptic curves[END_REF]. Moreover a large prime degree has to be chosen due to the birthday paradox and has to provide a sufficient level security. We shall target a minimum security level of 128 bits. In this way we shall use fields with an extension degree of at least 256 bits.

• Number of points:

The curve E shall have a cardinal of 2 c p where p is a large prime and c shall be small (1, 2 or 3).

• Number of points on the Twist: We have the same requirement on the number of points of the twist of the curve E in order to ensure resistance against some fault attacks on the original curve.

• j-invariant: The j-invariant of the BEC, 1/d 8 , shall generate the extension field.

• Avoiding small discriminants: The discriminant ∆ E = Tr(E) 2 -4q shall be divisible by a large prime exactly once.

• Avoiding pairing attack: The multiplicative order k of 2 m modulo the large prime of |E| shall be large. We have the same condition with the prime of the twist. k is called the embedding degree of the subgroup of the elliptic curve. NIST advises in (Information Technology Laboratory, 2013) to have an embedding degree greater than 20. However with the recent improvements of index calculus attacks [START_REF] Galbraith | Recent progress on the elliptic curve discrete logarithm problem[END_REF], we have to increase this requirement. Brainpool curves [START_REF] Brainpool | ECC Brainpool Standard Curves and Curves Generation v1[END_REF] have an embedding degree greater than l-1 100 , where l is the prime order of big subgroup of the elliptic curve. This requirement is somehow "overrated", but we choose to check it anyway, as it is very easy for an elliptic curve with a large subgroup of prime order to satisfy this criterion.

Based on these requirements, we can check, for each d, if a given curve is secure. In practice, we start with parameters d with a small Hamming Weight, first with binomials, trinomials, quadrinomials and finally with pentanomials. Almost all tested curves fail to the "number of points" criterion. Hence the principal operation for generating & testing a BEC is to run a point counting algorithm. We implemented the AGM method [START_REF] Cohen | Handbook of elliptic and hyperelliptic curve cryptography[END_REF] to perform this counting. Moreover, to avoid counting the number of points on the twist, we use the relationship between a curve E and its twist E tw such that |E| + |E tw | = 2q + 2 to determine the number of points on the twist.

Using the above implemented algorithms, we managed to build new BECs on alternative fields with optimal finite field arithmetics and respecting all the security requirements. We implemented this algorithm in C language and based on the FLINT library [START_REF] Flint | FLINT: Fast Library for Number Theory[END_REF]. For our first trial, the program was run over an Intel i5-6200 at 2.40GHz over a Linux virtual machine on Windows. After two full days of running, we built the following BEC on F 2 313 :

F 2 313 [t] = F 2 [X] (X 313 + X 121 + 1) d = t 38 + t 33 + t 28 + 1 |E| = 2 2 × (2 311 -23801035639178335780622897191462018186156241745) |E tw | = 2 × (2 312 + 47602071278356671561245794382924036372312483491)
This curve verifies all the above security requirements as illustrated below:

• 313 is a prime number and avoids the GHS attack on F 2 313 .

• The number of points of E in F 313 is 4p 1 where p 1 is a large prime.

• The number of points of the quadratic twist of E, E tw , is 2p 2 where p 2 is a large prime.

• The j-invariant 1/d 8 generates F 2 313 .

• (2 313 + 1 -4p 1 ) 2 -2 315 is divisible by a large prime and avoids small discriminants.

(2 313 + 1 -4p 1 ) 2 -2 315 -1887935404880658993146706844859746511

• The multiplicative order of 2 313 modulo p 1 is not small, exactly (p 1 -1)/2. The multiplicative order of 2 313 modulo p 2 is not small, exactly p 2 -1.

Likewise other curves with different sizes have been generated as given in Appendix.

Optimal group generator

Now that we have set methods and tools for optimal (i.e. respecting our security and 32-bit oriented implementation constraints) field representations and parameter d for a BEC, we have now to choose a generator of the elliptic curve sub-group used to perform an exponentiation. In this case, the security of the protocol is based on the difficulty of the DLP in this sub-group. In fact DLP is hard to solve onto elliptic curves group. We selected a curve with a large sub-group, i.e. with a number of points divisible by a large prime. In our case we have curves with a number of points n = 2 c p with c = 1, 2 or 3 and p a large prime. So, we look for a point of the curve of order p.

We can select a random point G of order p but, in this case, the w-differential coordinate of G is a random scalar. However in the "Add and Double" formula of the BEC (Algorithm 2) we have a multiplication by the scalar 1 w 1 . Due to the property of the Montgomery Ladder, 1 w 1 is invariant at each step of the ladder. In this case, if we can find a generator G with a 1 w(G) having a small Hamming Weight, we can expect to reduce the number of multiplications for the addition and doubling steps from 5 to 4 multiplications.

If we consider w = x + y, then from the equation 1 of a BEC with a parameter d, we have a new equation of BEC:

d(w + w 2 ) = x 4 + (1 + w + w )x 2 + (w + w 2 )x (2)
So, for a chosen w we can find an x coordinate from the equation 2 and a y coordinate from the initial equation 1. We choose an inverse w-coordinate for G with a small Hamming Weight.

In the previous section we took as example the BEC over F 2 313 with modulus x 313 + x 121 + 1 and d = t 38 + t 33 + t 28 + 1. If we use this method to build a friendly generator for this curve we find a generator G with:

G x =0x15c67e3024c7c27466e72a3391256e9a729fc15 8092053d89087c0f38408b214b0ade57363ea938 G y =0x15c67e3024c7c27446e72a3391256e9a529fc15 8092053d8b087c0f38408b214b0ade57363ea938 G 1/w = 0x10000000000000001 1
w is reduced to t 64 +1. So, the multiplication by w 1 in the "Add and Double" formula, Algorithm 2, is simply done using 9 XOR operations.

Considering physical attacks

We are working on a software implementation of these new curves for a 32 bit RISC V architecture. The implementation has been done keeping in mind resistance to side channel and fault attacks, based on several factors: the intrinsic choices made during the curves' and parameters' generations, the algorithms used for implementing some of the basic operations and the additional countermeasures from the literature. For example the completeness proprieties of BECs and the use of the Montgomery ladder offer a first level of security against side channel attacks. Moreover we have hard-coded the d parameter of the curve and the 1 w parameter of the point generator to avoid any attack of invalid curve analysis and invalid point analysis. Against attacks like Correlation Power Analysis (CPA), Template Attacks, Horizontal Attacks or Differential Fault Attacks (DFA) we added countermeasures like randomization of the coordinates and of the scalar (Fan et al., 2010a). Table 2 summarizes the theoretical advantages of the BEC model with respect to Side Channel Attacks (SCA) and Fault Attacks (FA).

Performances

We implemented this new set of elliptic curves on a RISC V core running at 100 MHz. No special instruction has been added to the RISC V, which means that we have a software only implementation. We have decided to implement the differential w-coordinates system, as described above, with two cases, one with the optimized generator and the other with a random generator. We added the projective coordinates system for an interesting comparison and to show how the w-coordinates could be useful for ECC implementations. For the finite field multiplication in GF(2 n ) we use a López-Dahab multiplication [START_REF] Aranha | [END_REF]. The exponentiation computation is performed by using a Powering Montgomery Ladder (Joye and [START_REF] Yen | The Montgomery powering ladder[END_REF]. All functions for the finite field arithmetics are written in Assembly Language and functions for ECC arithmetics in C language. The performances are given in Table 3 and the parameters of these new set of BECs are provided in the appendix. These performances are for a computation of kP without countermeasures, even if as explained in Section 5 and in Table 2, we are theoretically secure against against a set of side channel attacks and fault attacks. First of all, the gap between projective coordinates and w-coordinates is huge, of the order of 355% between both coordinates system. The gain in performance between an optimized generator and a random generator in w-coordinate is around 19%. In terms of absolute time, this gain is low for small curves but for large curves this gain could be very important, almost 100ms for the largest curve. In

Physical Attacks

Intrinsic Resistance

Remaining vulnerability

Due to choice of parameters of BEC Due to implementation done Additional countermeasure implemented Timming Attack [START_REF] Kocher | Timing Attacks on Implementations of Diffie-Hellman, RSA, DSS, and Other Systems[END_REF] Unified arithmetics Use of Montgomery ladder -SPA [START_REF] Fan | State-of-the-art of secure ECC implementations: a survey on known sidechannel attacks and countermeasures[END_REF] Constant time programming DPA [START_REF] Kocher | Differential Power Analysis[END_REF] --

Randomization

Template Attack [START_REF] Herbst | Using templates to attack masked montgomery ladder implementations of modular exponentiation[END_REF] of coordinates [START_REF] Coron | Resistance against Differential Power Analysis for Elliptic Curve Cryptosystems[END_REF] Relative Doubling Attack [START_REF] Yen | Relative doubling attack against montgomery ladder[END_REF] --Blinded scalar [START_REF] Coron | Resistance against Differential Power Analysis for Elliptic Curve Cryptosystems[END_REF] RPA/ZPA [START_REF] Akishita | Zero-value point attacks on elliptic curve cryptosystem[END_REF] w-coordinates arithmetics Direct implementation of the generator -Carry-based Attack [START_REF] Fouque | The Carry Leakage on the Randomized Exponent Countermeasure[END_REF] Binary curves chosen --Horizontal Attack [START_REF] Clavier | Horizontal correlation analysis on exponentiation[END_REF] --Attack model to be tested Safe error [START_REF] Fan | State-of-the-art of secure ECC implementations: a survey on known sidechannel attacks and countermeasures[END_REF] -Use of Montgomery ladder -Invalid point analysis [START_REF] Biehl | Differential fault attacks on elliptic curve cryptosystems[END_REF] --

Verify that point in on the curve Invalid curve analysis [START_REF] Ciet | Elliptic curve cryptosystems in the presence of permanent and transient faults[END_REF] Curves' parameters on Twist for eg Direct implementation of curve parameter -Twist Attack (Fouque et al., 2008a) Curves' parameter for Twist --Differential Fault Attack [START_REF] Biehl | Differential fault attacks on elliptic curve cryptosystems[END_REF] --

Blinded scalar [START_REF] Coron | Resistance against Differential Power Analysis for Elliptic Curve Cryptosystems[END_REF] Table 2: SCA and Fault Attacks against BEC-based ECC our opinion, the next step for better performance is to add a carry-less multiplication instruction to the set of the RISC V instructions. In fact, the finite field multiplication is the bottleneck for the computation of kP.

It is hard to compare these results with other results in the literature. To the best of our knowledge, the performances presented in this paper are the first public results for a RISC V architecture. We can compare with performances of MbedTLS on an ARM M3 running at 96MHZ presented in [START_REF] Tschofenig | Performance of State-of-the-Art Cryptography on ARM-based Microprocessors[END_REF]. For a security level of 128 bits, an ECDSA signature (where the main operation is the kP operation) takes 122 ms with the standard secp256r1 curve, which is to be compared to 46 ms (for the kP operation) in our case. Note moreover that the performances in [START_REF] Tschofenig | Performance of State-of-the-Art Cryptography on ARM-based Microprocessors[END_REF] are given for a sliding window implementation, which is very efficient in terms of timing performances, but still vulnerable to side channel and fault attacks.

Conclusion

In this paper, we described how we generated a new set of Binary Edwards Curves which might be efficient, secure and low power alternatives for implementing elliptic curves cryptography in IoT devices. The optimization of the w-coordinate of the point gen-erator allows a speed up if these coordinates are used with the Montgomery Ladder. This new set of curves covers different security levels to offer the possibility of making the best trade-off between security and execution time. These curves are built with a 32-bit architecture in mind but they are also compliant with smaller architectures (8 and 16 bits) and larger ones (64 bits). Using Binary Edwards Curves might allow to have some intrinsic protections against a set of side channel and fault attacks. Future work will consist in validating the expected "physical security" properties in practice. As the BEC model allows a complete compliance with the Weierstrass model, we shall also study how the use of BEC curves can be compliant, in practice, with already deployed implementations of ECDSA on Weierstrass curves. We shall also work on the implementation of a dedicated hardware accelerator in order to take advantage of the carry-less operations. For each curve, we selected the one with the best friendly generator among a set of generated curves.

m = 257 f = x 257 + x 65 + 1 d = t 65 + t 31 + t 14 + 1 G x = 16b46e24aa4b12ab2289fcd3417615387810f0 83f43419d8cae38ad9ac640d960 G y = 16b46e24aa4b12aba289fcd3417615383810f0 83f43419d8eae38ad9ac640d968 G 1/w = 1000000000000000000000000000000000000 000000000000 n = 231584178474632390847141970017375815706 332616967362709021140632923291797618908 p = 578960446186580977117854925043439539265 83154241840677255285158230822949404727 h = 4 m = 313 f = x 313 + x 121 + 1 d = t 38 + t 33 + t 28 + 1 G x = 15c67e3024c7c27466e72a3391256e9a729fc 158092053d89087c0f38408b214b0ade5736 3ea938 G y = 15c67e3024c7c27446e72a3391256e9a529fc 158092053d8b087c0f38408b214b0ade5736 3ea938 G 1/w = 10000000000000001 n = 166873987181321100187111070794496258 953336290808161456226545492179886000 18895406224309766337212 p = 417184967953302750467777676986240647 383340727020403640566363730449715000 4723851556077441584303 h = 4 m = 431 f = x 431 + x 303 + x 239 + x 111 + 1 d = t 83 + t 66 + t 17 + 1 G x = 4e1765c1f2f6140db17d5ef2f14c59a38a93e5 b65ba9acca547bf2cc34f3d55bd85ccf4daeaf 7ca1becaa8ee877b01f8d8acae12b210 G y = 4e1765c1f2f6140d317d5ef2f14c59a30a93e5 b65ba9acca547bf2cc34f3d55b585ccf4daeaf 7ca13ecaa8ee877b01f8d8acae12b210 G 1/
Computations were done on a cluster of 80 cores running Scientific Linux 6 during one and a half month generating 48 new curves of different sizes.

Table 1 :

 1 Comparison (in terms of number of logical operations) between different polynomials as modulus (the ones in bold characters are ours).where d is not a square in F p . In[START_REF] Bernstein | Binary edwards curves[END_REF] [START_REF] Bernstein | Batch binary Edwards[END_REF], Bernstein et al. constructed an equivalent to Edwards curves onto binary fields. Definition 1. Binary Edwards Curves (BEC): Let K be a field of characteristic 2 and let d 1 , d 2 ∈ K with d 1 = 0 and d 2 = d 2 1 + d 1 . The Binary Edwards Curve, with coefficients d 1 and d 2 , is the affine equation:

	Security level	Polynomials	# logical operations # 32-bit registers
	80	x 163 + x 7 + x 6 + x 3 + 1	74	6
	112	x 223 + x 159 + 1	38	7
		x 233 + x 74 + 1	62	8
	128	x 257 + x 65 + 1	44	9
		x 283 + x 12 + x 7 + x 5 + 1	152	9
		x 313 + x 121 + 1	53	10
	192	x 409 + x 87 + 1	108	13
		x 431 + x 303 + x 239 + x 111 + 1	99	14
	224	x 479 + x 255 + 1	78	15
		x 487 + x 295 + x 167 + x 39 + 1	114	16
	256	x 521 + x 489 + 1	79	17
		x 569 + x 441 + x 313 + x 121 + 1	135	18
		x 571 + x 10 + x 5 + x 2 + 1	296	18

  487 + x 295 + x 167 + x 39 + 1 d = t 69 + t 33 + t 15 + 1

			m =	487
			f = x G x = 339b843c53c409543f396d39e57efde813f06
				e3099735004b999b15776a75a4c3a22dcaf1e
				91e261fe479b89a64d65103928195d727bd3
				d157735b2071
			G y =	339b843c53c40954bf396d39e57efde893f0
				6e30997350043999b15776a75a4c3a22dcaf
				1e91e261fe479b89a64d6510b928195d727b
				d3d1d7735b2071
			G 1/w =	10000000000000001
			n =	3995838144404470056168444454135252871
				3582056226111630730997209083204758260
				1818621382266938545215444668253454986
				755337077661569719769645413325977844
			p =	9989595361011175140421111135338132178
				3955140565279076827493022708011895650
				4546553455667346363038611670633637466
				88834269415392429942411353331494461
			h =	4
			m =	521
			f =	x 521 + x 489 + 1
			d =	t 66 + t 29 + t 28 + 1
			G x =	16b369b497b805e6199a342909aa4608cdc
				ecb10e0988ba73eb1f1186039c8b1f6d2a9
				db39b1302d29d9d449b9aa459cc5d6bbb4e
				33a1eb8fcc056ce724cde5aaa8
			G y =	16b369b4b7b805e6199a342909aa4608cdc
				ecb10e0988ba73eb1f1186039c8b1f6d2a9
				db39b1302d29d9d449b9aa459cc5d6bbb4e
				33a1eb8fcc056ce724cde5aaa8
			G 1/w =	100000001
			n =	686479766013060971498190079908139321
				726943530014330540939446345918554318
				339765708943950328983162243962653943
				419778613112216264068985797800513224
				0328602782204
			p =	171619941503265242874547519977034830
				431735882503582635234861586479638579
				584941427235987582245790560990663485
				854944653278054066017246449450128306
				0082150695551
	w =	10000000000000001	h =	4
	n =	554533938824162971915682836828616740	m =	569
		687287415075163315034095916131188222	f =	x 569 + x 441 + x 313 + x 121 + 1
		853620417205051641025857549800325003	d =	t 56 + t 45 + t 41 + 1
		9698819164222488620692	G x =	195b22b2864ee08dd456bab1a95cdd8c7e3
	p =	138633484706040742978920709207154185		fd330fddf630f9c3bb5c33f062b341c919c
		171821853768790828758523979032797055		6bb4cbf1d4335a344ed023b319585ea0e16
		713405104301262910256464387450081250		f03453cc5ba9a86a4b28b16e1c72ad75f11
		9924704791055622155173		41f
	h =	4	G y =	195b22b2864ee08df456bab1a95cdd8c5e3
				fd330fddf630f9c3bb5c33f062b341c919c
				6bb4cbf1d4135a344ed023b319785ea0e16
	m = f = d = G x = G y = G 1/w = n =	479 x 479 + x 255 + 1 t 73 + t 29 + t 3 + 1 7bdd9f19e11e888e80d7c093092d208b4fe99 6e8fcbdffa28cc90173ece2c43673f1372e975 ba9dcd3a06332abf15dbe9b679f6c63e30b88 4ab93272 3bdd9f19a11e888e40d7c093492d208b8fe99 6e8bcbdffa24cc90173ece2c436f3f1372e175 ba9dc53a063322bf15dbe1b679f6ce3e30b88 cab93272 10000000000000001 15608742751579961156907986148965831 52874299071332485575429578479812685	G 1/w = n = p = h =	f03453ce5ba9a86a4b28b16e1c72ad75f11 41f 10000000000000001 193226876150862917234767594546599367 214946366485321749932861762572575957 114478021226803243121352042256909941 058335965554184749327731790508484922 4650340025220880067199309308 483067190377157293086918986366498418 037365916213304374832154406431439892 786195053067008107803380105642274852 645839913885461873319329476271212306 162585006305220016799827327 4
		86941544801971795445881886763046934
		69803241139597788961643097959459945
		58356
	p =	390218568789499028922699653724145788
		218574767833121393857394619953171467
		353862004929488614704716907617336745
		081028489944724041077448986498639589
	h =	4

L q (a) = e c log(q) a log log(q) 1-a

APPENDIX

In this annex we propose a new set of binary elliptic curves with the finite field representation, BEC parameter and the Weierstrass form, a "friendly" (i.e. adapted to 32-bit architectures) generator G and coordinates in Weierstrass form. Each curve was generated using the previously described methods and they are composed of the following parameters:

• m is the size of the extension of F 2 .

• f is the modulus defining the representation of the finite field.

• d is the parameter of the curve in BEC form, equation 1.

• G x , G y and G 1/w are the coordinates of a friendly generator in BEC form with G 1/w the wcoordinate such as G 1/w = 1 w where w = G x + G y .

• n is the number of points of the curve.

• p is the prime order of the sub-group generated by G

• h is the cofactor.

We make the choice of generating a curve over 223 bits which does not reach the 128 bits security level due to the very interesting performance properties of such curves: in a very constrained environment on a very low power devices, a trade off between security and efficiency could be done.


Neila El Heni
email: neila.elheni@telecom-bretagne.eu

Xavier Lagrange

N El

Analytical evaluation of multicast packet delivery and user-clustering schemes in high-speed cellular networks

Keywords: Clustering, Fairness, HSDPA, Multicast transmission, Rayleigh fading

Transmission on data-oriented radio interfaces of cellular networks has been primarily designed for unicast applications. Nevertheless, unicast may not optimize the resource usage when the same content has to be transmitted to several users in the same cell. In this context, multicast seems to be an efficient means to convey data. In this paper, we develop an analytical model that allows the computation of the mean bitrate for both multicast and multiple-unicast transmission schemes. Furthermore, we propose a multicast transmission scheme called the equal-bitrate (EB) algorithm that allocates bandwidth to mobiles according to their instantaneous channel quality. We compare it to adaptations of the well-known max-signal-to-noise ratio and round robin to multicast. We propose to group users into clusters. The clustering method combines multicast and unicast transmission schemes according to the user's average channel conditions. We use the analytical model to evaluate the proposed solutions. We compare the resulting performance against pure multicast and multiple-unicast approaches. We show that the EB algorithm offers a good trade-off between throughput and fairness. Also, we show that mixed clustering achieves good performance compared to conventional clustering methods.

Introduction

Third-generation (3G) data networks have been introduced for the support of multimedia services over wireless links. They have evolved to CDMA2000 1xEvolution Data Optimized [START_REF] Bender | CDMA/HDR: a bandwidth efficient high speed wireless data service for nomadic users[END_REF] and High-Speed Downlink Packet Access (HSDPA) [START_REF] Gpp | High speed downlink packet access (HSDPA): overall description[END_REF][START_REF] Xavier | Principes et évolutions de l'UMTS[END_REF], commonly known as 3.5G systems. The so-called 3.5 systems are enabling the delivery of very demanding applications in terms of network resource. These include video streaming such as Mobile TV [START_REF] Mobitv | Mobile television and radio service provider[END_REF], the VCast service [START_REF] Segan | PC Magazine Verizon's vcast: video over evdo phones[END_REF] from Verizon, multiparty games or group-chat application, etc.

Although 3.5G networks offer high capacity, the expected demand will certainly overcome the available resources. At present, multimedia services are delivered over separate connections on the radio interface for each recipient. Packets are duplicated if there are multiple users who want to receive the same data. In this paper, such an approach is called multiple-unicast. Because spectrum is a limited and expensive resource, the radio interface can support only a handful of highbitrate users simultaneously, and the wireless link may easily become a bottleneck. In this context, multicast seems a promising alternative to conventional multipleunicast. The main benefit of multicast is bandwidth saving, as it prevents clogging up the air interface with redundant transmissions. Also, the overhead on the sender side is reduced since it does not have to maintain addresses of each destination. As commonly accepted, we define a multicast group as a set of terminals listening to the same service.

Multimedia Broadcast/Multicast Service [START_REF] Gpp | Multimedia broadcast/multicast service (MBMS); architecture and functional description[END_REF][START_REF] Gpp | Introduction of the multimedia broadcast multicast service (MBMS) in the radio access network (RAN); stage 2[END_REF] is currently specified in the 3GPP recommendation to address multicast and broadcast services in UMTS/ WCDMA networks. Equivalently, BroadCast and Mul-tiCast Services [8] are defined by 3GPP2 standard to support multimedia content multicast in CDMA2000 networks. However, the focus is on the access and core network rather then on the radio interface. It is implicitly assumed that transmissions are duplicated by the base station (BS) if several user equipments (UEs) in a cell belong to the same multicast group. Multicast over 3G/3.5G radio interfaces is only used with low bitrate channels, which does not scale well when the number of multimedia services increases.

A common point between these networks is that they use a shared channel to serve UEs with data packets on a time-multiplexing basis. In this framework, a radio resource consists in a transmission time interval (TTI), and the scheduler's task is to share time between users of the system.

Multiple-unicast vs multicast issues

During a service session, users may experience different channel conditions from one TTI to another. At each TTI, the packet scheduler uses the reported signalto-noise ratio (SNR) to determine the instantaneous bitrate capacity of each user. Then, depending on its scheduling strategy, it chooses the user to serve with the suitable bitrate. This scheme is illustrated in Fig. 1. Given that the TTI is constant, determining the instantaneous bitrate capacity amounts to determining the packet size, also called the transport block size (TBS), like in HSDPA.

In multiple-unicast, user-specific adaptation of the radio parameters (e.g., to track fast fading) can be applied. As transmissions to users are considered separately, all classical scheduling mechanisms may be used, for instance, Max C/I [START_REF] Borst | User-level performance of channel-aware scheduling algorithms in wireless data networks[END_REF], to maximize the system throughput. However, fairness has to be guaranteed between users who receive the same data. The main drawback with multiple-unicast is that a TTI is restricted to a single user. Alternatively, a TTI may be allocated to all the users listening to the same service.

We call this scheme pure multicast. It is also possible to divide the multicast group into different clusters and to serve them at different times. As such, a cluster includes one or several users who listen to the same service transmission in the same TTI. If N represents the number of users in the multicast group, multiple-unicast is then equivalent to N single-user clusters, whereas pure multicast amounts to a single N-user cluster. A mixed clustering includes, for example, N single-user clusters and a single (N -N )-user cluster. Figure 2 gives the example of two clusters.

Multicast schedulers not only have to select the cluster to serve at each TTI, but they also have the additional difficulty of determining the transmission bitrate that is common to the selected users despite their various channel conditions. Indeed, bitrate selection is a challenging task in multicast scheduling. If the transmission bitrate exceeds the bitrate capacity of a UE, this mobile will not be able to decode the transmitted packet, which is consequently lost. This is particularly true for mobiles at the cell edge. While the usage of retransmission algorithms [START_REF] Vartiainen | Point-to-multipoint multimedia broadcast multicast service (MBMS) performance over HSDPA[END_REF] is possible to counter losses, subsequent bandwidth consumption and signaling may be costly. In all cases, the useful multicast bitrate is constrained by the UE having the worst channel, as this determines the service coverage. Hence, the gain of multicast compared to multiple-unicast is not so obvious and needs to be studied. Indeed, it is difficult to predict whether it is better to duplicate data with a relatively high bitrate or to send it only once but with the lowest bitrate. Very few studies address the multicast schedulers' design.

The study [START_REF] Won | Multicast scheduling in cellular data networks[END_REF] proposes multicast proportional fairness and intergroup proportional fairness schedulers that are adaptations of the PF scheduler. However, packet loss with these algorithms may be frequent, wasting bandwidth because of retransmissions.

Paper scope and outline

In this study, we propose a new algorithm for multicast transmission over high-speed cellular systems and define a user clustering strategy that combines pure multicast and multiple-unicast schemes. While existing work on multicast generally addresses streaming services, our study is agnostic of the service type, which could be file transfer, video streaming, etc. We consider the transmission of one service. Scheduling multiple services amounts to managing priority between these services according to their quality-of-service (QoS) requirements. These issues have already been extensively developed in literature, e.g., [START_REF] Lundevall | Streaming applications over HSDPA in mixed service scenarios[END_REF][START_REF] Kazmi | Scheduling algorithms for HS-DSCH in a WCDMA mixed traffic scenario[END_REF][START_REF] Braga | Packet scheduling for VoIP over HSDPA in mixed traffic scenarios[END_REF]. This paper is organized as follows: In Section 2, the system model and assumptions are given. Section 3 explains the proposed equal-bitrate (EB) algorithm and presents the model for throughput calculations of max-SNR and round robin (RR). Section 4 places our model in the framework of a theoretical generic system. Particular case studies are developed in Section 5. In Section 6, we introduce state-of-the-art clustering schemes before defining a new clustering strategy. Section 7 analyzes the performance of the proposed solutions. Conclusions are drawn in Section 8.

Model description

General assumptions

We consider a system with one serving cell and N users in the same multicast group. Dealing with multicast in different cells is equivalent to studying multicast in fixed networks, as data have to be transmitted over wired links to different BSs. This has been developed in [START_REF] Rummler | Modeling and analysis of an efficient multicast mechanism for UMTS[END_REF]. Here, we consider downlink transmission like in a HSDPA system with only one code. Large-scale mobility aspects are not considered; signaling for clustering dynamics is, then, not developed here but may be considered as a future work.

Users are randomly distributed in the cell. To apply our multicast transmission algorithm, clusters should already be defined according to the proposed clustering strategy. Let G be the number of clusters and S j be the size of cluster j. Let γ i, j be the instantaneous SNR of user device i, which is a member of cluster j, and γ i, j , its average SNR. Note that indices i or j may be sometimes omitted for simplicity. Due to channel variations, γ i, j are identically and independently distributed variables that change randomly from one TTI to another but are assumed to be constant during a TTI. Without loss of generality, we assume that the user having the highest average SNR is indexed by 1 and belongs to cluster 1. Let k i, j (i = 1..S j , j = 1..G) be the SNR ratios such that k i, j = γ 1,1 γ i, j . By definition, k 1,1 = 1 and k i, j ≥ 1 for i and j higher than 1. Assuming the BS is serving user i of cluster j, we define β i, j as the largest TBS supported by that user. Let g be the function that relates β i, j to the reported γ i, j of the served user; hence, β i, j = g(γ i, j ). It is easy to see that g is a strictly increasing function. Its associated inverse function is denoted as h, i.e., γ i, j = h(β i, j ).

The main notations used in this paper are listed in Table 1. We also define P X (x) as the cumulative distribution function (CDF) of a random variable X. Hence, P γ i, j (x) = Pr(γ i, j ≤ x). Similarly, p X (x) denotes the probability distribution function (PDF) of X. Largest TBS supported by UE i , which is member of cluster j β j TBS of cluster j γ i, j

Instantaneous SNR of user i, which is member of cluster j γ i, j

Average SNR of user i, which is member of cluster j γ j SNR within cluster j

Propagation model

The average SNR may be computed by using a conventional propagation model. Let P i be the transmit power to user i. The received power denoted as P r is, then,

P r = P i h i χ i , (1)
where h i is the path gain including shadowing and distance loss for user i and χ i is the fast fading between user i and the BS. Variable χ i is a random variable that represents Rayleigh fast fading. It therefore has an exponential distribution. The path gain for user i at distance r i from the BS is described by

h i = 10 A 0 10 r β i ( 2
)
where β is the pathloss exponent and A 0 is the distanceloss at 1 m (with a BS antenna height of 30 m, a UE antenna height of 1.5 m, and a carrier frequency of 1,950 MHz).

The signal-to-interference ratio received by user i is [START_REF] Heni | Multicast vs multiple-unicast scheduling in high-speed cellular networks[END_REF]

γ i = P i h i χ i α(P max -P i )h i χ i + I ext (3)
where P max is the total transmit power of the cell, α is the orthogonality factor, and I ext represents the intercell interference. Note that I ext is assimilated to a constant similarly to the model given by the "enhanced UMTS radio access network extension for ns-2" (Eurane) simulator [START_REF] Eurane | Eurane website[END_REF]. By using further modifications, we can obtain the average SNR of user i

γ i = γ sup 0 exp -I ext 10 A 0 10 r β i x P i -α(P max -P i )x dx, (4)
where γ sup is defined by

γ sup = P i α(P max -P i ) . ( 5
)

Multicast scheduling

Multicast scheduling is performed in two steps. First, the scheduler determines the convenient transmission bitrate for each cluster. The intracluster bitrate allocation strategy is conservative. We then have

γ j = min i=1..S j γ i, j (t). ( 6
)
The CDF of γ j is equal to

P γ j (x) = 1 - S j i=1 1 -P γ i, j (x) . ( 7
)
Once the bitrate of each cluster is determined, the scheduler chooses the cluster to serve. In order to maximize the global throughput, a natural solution is to use max-SNR scheduling that privileges the cluster having the highest bitrate capacity. However, the scheduler must guarantee fairness between clusters. RR scheduler is fair with respect to system resources (i.e., time slots), but in general, the bandwidth is not shared equally among the UEs because they have different bitrate capacities. The weighted-fair queuing RR (WFQ-RR) is another variant of RR; it uses weights to achieve the same throughput for all clusters; however, WFQ-RR is agnostic of the channel conditions. In this context, we introduce the EB scheduler that realizes the same average bitrate for all the clusters while using link adaptation to address channel variations. In the next sections, we introduce our EB algorithm, then we present briefly the computation models for max-SNR and WFQ-RR.

EB scheduler

As all the clusters belong to the same multicast group, and in order to achieve a maximum fairness between them, EB targets the same average bitrate for all the clusters. For this purpose, we define fairness factors {F j } j=1..G such that a cluster having a higher instantaneous bitrate capacity is served with a lower probability, i.e., the time resource is not uniformly shared between clusters. For clusters with low average channel quality, the bitrate capacity is low. Hence, in order to increase the average throughput for these clusters, their corresponding fairness factors are relatively high such that they are served more frequently. Note that, for clusters having the same average channel quality, their fairness factors are equal. At every TTI, cluster j is served if the product of its instantaneous SNR and its corresponding fairness factor F j is the highest among all the clusters. Hence, cluster j is served if and only if

γ j F j = max l=1..G (γ l F l ). ( 8
)
If condition 8 is fulfilled, the TBS of cluster j, denoted as β j , is then equal to g(γ j ). Otherwise, β j = 0. Let β j (γ j ) be the mean TBS of cluster j given that γ j is known. We have

β j (γ j ) = g(γ j ) G l=1,l = j P γ min,l γ j F j F l . ( 9
)
Note that if there is only one cluster, the mean TBS is equal to g(γ j ). Equation 9 is still valid if the product of a null number of factors is considered equal to 1. The mean bitrate of cluster j denoted as R j is obtained by considering all the γ j values

R j = 1 D TTI ∞ 0 p γ j (x)β j (x)dx, ( 10
)
where D TTI is the TTI duration. Substituting Eq. 9 in Eq. 10, we obtain

R j = 1 D TTI ∞ 0 p γ j (x)g(x) G l=1,l = j P γ min,l xF j F l dx. ( 11
)
Equation 11 provides both a general formula for the average bitrate per cluster j and a method to determine the values of {F j } j=2..G . It can be seen from Eq. 8 that only the ratio F l F j impacts the algorithm. Hence, F 1 may be arbitrarily set to 1. As we set the average bitrate given in Eq. 11 to the same value for all the clusters, we then have a set of G -

1 equations R j = R 1 ( j = 2..G) with G -1 unknown values, namely, F j ( j = 2..G).
In an operational system, the scheduler may start by allocating the same fairness factor for all clusters and then increases this parameter gradually in clusters for which the buffer at the BS is close to saturation. A steady state is reached when all buffers have equivalent loads. Alternately, it is possible to use a lookup table at the BS to map the average SNRs to their corresponding fairness factors. Of course, that supposes a prior knowledge of the average channel conditions for each user. In the case of HSDPA, measurement reports and uplink feedbacks of the channel quality indicator may be used by the BS to estimate the average SNR.

Max-SNR scheduler

In order to evaluate our proposed EB scheduler, it would be interesting to compare it with reference schedulers like max-SNR or RR. Max-SNR is the most greedy scheduler with respect to throughput. At each TTI, it selects the cluster with the maximum SNR. Hence, throughput formulation can be deduced from EB scheduler for {F j } j=1..G = 1.

R j = 1 D TTI ∞ 0 p γ j (x)g(x) N l=1,l = j P γ l (x)dx. ( 12
)
3.3 WFQ-RR WFQ-RR scheduler serves users in a RR scheme without caring about their radio channel variations. It tar-gets the same throughput for all clusters. Let R j be the average bitrate capacity of cluster j. In other words, R j represents the throughput of cluster j if it were served all the time. Hence,

R j = 1 D TTI E[g(γ j )], (13)
which is developed as follows:

R j = 1 D TTI ∞ 0 1 -P γ j (h(x)) dx. ( 14
)
Combining Eqs. 7 and 14, we deduce the average bitrate capacity of cluster j

R j = 1 D TTI ∞ 0 S j i=1 1 -P γ i, j (h(x)) dx. ( 15
)
In order to guarantee fairness among clusters with different channel conditions, the scheduler allocates to each cluster j a bandwidth fraction ρ j such that

⎧ ⎪ ⎪ ⎨ ⎪ ⎪ ⎩ ρ i R i = ρ j R j ∀i, j = 1..G G j=1 ρ j = 1 (16)
After a few elementary computations, the resolution of Eq. 16 gives that

ρ j = G i = j R i G k=1 ⎛ ⎝ G i =k R i ⎞ ⎠ . ( 17
)
Then, the effective average bitrate per cluster is given by

R j = ρ j R j . ( 18
)
4 Computation for a generic system

Generic system considerations

We consider a generic system based on the Shannon formula [START_REF] Shannon | A mathematical theory of communication[END_REF] that gives the maximum reachable system capacity. Shannon assumes a perfect error correcting system. It allows the computation of the maximum bitrate depending on the SNRs. If W is the available bandwidth, the maximum bitrate of user i is then

R max (bps) = W log 2 (1 + γ i ). ( 19
)
The TBS is then derived as follows:

g(γ i ) = W D TTI log 2 (1 + γ i ). ( 20
)
Next, we can deduce the associated inverse function

h(x) = 2 x W D TTI -1. (21)
In [START_REF] Knopp | Information capacity and power control in single-cell multiuser communications[END_REF], Knopp and Humblet proposed a reference radio channel model based on an exponential distribution for γ .

p γ (x) = 1 γ . exp -x γ if x > 0 0 otherwise, ( 22
)
where γ is the average received SNR. The associated CDF is

P γ (x) = 1 -exp -x γ . ( 23
)
Note that the SNR given in Eq. 3 does not follow an exponential distribution. However, we can verify that, when γ sup tends to infinity, the SNR has an exponential distribution. Hence, the approximation of using both Eq. 23 for the CDF of γ and Eq. 4 for its mean value makes sense.

Assuming that the SNR of each user follows an exponential distribution, we can deduce the CDF of γ j from Eqs. 7 and 23 as follows:

P γ j (x) = 1 -exp ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ -x S j i=1 k i, j γ 1,1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ . ( 24
)
Its corresponding PDF is derived as follows:

p γ j (x) = S j i=1 k i, j γ 1,1 exp ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ -x S j i=1 k i, j γ 1,1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ . ( 25
)
4.2 EB throughput Substituting Eqs. 20, 24, and 25 in Eq. 11, we reformulate the EB average bitrate per cluster j as follows:

R j = W S j i=1 k i, j γ 1,1 ∞ 0 ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ exp ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ -x S j i=1 k i, j γ 1,1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ log 2 (1 + x) × G l=1,l = j ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 -exp ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ - xF j F l S l i=1 k i,l γ 1,1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ dx. ( 26
)
Equation 26 gives the average bitrate per cluster j considering Shannon and an exponential distribution for SNR.

As multiple-unicast case corresponds to N singleuser clusters, the average bitrate is derived from Eq. 26 for G = N and S j = 1 (∀ j = 1..G). We denote it as R ucast and we formulate it hereafter

R ucast = W ∞ 0 1 γ j . exp -x γ j log 2 (1 + x) × N l=1,l = j 1 -exp - xF j F l γ j ⎤ ⎦ dx (27) with F j verifying R j = R 1 (∀ j = 1..N).
In the framework of pure multicast, We denote the average bitrate as R mcast . It is derived after a few computations including integration by parts and variable changing.

R mcast = 1 D TTI . ∞ 0 exp -2 x W.D TTI -1 N i=1 k i γ 1 dx. ( 28
)
Note that cluster indice j is omitted from k and γ parameters as there is only one cluster. Finally, formulations for max-SNR are deduced from EB expressions with fairness factors set to 1.

WFQ-RR throughput

Considering the (K&H) model given in Eq. 23 and the inverse of shannon formula given in Eq. 21, we reformulate Eq. 15 as follows:

R j = 1 D TTI ∞ 0 S j i=1 exp -2 x W D TTI + 1 γ i, j dx, ( 29
)
which is equivalent to

R j = 1 D TTI ∞ 0 exp ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ -2 x W D TTI + 1 S j i=1 k i, j γ 1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ dx (30)
From Eqs. 18 and 29, the average bitrate per cluster is

R j = ρ j D TTI ∞ 0 exp ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ -2 x W D TTI + 1 S j i=1 k i, j γ 1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ dx (31)

Case studies

In this section, we target performance evaluation for some particular cases like equivalent users or uniform user distribution. For the sake of simplicity and due to space limitation, we restrict ourselves to the EB algorithm, first because it makes part of the contributions of this work, and second because it seems to offer a good trade-off between system throughput and user fairness. In fact, unlike WFQ-RR, this scheduler uses link adaptation to track the channel variations, and, unlike max-SNR, it realizes throughput fairness between clusters. The objective is to study, in simple examples, the behavior of our scheduler. Also, this part orients us toward a new efficient clustering scheme. The comparison with the other schedulers is addressed in a more general case in Section 6.

Serving equivalent users

Considering N users with the same average SNR, we compare the single N-user and the N single-user clustering schemes. We have compared these schemes in [START_REF] Heni | Multicast vs multiple-unicast scheduling in high-speed cellular networks[END_REF]; however, we take back the obtained results to justify our new clustering strategy. Let γ be the common average SNR value for all users. In this context, fairness factors F j and SNR ratios k j , ( j=1..N) are equal to 1. Hence, the average unicast bitrate is deduced from Eq. 27 as follows:

R ucast = W γ ∞ 0 exp -x γ log 2 (1 + x) × 1 -exp -x γ N-1 dx (32)
and Eq. 28 becomes

R mcast = 1 D TTI . ∞ 0 exp -(2 x W.D TTI -1)N γ dx. ( 33
)
In the framework of multiple-unicast, our EB algorithm is equivalent to the opportunistic max C/I scheduler. As users have the same average SNR, they will receive the same bandwidth fraction on average and fairness is not affected. We evaluate the multicast gain, which represents the ratio between the pure multicast bitrate and the multiple-unicast bitrate. We denote this gain as Γ m . Hence,

Γ m = R mcast R ucast . ( 34
)
Simulation parameters are listed in Table 2; they are compatible with Eurane simulator [START_REF] Eurane | Eurane website[END_REF] and coherent with our propagation model to be applied in typical urban environments. Two regions can be depicted. Figure 3 presents the multicast gain Γ m vs the average SNR for various N. For high average SNRs (above around 3.7 dB), it is better to use multicast, for which the gain increases with N. In this case, bandwidth saving by transmitting packets only one time to several users has a higher impact than the bitrate limitation by the minimum SNR. It is also observed that the larger the average SNR becomes, the greater the span in system gain becomes between different values of N. However, when average SNRs are low (below 3.7 dB), multiple-unicast becomes more attractive. In this case, it is better to take advantage from opportunistic scheduling and to select the UE that has the best SNR even if the same data are duplicated over several times to the different users. Moreover, the greater N is, the smaller the multicast gain is. Performance evaluation for HSDPA is coherent with the generic case. Multicast can largely outperform multiple-unicast for acceptable average channel quality. When the mean SNR is high enough, minimum and maximum bitrate capacities within a multicast group are equivalent, and then, as multicast serves all terminals at one TTI instead of N TTIs, the multicast gain limit is

lim γ →∞ Γ m = N. ( 35
)
This gain value is also obtained for multicast in fixed networks.

Study of a uniform user distribution

In this section, we consider five users (U E i=1..5 ) that are uniformly distributed over a cell and sorted by their distance from the BS in an increasing order of their indices. As such, U E 1 , . . . , U E 5 are distant from the BS by 223, 316, 387, 447, and 500 m, respectively. Their average SNRs are deduced from their respective distances from the BS as shown in Eq. 4. Other simulation parameters are listed in Table 2. The objective is to deduce in a simple example the best clustering strategy to be coupled with the EB scheduler, and also to verify coherence with results of Section 5.1. Table 3 gives the most interesting tested configurations considering different clustering schemes. Note that F j values are computed by numerical resolution of an equation system, as has been explained at the end of Section 3. Configuration 1 corresponds to pure multicast where all users belong to cluster C 1 . Configurations 2, 3, and 4 consider two clusters. For instance, in configuration 2, users with high γ values and indexed by 1, 2, and 3 are grouped in the same cluster (C 1 ) and users 4 and 5 in a second cluster (C 2 ). Configuration 5 considers three clusters. Finally, configuration 6 illustrates the multipleunicast case with five clusters.

We can verify in configurations 1 and 2 that grouping users with low channel quality (U E 4 and U E 5 ) offers the worst performance (1.69 and 1.87 Mbps). Multipleunicast is then better for these users. On the other hand, configuration 5, which isolates users with low average SNRs while grouping together terminals with good channel quality, offers the best throughput (3.11 Mbps). This confirms the idea of using multicast for users near the BS. It is noteworthy that, instead of serving users having low average SNRs in a unicast mode, it seems interesting to cluster them with users close to the BS. In fact, as the bitrate allocation strategy to each cluster is conservative, the bitrate value is influenced by the lowest capacity within the cluster. Hence, adding a good user will have no tangible impact on the throughput. This can be verified when we look at configurations 3 and 4 in Table 3. In the former, the worst user is isolated in a separate group; in the latter, the user with the highest SNR is grouped with the worst user. The resulting performance is the same.

User clustering

Multiple-unicast and pure multicast are two extreme clustering strategies. In a previous work [START_REF] Heni | Multicast vs multiple-unicast scheduling in high-speed cellular networks[END_REF], we compared these strategies for the particular case of equivalent users, i.e., users with the same average SNR. In this study, we consider a mixed clustering strategy in a more general case where users are anywhere in the cell.

(γ 1 = 24, γ 2 = 9, γ 3 = 5, γ 4 = 0.8 and γ 5 = -3 dB) Configuration Results ID G {U E 1 ,U E 2 ,U E 3 ,U E 4 ,U E 5 } {F j } j=1..G Rate (Mbps) 1 1 C 1 , C 1 , C 1 C 1 , C 1 {1} 1.69 2 2 C 1 , C 1 , C 1 , C 2 , C 2 {1, 28.4} 1.87 3 2 C 1 , C 1 , C 1 , C 1 , C 2 {1, 2.3} 2.06 4 2 C 1 , C 2 , C 2 , C 2 , C 1 {1, 1.92} 2.06 5 3 C 1 , C 1 , C 1 , C 2 , C 3 {1, 2.1, 11} 3.11 6 5 C 1 , C 2 , C 3 , C 4 , C 5 {1, 59
, 187, 680, 2980} 1.96

Previous work on multicast clustering

Existing multicast schemes introduce two clustering approaches. First, we find basic clustering that ignores QoS aspects. It includes service-oriented clustering [START_REF] Eusebio | Management scenarios for multicast groups in enhanced-UMTS[END_REF], where clusters are created depending on the service availability. This means that a user subscribing to a specific multicast service is integrated into the cluster where this service is provided, independently of the user's location. If only one service is considered, this scheme is equivalent to pure multicast with one cluster. Another variant of the service-oriented clustering is given in [START_REF] Bakker | A method of providing a multimedia service in a digital wireless communication network[END_REF]. Instead of automatically using multicast for users subscribing to the same service, authors propose to first check the number of these users and then decide the transmission mode. If this number is below a predefined threshold, multiple-unicast is adopted. Otherwise, pure multicast is used. Considering only one service, we clearly see that this clustering scheme is equivalent to either multiple-unicast or pure multicast.

As basic clustering schemes do not consider the user channel quality, the deviation between the best and the worst user within a cluster may be so important that the allocation of a common bitrate inside a cluster will lack accuracy. In fact, if the chosen bitrate is high, there will be losses for users having bad radio conditions. Conversely, if the bitrate is too low, terminals with good radio quality will be penalized as their channels will be under-exploited. A second clustering approach is proposed in [START_REF] Eusebio | Management scenarios for multicast groups in enhanced-UMTS[END_REF] to alleviate the aforementioned problems. It involves a QoS-oriented approach where users are grouped according to their received quality. Conceptually, UEs near the BS should belong to a cluster with better quality and faraway UEs belong to other clusters with lower quality. The division of UEs into quality groups ensures that the deviation between the best and the worst reception within a cluster is less significant. However, according to the previous section, grouping users with low SNRs decreases the system performance. For this reason, we propose a new QoSbased clustering scheme, namely, mixed clustering.

Proposed clustering strategy

We have seen in Section 5 that multicast is attractive only for high average SNRs. Our proposed clustering scheme is then deduced as illustrated in Fig. 4.

-An average SNR threshold is fixed so that the system can differentiate users. The average SNR is declared as low if it is below a threshold value denoted as γ thres . Let N low be the number of users having a low average SNR. As it combines pure multicast and multiple-unicast, our clustering scheme is called mixed clustering.

Evaluation results

In order to evaluate our mixed clustering approach, we compare it to multiple-unicast with N single-user clusters and pure multicast with a single N-user cluster. We take γ thres equal to 3.7 dB, which is the threshold value depicted in Fig. 3. Section 7.3 shows that the tolerance margin for γ thres is quite large. We consider randomly located users. Several iterations are performed until acceptable confidence intervals are achieved. Each iteration corresponds to a given user distribution. Then, the average performance is computed over all the iterations including the 95% confidence intervals.

Performance evaluation in terms of average bitrate

We compare the throughput performance of the studied schedulers associated to different clustering schemes, namely, multiple-unicast, multicast, and mixed. For instance, considering the EB scheduler, we record in each iteration the average bitrate given by Eq. 26 for mixed clustering, by Eq. 27 for multipleunicast, and by Eq. 28 for pure multicast. The throughput performance averaged over all the iterations for EB scheduler in comparison with max-SNR and WFQ-RR schedulers is shown in Fig. 5. We can state that, for all the schedulers, mixed clustering offers the highest performance. We can also verify that curves of pure multicast are the same for all schedulers. In fact, as pure multicast throughput deals with a single cluster, it is independent of the cluster selection scheme. Instead, it depends on the link adaptation scheme, which is conservative for all our schedulers. Dealing with the other clustering schemes, we see that max-SNR offers the highest throughputs, while the lowest performance is obtained with WFQ-RR. In fact, the latter does not take profit from link adaptation.

Figure 6 shows the ratio between the EB throughput with mixed clustering and EB throughput with conventional clustering, namely, pure multicast and multipleunicast. This gain is computed for a fixed number of terminals and different user distributions. In each distribution, the average distance between UEs and the BS is represented on the x axes. When all the users have high SNRs, our clustering solution is equivalent to pure multicast, the gain against pure multicast being equal to 1. On the other hand, the case when all the users have low average SNRs comes down to multipleunicast as the gain against multiple-unicast is equal to 1. In the other cases where users have heterogeneous channel qualities, the EB algorithm combined with the new clustering scheme generally offers the best performance, as the gain values in Fig. 6 are always higher than 1.

Figure 7 shows the number of clusters vs the average distance between users and the BS. For low average Fig. 7 Number of clusters vs the average distance from the BS distances from the BS, the number of clusters is equal to 1; hence, we verify that our strategy corresponds to pure multicast. When users get closer to the cell border, the number of clusters increases and our strategy tends to multiple-unicast.

Fairness evaluation

Fairness is an important issue for the scheduling design. A fair scheduler prevents the automatic resource allocation to users with a good channel quality. In order to evaluate fairness for each scheduler, we use the coefficient of variation (CV), which is a normalized measure of the probability distribution dispersion. It is defined as the ratio of the standard deviation to the mean.

As EB and WFQ target the same throughput for all clusters, their CV is equal to zero for all the clustering schemes. Conversely, max-SNR scheduler results in high CV. In Fig. 8, we estimate the mean CV values for the different schedulers with multiple-unicast and mixed clustering. Of course, the CV is null for pure multicast. The CV is lower with mixed clustering compared to multiple-unicast, and it increases for a higher number of users. We can clearly see that, in spite of the good performance of max-SNR in terms of throughput and transmission duration, it is unfair to users with a bad channel quality. The CV reaches 28% for 15 users served in a unicast scheme.

Impact of the average SNR threshold γ thres

In order to see the impact of γ thres on the new clustering scheme, we perform additional simulations for different values of this parameter. The gain is averaged over 100 iterations. Let Q ucast be the gain of the new clustering method compared to multiple-unicast. Hence,

Q ucast = R mixed -R ucast R ucast , ( 36
)
where R mixed and R ucast are the average bitrate performances of the new clustering strategy and multiple- unicast, respectively. Similarly, we define Q mcast as the gain compared to pure multicast. It is noteworthy that the performance of multiple-unicast and pure multicast is independent of γ thres . Only the bitrate performance of the new clustering scheme is influenced by γ thres values.

Results are given in Table 4. The proposed clustering method provides the highest gain for γ thres around 3.7 dB, which is coherent with the value obtained in Fig. 3. For low γ thres values (-5 dB), we see that the new clustering scheme is similar to pure multicast as Q mcast is equal to 0%. In fact, in this case, all users are considered with high average channel quality and are therefore clustered together. For high γ thres values (15 dB), the new clustering strategy is similar to multiple-unicast, as Q ucast is equal to 0%. In fact, in this case, all users tend to be considered with low average SNRs and are therefore separated into different clusters. According to Fig. 5, the global performance of EB with multiple-unicast is lower than the one with pure multicast. Hence, Q mcast is negative when γ thres is very high.

Except for the highest value of γ thres , the proposed clustering method always provides a gain compared to both multicast and multiple-unicast, whatever γ thres is. Among moderate values of γ thres , the threshold of 3.7 dB provides the highest gain. It is noteworthy that the gain with neighboring threshold values (like 2 and 5 dB) remains interesting. This allows an important tolerance margin for the choice of γ thres in which the use of our clustering scheme is favorable.

Conclusion

This study deals with multicast transmission over highspeed cellular networks. In this framework, we propose the EB algorithm, which allows fairness among clusters while offering good system throughput. We have developed an analytical model for the average bitrate calculations in order to evaluate the resulting scheduling performance and compared it with two reference schedulers, namely, max-SNR and WFQ-RR. To ensure an efficient usage of our algorithm, we have proposed to use a mixed clustering strategy. We have shown that our solution outperforms both multipleunicast and pure multicast strategies. Although the gain is influenced by the threshold value for the average SNR, a wide tolerance margin for this parameter exists in which it is always favorable to use the new clustering scheme.

The proposed solution is agnostic of the service type and may be refined by considering both radio metrics and QoS metrics. However, it is quite complicated to model these QoS metrics, and resorting to simulations seems to be necessary. This could be considered in a future work.

Fig. 1

 1 Fig. 1 Conventional unicast scheduling

Fig. 2

 2 Fig. 2 Multicast transmission with an example of clustering

Fig. 3

 3 Fig. 3 Multicast gain Γ m vs γ for different N

Fig. 4

 4 Fig. 4 Proposed clustering strategy

Fig. 5 Fig. 6

 56 Fig. 5 Throughput (Mb/s) for the different scheduling and clustering strategies with 95% confidence intervals

Fig

  Fig. 8 CV for the different schedulers and clustering schemes

Table 1

 1 List of the main notations

	Notation	Meaning
	g	Function that associates a TBS
		to the reported SNR
	G	Number of clusters
	i, j	Indices of a UE and a cluster, respectively
	k i, j	SNR ratios such that k i, j =	γ 1,1 γ i, j
	F j	Fairness factor of cluster j
	N	Number of users in the cell (all belonging
		to the same multicast group)
	Q mcast	Gain of the new clustering method
		compared to pure multicast
	Q ucast	Gain of the new clustering method
		compared to multiple-unicast
	R j	Mean bitrate of cluster j
	S j	Size of cluster j
	β i, j

Table 2

 2 Simulation parameters

	Frame period (D TTI )	2 (ms)
	BS Transmission power (P max )	3 8 ( d B m )
	Power of Intra-cell interferants (P max -P i )	3 0 ( d B m )
	Inter-cell interference (I ext )	-100 (dBm)
	β	3.52
	A 0	31.8 (dB)
	W	5 (MHz)

Table 3

 3 Clustering schemes' performance for five uniformly distributed users

Table 4

 4 Gain of the new strategy compared to pure multicast and multiple-unicast for various γ thres values

	γ thres (dB)		-5	2	3.7	5	15
	Q mcast	5 UEs	0%	13%	19%	4%	-18%
	10 UEs	0%	18%	23%	10%	-26%
	15 UEs	0%	23%	30%	14%	-35%
	Q ucast	5 UEs	23%	40%	46%	28%	0%
	10 UEs	36%	56%	68%	34%	0%
	15 UEs	54%	71%	100%	46%	0%