Marin Bougeret
email: marin.bougeret@lirmm.fr

György Dósa

Noam Goldberg
email: noam.goldberg@biu.ac.il

Michael Poss
email: michael.poss@lirmm.fr

G Y Örgy

D Ósa

Constant-Ratio Approximation for Robust

Keywords: Bin-packing, robust optimization, approximation algorithms, Next-fit-decreasing, dynamic programming

1. Introduction. This paper studies approximation algorithms for a generalization of the bin packing problem that requires its solutions to be robust with respect to uncertain input data. In general, approximation algorithm bound results may involve asymptotic or absolute approximation ratio guarantees whose definition can be stated as follows: An algorithm has an asymptotic approximation ratio ρ for a minimization problem if and only if there exists a constant c such that it outputs a solution of value that is at most ρOPT(I) + c, for every instance input I. The definition of an absolute approximation ratio ρ is similar with the same upper bound only that c = 0.

Bin packing is the problem of assigning a given set of n items, each item of a specified size, to the smallest number of unit capacity bins. The problem has been the subject of study in an extensive body of research initiated by several publications in the 1970s including the work of Johnson et al. [START_REF] Johnson | Worstcase performance bounds for simple one-dimensional packing algorithms[END_REF]. The problem is N P-hard and in fact a straightforward reduction from the partition decision problem implies that it is N P-hard to determine whether a bin packing instance has a solution using only two bins. This also shows that the problem cannot be approximated within a factor less than 3/2 unless P = N P. An absolute approximation factor guarantee of 3/2 has been proven for the first-fit decreasing algorithm [START_REF] Simchi-Levi | New worst-case results for the bin-packing problem[END_REF]. Several results concern asymptotic setting of classical bin packing. It is known since the problem has been first analyzed that the first-fit decreasing algorithm has an asymptotic approximation ratio of 11 9 [START_REF] Johnson | Worstcase performance bounds for simple one-dimensional packing algorithms[END_REF]. This result has also been extended to analyze the constant for first-fit decreasing, in particular it turns out that 11 9 OPT(I) + 6 9 is a tight upper bound (so c = 6 9 ) [START_REF] Dósa | The tight bound of first fit decreasing bin-packing algorithm is FFD(I) <= 11/9 OPT(I) + 6/9[END_REF]. A well known result is the asymptotic fully polynomial-time approximation scheme for bin packing established by [START_REF] Karmarkar | An efficient approximation scheme for the one-dimensional bin-packing problem[END_REF], following the approximation scheme first proposed by [START_REF] Vega | Bin packing can be solved within 1+epsilon in linear time[END_REF].

Another distinct line of studies, which may not include the focus of the current paper but has received significant attention in the literature, is that of bin packing online algorithms, a setting in which the instance is not given in advance but items are revealed and packed one at a time as they arrive (note that first-fit decreasing is not an online algorithm). In this context the performance measure, which bounds the ratio of the algorithm's solution objective value to the optimal offline solution objective value, is called the competitive ratio. The best asymptotic and absolute online competitive ratios of 1.578 and 5/3, respectively, have been shown in [START_REF] Balogh | A New and Improved Algorithm for Online Bin Packing[END_REF] and [START_REF] Balogh | The optimal absolute ratio for online bin packing[END_REF], respectively, while the best lower bound for the asymptotic competitive ratio of around 1.54278 has been provided in [START_REF] Balogh | A new lower bound for classic online bin packing[END_REF]. Approximation results in the offline and online settings have also been developed for different extensions and generalizations of the classical bin packing problem; for example see [START_REF] Christensen | Approximation and online algorithms for multidimensional bin packing: A survey[END_REF][START_REF] Goldberg | Online packing of arbitrary sized items into designated and multipurpose bins[END_REF] and references therein.

In many applications, the sizes of the items to be packed are not fully known at the time that the packing is carried out. In cargo shipping, for example, the actual weight of a container may deviate from its declared weight or its measurements may be inaccurate. Bin packing has been used to model the assignment of elective surgeries to operating room in hospitals [START_REF] Dexter | Which algorithm for scheduling add-on elective cases maximizes operating room utilization? use of bin packing algorithms and fuzzy constraints in operating room management[END_REF]. Here a bin is a shift of a properly equipped and staffed operating room for performing a certain type of elective surgeries. The room scheduler has to fit in the bins as many cases (patients) as possible. In this setting clearly the length of time of performing each surgery is subject to uncertainty for example in the event of complications. Bin packing variants has also been used to model the scheduling of internet advertising; see for example [START_REF] Adler | Scheduling space-sharing for internet advertising[END_REF]. When multiple commercials are packed into advertising breaks of online video services and each commercial may be skipped by the viewer then the length of each ad is subject to uncertainty.

One way to model the uncertainty that falls into the framework of robust optimization is to assume that the sizes are uncertain parameters taking any value in a given set U ⊂ R n , where each a ∈ U represents a possible scenario. Modern robust optimization considers such sets U in continuous as well as combinatorial optimization, and most common are discrete, polyhedral, and ellipsoid uncertainty sets; see [START_REF] Ben-Tal | Robust Optimization[END_REF] for a comprehensive reference on continuous robust optimization and [START_REF] Buchheim | Robust combinatorial optimization under convex and discrete cost uncertainty[END_REF] as a survey of robust combinatorial optimization. Incorporating U in bin packing leads to the following problem (note that the description of U may not be explicit to avoid an exponential length in n).

Robust bin packing (RBP)

Input: U ⊂ R n Output: A solution that is a partition of [n] into k bins b 1 , . . . , b k such that max a∈U i∈bj a i ≤ 1 for each j ∈ [k] Minimize: k
In this paper we focus on two specific uncertainty sets, defined by simple budget constraints. One of these widely used uncertainty sets [START_REF] Bertsimas | Robust discrete optimization and network flows[END_REF], U Γ , supposes that the size of each item is either its given nominal size āi , or its peak value āi + âi . Furthermore, in any scenario, at most Γ ∈ N of the items may take their peak value simultaneously.

Formally, U Γ can be defined as

U Γ = {a|∀i ∈ [n], a i ∈ {ā i , āi + âi } and i∈[n]
(a i -āi )/â i ≤ Γ}. 1 We also consider a second uncertainty set (used in [START_REF] Gounaris | The robust capacitated vehicle routing problem under demand uncertainty[END_REF][START_REF] Pessoa | Branch-and-cut-and-price for the robust capacitated vehicle routing problem with knapsack uncertainty[END_REF][START_REF] Tadayon | Algorithms and complexity analysis for robust single-machine scheduling problems[END_REF], among others), characterized again by ā and â, as well as the number Ω ∈ [0, 1] stating how much deviation can be spread among all sizes, formally

U Ω = {a ∈ × i∈[n] [ā i , āi +â i ] | i∈[n] (a i -ā i ) ≤ Ω}.
The main purpose of this paper is to prove the following two theorems.

Theorem 1.1. There exists a 2-approximation algorithm for RBP with uncertainty set U Ω . Theorem 1.2. There exists a min(2Γ, 4.5)-approximation algorithm for RBP with uncertainty set U Γ .

1.1. Literature review. Uncertainty set U Γ has been widely used in robust combinatorial optimization with a constant number of constraints because the set essentially preserves the complexity properties of the nominal problem, and to some extent, its approximability properties as well. While the results for this model were first developed for min-max problems and uncertain objective functions [START_REF] Bertsimas | Robust discrete optimization and network flows[END_REF], they were later extended to uncertain constraints independently in [START_REF] Álvarez-Miranda | A note on the Bertsimas & Sim algorithm for robust combinatorial optimization problems[END_REF][START_REF] Goetzmann | Optimization over integers with robustness in cost and few constraints[END_REF]. Our second uncertainty set (used in [START_REF] Gounaris | The robust capacitated vehicle routing problem under demand uncertainty[END_REF][START_REF] Pessoa | Branch-and-cut-and-price for the robust capacitated vehicle routing problem with knapsack uncertainty[END_REF][START_REF] Tadayon | Algorithms and complexity analysis for robust single-machine scheduling problems[END_REF], among others), U Ω benefits from similar positive results as U Γ ; see [START_REF] Poss | Robust combinatorial optimization with knapsack uncertainty[END_REF].

The above positive complexity results (e.g., [START_REF] Goetzmann | Optimization over integers with robustness in cost and few constraints[END_REF][START_REF] Poss | Robust combinatorial optimization with knapsack uncertainty[END_REF]) imply, for instance, that under mild assumptions there exists a fully-polynomial time approximation scheme (FP-TAS) for the robust knapsack problem with uncertain profits and uncertain weights in U Ω or U Γ . Interestingly, these positive results may not be applied to most scheduling problems due to the presence of nonlinearities of the objective function in the uncertain data parameters. Also, these general results cannot be applied to the bin packing problem because it involves a non-constant numbers of robust constraints. Further, in the current paper it is proved that approximation results of standard bin packing in the form of asymptotic approximation scheme, cannot carry over (unless P = N P) to the robust counterpart with U Γ . While in previous papers [START_REF] Bougeret | Robust scheduling with budgeted uncertainty[END_REF][START_REF] Bougeret | Approximation results for makespan minimization with budgeted uncertainty[END_REF] (with authors in common) we provided approximability results on robust scheduling, no such results have yet been proposed for the bin packing problem, the only previous work focusing on numerical algorithms [START_REF] Song | The robust machine availability problem-bin packing under uncertainty[END_REF]. The purpose of this paper is to fill these gaps, as we present constant-factor approximation algorithms the bin packing problem, both for U Ω and U Γ .

1.2. Notation, problems definitions, and next-fit algorithm. In this paper we consider two special cases of RBP. In the first one, ΓRBP, the input is

I = (a, â, Γ) ∈ [0, 1] n × [0, 1] n × N where n ∈ N, and U = U Γ . In the second one, ΩRBP, the input is I = (a, â, Ω) ∈ [0, 1] n × [0, 1] n × [0, 1] where n ∈ N and U = U Ω .
Let us now define some notation that may be required for formally stating the ΓRBP and ΩRBP problems. Given n ∈ N, sets {0, 1, . . . , n} and {1, . . . , n} are respectively denoted [n] 0 and [n]. Set {i, . . . , j} is denoted by i, j . Given a vector v ∈ [0, 1] n and a subset X ⊆ [n], we define v(X) = i∈X v i . Given two vectors a ∈ [0, 1] n , â ∈ [0, 1] n and a subset of items X ⊆ [n], we define âΩ (X) = min{â(X), Ω}, 1 U Γ is often defined alternatively in the literature, as the polytope {a ∈

× i∈[n] [ā i , āi + âi ] | i∈[n] (a i -āi )/â i ≤ Γ}.
For the bin packing problem, one readily verifies using classical arguments that the two definitions lead to the same optimization problem. Γ(X) as the set of Γ items in X with largest â values (ties broken by taking smallest indices), or Γ(X) = X if |X| < Γ, and âΓ (X) = â(Γ(X)). Accordingly, we define the fill of a bin b ⊆ [n] as f Γ (b) = a(b) + âΓ (b) for set U Γ , and f Ω (b) = ā(b) + âΩ (b) for set U Ω . The fill of a bin for a general uncertainty set U is denoted as f U (b) = max a∈U a(b).

Consider the following example. We are given an ordered set of pairs (ā i , âi ), X = {(0.3, 0.2), (0.4, 0.2), (0.3, 0.1), (0.2, 0.5)} with Γ = 2 and Ω = 0.3. Thus, Γ(X) = {(0.3, 0.2), (0.2, 0.5)}, ā(X) = 1.2, âΓ (X) = 0.7, and f Γ (X) = 1.9. Similarly, âΩ (X) = 0.3 and f Ω (X) = 1.5. Now, observe that max a∈U i∈bj a i ≤ 1 (the constraint required in RBP) is equivalent to f U (b) ≤ 1, and thus to f Γ (b j ) ≤ 1 for ΓRBP and f Ω (b j ) ≤ 1 for ΩRBP. For example in ΓRBP, f Γ (b j ) ≤ 1 simply means that the total nominal (a) size of the items plus the deviating size (â) of the Γ largest (in â values) items must not exceed one. Then, the two optimization problems studied in this paper can be equivalently formulated as follows.

Γ-robust bin packing (ΓRBP)

Input: I = (a, â, Γ) ∈ [0, 1] n × [0, 1] n × N. Output: A solution that is a partition of [n] into k bins b 1 , . . . , b k such that f Γ (b j ) ≤ 1 for each j ∈ [k] Minimize: k Ω-robust bin packing (ΩRBP) Input: I = (n, a, â, Ω) where n ∈ N, a ∈ [0, 1] n , â ∈ [0, 1] n , and Ω ∈ [0, 1]. Output: A solution that is a partition of [n] into k bins b 1 , . . . , b k such that f Ω (b j ) ≤ 1 for each j ∈ [k] Minimize: k
The optimal solution value or cost of either problem is denoted by OPT(I) = k * (I may be omitted when the instance is clear from the context) and a corresponding optimal solution is denoted by s

* = {b * 1 , b * 2 , . . . , b * k * }.
We introduce in Algorithm 1 a variant of the next-fit-decreasing algorithm. In Section 2 we prove that there is no asymptotic approximation for ΓRBP with a factor less than 3/2 (unless P = N P). In Sections 3 and 4, we analyze approximation-factor guarantees for Next-Fit-Decreasing in the case of ΩRBP and ΓRBP, respectively. For ΩRBP, using ordering (3.1) (non-increasing ordering on âi āi ) the ratio is equal to 2. For ΓRBP, using ordering (4.1) (non-increasing ordering on âi ), the ratio is bounded by 2Γ. As Theorem 4.5 shows that neither ordering (3.1) or (4.1) leads to a constant ratio using Next-Fit-Decreasing, this raises the question of the existence of a constant approximation for ΓRBP. We answer the question in Section 5 by providing a dynamic programming algorithm (DP) giving a ratio of 4.5 for ΓRBP and any Γ ∈ N, which is our main result.

2. Innaproximability. First we establish that as opposed to classical bin packing that has asymptotic approximation schemes [START_REF] Vega | Bin packing can be solved within 1+epsilon in linear time[END_REF][START_REF] Karmarkar | An efficient approximation scheme for the one-dimensional bin-packing problem[END_REF], such approximation schemes are not possible (unless P = N P) for ΓRBP.

Lemma 2.1. Suppose Γ = 1. For any K ≥ 1, it is N P-complete to decide whether ΓRBP has a solution using 2K bins or whether at least 3K bins are required.

Proof. Consider an instance I, of the N P-Complete 2-partition decision problem [START_REF] Garey | Computers and Intractability: A Guide to the Theory of NP-Completeness[END_REF], given by a set of distinct positive rational numbers whose sum is 2 (in particular each instance given by a set of N integers can be normalized to get such a set of rational numbers), so I = {s 1 , s 2 , ..., s N }, where s 1 +s 2 +...+s N = 2 and 0

≤ s i ≤ 1 for i ∈ [N ].
In particular, it is N P -complete to decide the question whether the given set of numbers can be partitioned into two subsets having a sum of one in each subset (that is to decide whether I is a "YES-instance"). Let s min = min i∈[N ] {s i } and let 0 < δ < s min (assume that δ is a rational number). Now choose some (arbitrarily large) integer K ≥ 1, and create K sets (each of size N ) of bin packing items as follows. For k ∈ [K], let ε k = δ k . Elements of the K copies of the N numbers will be denoted by index set pairs in [N ] × [K]. The nominal and deviation sizes of the ΓRBP instance items are defined for i ∈ [N ] and k ∈ [K] as

a i,k = ε k • s i and âi,k = 1 -ε k , respectively.
Note that for each k, i∈[N ] a i,k = 2ε k , and

a i,k ≤ ε k for each i ∈ [N ]. Also note, that by construction, for all i ∈ [N ], a i,k > ε k+1 , for 1 ≤ k < K.
First it can be observed that if I is a YES-instance of the 2-partition problem, then the constructed ΓRBP instance 2K bins: K pairs of bins, where in each of the two bins of pair k ∈ [K], the total nominal is ε k and the deviation is 1k .

Suppose now that I is a NO-instance of 2-partition. Let b 1 , . . . , b L be a ΓRBP solution for some L > 1. To show that L ≥ 3K, first define for each k ∈ [K],

I k = {(i, k), i ∈ [N ]
} to be referred to as the set of items of type k, and observe that no bin can contain items of two different types. Indeed, suppose by contradiction that for some k, k

∈ [K] satisfying k = k (say k > k ), and i, i ∈ [N ], item (i, k) is packed with item (i , k ). This implies a i,k + a i ,k + max{â i,k , âi ,k } ≥ a i,k + âi ,k > ε k + (1 -ε k ), a contradiction. Thus, none of the bins b 1 , . . . , b l contain items of two different types, that is for each l ∈ [L], there exists some k ∈ [K] such that b l ⊆ I k . Since I is a NO-instance, each I k is packed in at least 3 distinct bins out of b 1 , . . . , b L , implying L ≥ 3K.
The following theorem is immediate from Lemma 2.1.

Theorem 2.2. Unless P = N P there is no asymptotic approximation for ΓRBP of ratio ρ < 3/2.

Note that this result is in contrast with the approximation schemes that exist for classical bin packing [START_REF] Karmarkar | An efficient approximation scheme for the one-dimensional bin-packing problem[END_REF][START_REF] Vega | Bin packing can be solved within 1+epsilon in linear time[END_REF]. Furthermore, recall that in contrast to Theorem 2.2, for the classical bin packing, the first-fit decreasing algorithm packs the items of any input I into at most 11 9 OP T (I) + 6/9 bins [START_REF] Dósa | The tight bound of first fit decreasing bin-packing algorithm is FFD(I) <= 11/9 OPT(I) + 6/9[END_REF].

3. Next-fit-decreasing for ΩRBP. Unlike the classical bin packing problem, executing Next-Fit-Decreasing on arbitrarily ordered items can lead to arbitrarily bad solutions. For example, given with 0 < ≤ 1 2n , consider an instance with Ω = 1 -, and items ((2 , 0), (0, 1 -), . . . , (2 , 0), (0, 1 -)), where item i ∈ [n] is denoted by the pair (ā i , âi ). Using this ordering, Next-Fit-Decreasing will create n/2 bins b j with f Ω (b j ) > 1 for any j ∈ [n] (which will be turned into n bins {b 1 j , b 2 j }), whereas the optimal solution uses 2 bins. This example also illustrates that, unlike in the standard bin packing, the total size argument no longer apply to the robust counterpart as having f Ω (b j ) > 1 for any j ∈ [n] does not imply a large (depending on n) lower bound on the optimum.

Next, we consider an ordering of the items such that

(3.1) â1 /ā 1 ≥ • • • ≥ ân /ā n ,
and recall that k is the number of bins opened in Step 1 of Next-Fit-Decreasing.

* = {j ∈ [k * ] | â(b * j ) > Ω}. Let A = i∈[n] (ā i + âi ) = j∈[k ] (ā(b j ) + â(b j )) = j∈[k * ] ā(b * j ) + â(b * j ) . Let G denote the first |G * | bins opened in Step 1 of Next-Fit-Decreasing. If k ∈ G then clearly k ≤ k * . Otherwise, it can be observed that for each l ∈ G, ā(b l ) > 1 -Ω (as ā(b ) + âΩ (b ) > 1 and âΩ (b ) ≤ Ω) and 1 -Ω ≥ max j∈G * ā(b * j ) (as f Ω (b * j ) ≤ 1)
. Thus, j∈G ā(b j ) > j∈G * ā(b * j ) and so by the assumed ordering (3.1) of the items, following a standard knapsack argument, j∈G â(b j ) 4. Next-fit-decreasing for ΓRBP. From now on, we focus on problem ΓRBP. Remark first that using an arbitrary ordering leads to arbitrarily bad solutions, considering Γ = 1 and the same items ((2 , 0), (0, 1 -), . . . , (2 , 0), (0, 1 -)) as in the previous section. Thus, we consider here an ordering of the items such that

> j∈G * â(b * j ). Letting Ḡ = [k ] \ G and Ḡ * = [k * ] \ G * , it follows that j∈ Ḡ (ā(b j ) + â(b j )) = A - j∈G (ā(b j ) + â(b j )) ≤ A - j∈G * a(b * j ) + â(b * j ) = j∈ Ḡ * ā(b * j ) + â(b * j ) (equality may hold throughout if G * = ∅). Further, for each j ∈ Ḡ \ {k }, ā(b j ) + â(b j ) ≥ f Ω (b j ) > 1 and for each j ∈ Ḡ * , ā(b * j ) + â(b * j ) ≤ 1. Therefore, | Ḡ| ≤ j∈ Ḡ (ā(b j ) + â(b j )) ≤ j∈ Ḡ * ā(b * j ) + â(b * j ) ≤ | Ḡ * |
(4.1) â1 ≥ • • • ≥ ân .
The main result of this section is the following. We also note that this result has been improved compared with the result that appears in the preliminary extended abstract version of this paper [START_REF] Basu Roy | Approximating robust bin-packing with budgeted uncertainty[END_REF]. Let M ⊆ s be the set of bins that contain only items that either do not deviate or are the smallest deviations, i * j for some j

∈ [k * ]. So for each b ∈ M , b ⊆ [n] \ k * j=1 Γ(b * j ) ∪ {i * 1 , . . . , i * k * }. In what follows, we bound |M | and |s \ M | by multiples of k * . Lemma 4.3. |s \ M | ≤ (Γ -1)k * Proof. As each bin b ∈ s \ M contains at least one item from k * j=1 Γ(b * j ) \ {i * 1 , . . . , i * k * }, and as | k * j=1 Γ(b * j ) \ {i * 1 , . . . , i * k * }| ≤ (Γ -1)k * , the result is immediate.
The case of M , which is slightly more involved, is addressed in the following lemma.

* * < k . Then, since k j=1 b j ⊆ k * * j=1 b * j , j∈[k ] ā(b j ) ≤ j∈[k * * ] ā(b * j ) ≤ k * * - j∈[k * * ] âΓ (b * j ).
We now show that j∈ ā(b j ) > 1 -

j∈[k] âΓ (b j ) ≥ 1 - j∈[k] âΓ (b * j ) ≥ j∈[k] ā(b * j ) ≥ j∈[k] ā(b * j \ Γ(b * j ) ∪ {i * j }).
The above inequality and the ordering of the items imply that âmax ≡ max

i∈ k * * j=k+1 (b * j \Γ(b * j )∪{i * j }) âi ≥ max i∈b k+1
âi .

Since âi * k+1 =â max it follows that âΓ (b * k+1 ) ≥ âΓ (b k+1 ). Together with the induction hypothesis it implies that j∈

[k+1] âΓ (b j ) ≤ j∈[k+1] âΓ (b * j ). Note that âΓ (b j ) ≤ 1 for j = 1, . . . , k ; otherwise if âΓ (b j ) > 1 for some j ∈ [k ],
then there is some i ∈ b j that has âi > 1/Γ. But then i ∈ b * l for some l ∈ [k * * ] and since i is either a nondeviating item and âΓ (b * l ) > 1, or i = i * l and it implies that b * l may contain only Γ -1 deviating items, thereby establishing a contradiction. Now, by definition of the algorithm,

j∈[k ] ā(b j ) > j∈[k -1] ā(b j ) > k -1 - j∈[k -1] âΓ (b j ) ≥ k * * - j∈[k * * ]
âΓ (b * j ). Proof of Theorem 4.1. After step 2 of the algorithm the total number of bins is at most 2Γk * , concluding the proof of Theorem 4.1.

To complete the analysis, we establish the following lower bound on the approximation ratio of Next-Fit-Decreasing.

( , 1/Γ -δ 1 ) (0, 1/Γ -δ 1 ) . . . (0, 1/Γ -δ 1 ) . . . . . . . . . . . . ( , 1/Γ -δ Γ ) (0, 1/Γ -δ Γ ) . . . (0, 1/Γ -δ Γ )
where ≤ 1/Γ and

δ 1 ≤ • • • ≤ δ Γ < /Γ. On the one hand, + Γ • (1/Γ -δ l ) > 1 for each l ∈ [Γ]
, so step 1 of Next-Fit-Decreasing outputs Γ bins. On the other hand, an optimal solution can pack all the items above except the ones in the first column into a single bin because Γ • 1/Γ -δ 1 ≤ 1. Further, the total weight of the first Γ/2 items of the first column sums up to Γ/2

• (1/Γ + ) - Γ/2 l=1 δ l ≤ 1 - Γ/2
l=1 δ l ≤ 1, and similarly for the last Γ/2 items, so an optimal solution may pack the first column using two bins.

This instance can be adapted to establish a lower bound for the approximation ratio of Next-Fit-Decreasing when items are ordered according to (3.1). We consider an example that yields a lower bound on the approximation ratio of Next-Fit-Decreasing in solving ΓRBP when the items are ordered according to (3.1). For some c ≥ Γ 2 , = The optimal solution opens 3 bins, 2 bins to store the items of the first column and another bin to store the rightmost Γ -1 columns. Although in this example all items i ∈ [n] are set to have ratios âi /ā i = c, the example can be extended in a straightforward manner with slight perturbations of the item sizes so that the ratios will be strictly decreasing for the items ordered from left-to-right and top-to-bottom in this matrix.

We conclude the section by emphasizing our results for the special cases of Γ = 1 and Γ = 2, where the straight-forward Next-Fit-Decreasing algorithm with ordering (4.1) obtains the best approximation guarantees using the analysis of the current paper. In the section that follows we consider constant-factor approximation factors that do not depend on the parameter Γ.

5.

A constant-factor approximation algorithm for ΓRBP. The algorithm presented in this section relies on three main ideas. First, we show in Section 5.1 that we can restrict ourselves to instances of ΓRBP with small items; that is,

a i ≤ 1 Γ and âi ≤ 1 Γ for each i ∈ [n]
. Specifically, we show how to convert any ρ-approximation algorithm for the latter special case of the problem with small items into a general (ρ + ρ bp )-approximation for ΓRBP, where ρ bp is an approximation factor guarantee of some algorithm for classical bin packing. First note that having only small items, any set of Γ/2 items can always be packed together into a bin. This fact is stated as the following observation.

Observation 1. Given an instance I to the ΓRBP satisfying âi ≤ 1/Γ and a i ≤ 1/Γ for each i ∈ [n], any subset X ⊆ [n] can be packed in at most 2|X| Γ bins.

Next, we introduce in Section 5.2 variants of ΓRBP where items are packed into bins as in standard bin packing but some of the items are placed in two designated special bins. The deviations will be ignored so that they are irrelevant in these two bins, while in each of the other regular bins only a single item can deviate. One of the special bins that will be referred to as the trash bin, cannot contain more than k(Γ -1) items, where k is the number of regular bins used in the solution. The trash contains and "mimics" the Γ -1 deviating items of each bin in an optimal solution. Following Observation 1, items in the trash can be packed into at most 2k additional bins. The problem with trash remains hard because of the capacity of the regular bins, so we focus on almost feasible solutions, which are allowed to exceed each regular bin by one item. We show in that section how finding almost feasible solutions no worse than the optimal solution for the problem with trash leads to an approximation algorithm for ΓRBP with small items.

Finally, we present in Section 5.3 a dynamic programming (DP) algorithm that outputs an almost feasible solution using a number of bins that is at most that of the optimal solution of the original ΓRBP problem. Essentially, for each bin the DP algorithm guesses a particular item that deviates. Then it greedily packs the remaining items with largest nominal values into the trash.

We maintain throughout the section the assumption that the items are ordered according to (4.1).

5.1.

Robust bin packing with small items. We define ΓRBP with small values as the ΓRBP problem restricted to inputs where for any i ∈ [n], a i ≤ 1 Γ and âi ≤ 1 Γ . The following proposition motivates our focus on ΓRBP with small values.

Proposition 5.1. Any polynomial ρ-approximation algorithm for ΓRBP with input satisfying âi ≤ 1/Γ and a i ≤ 1/Γ for each i ∈ [n], can be turned into a polynomialtime (ρ + ρ bp )-approximation algorithm for ΓRBP. 2Proof. Given an instance I of ΓRBP, we define the small items S = {i ∈ [n] : a i ≤ 1/Γ and âi ≤ 1/Γ} and the large item as B = [n] \ S. We use the given ρapproximation algorithm to pack S into k S bins, so that k S ≤ ρOPT(I). Then, we observe that in any packing of B, each bin contains at most Γ items, so that all items deviate in these bins. Hence, the least number of bins needed to pack items in B is is given by a solution that is optimal to the standard bin packing problem for this same set of items B where the size of each item i ∈ B is a i = a i + âi . Let I denote this instance of standard bin packing So, OPT bp (I ) ≤ OPT(I) (where OPT bp denotes the optimal objective value of standard bin packing), and observe that any optimal solution of standard bin packing for items B is a solution that is feasible for the ΓRBP instance I using the same number of bins. Using a ρ bpapproximation algorithm for standard bin packing to pack B in k B bins, it follows that k B ≤ ρ bp OPT bp (I ) = ρ bp OPT(I) ≤ ρ bp OPT(I). Thus, the packing of B and S is a packing of

I satisfying k S + k B ≤ (ρ + ρ bp )OPT(I).
Notice that instances with small items are not easier to approximate by Next-Fit-Decreasing as illustrated by the instance defined by (4.2), in Section 4, whose items satisfy a i , âi ≤ 1/Γ for all i ∈ [n].

5.2. Bin-packing with trash. For any X ⊆ [n], we define ãΓ (X) = Γâ 1 (X) (ã Γ (X) is Γ times the largest deviating value of an item in X) and f (X) = a(X) + ãΓ (X). We introduce next a decision problem ΓRBP-T that is related to ΓRBP.

ΓRBP-T (Robust bin packing with trash) Input:

(I, k, t) where I is an instance of ΓRBP and k, t ∈ N. Output: 'Yes' if a solution exists, which is a partition of the set of items into k + 1 sets b 1 , . . . , b k and T (called the trash) such that:

• f (b j ) ≤ 1 for each j = 1, . . . , k • |T | ≤ t and 'No' otherwise.
Notice that although the input of ΓRBP-T is assumed to include only small items, it is possible to have an item i ∈ [n] such that f ({i}) > 1, implying that i ∈ T . The following two lemmas suggest how the decision problem ΓRBP-T may be used to determine an approximate solution of ΓRBP. Proof. Given an optimal solution with objective value k * of ΓRBP we create a solution to ΓRBP-T problem as follows. For some (arbitrary) j ∈ [k * ], let b * j be a bin of the considered optimum. Let N j = b * j \ Γ(b * j ) (the non-deviating items of b * j ). For j = 1, . . . , k * , let

X j = max(Γ(b * j )) |Γ(b * j )| = Γ ∅ otherwise.
We define b j = N j ∪ X j , and adjoin items of The next lemma establishes that yes-instances of ΓRBP-T can be used to construct solutions of ΓRBP using at most a number bins that is a constant factor of the number of bins used by ΓRBP-T. Lemma 5.3. For any instance I of ΓRBP satisfying âi ≤ 1/Γ and a i ≤ 1/Γ for each i ∈ [n], and integer k, given a yes-instance of (I, k, Γk) of ΓRBP-T, we can compute in polynomial time a solution of 3k bins for I.

Proof. Given a solution b 1 , . . . , b k , T for (I, k, Γk) of ΓRBP-T the bins remain feasible in ΓRBP as for each j ∈

[k], f Γ (b j ) = a(b j ) + âΓ (b j ) ≤ a(b j ) + ã(b j ) = f (b j ).
Then, Observation 1 implies that the trash T can be packed into kΓ/(Γ/2) ≤ 2k additional bins.

Notice that the trash T in the instance constructed in this lemma has |T | = Γk, which exceeds the |T | = (Γ -1)k in Lemma 5.2. The additional k "slots" are necessary for deciding ΓRBP-T as will be illustrated in the analysis that follows.

In order to develop an algorithm and in particular a dynamic program (DP) for deciding ΓRBP-T it is convenient to consider an optimization variant of ΓRBP-T, to be called G-ΓRBP-T for "generalized robust bin packing with trash". This variant is defined in the following for a fixed instance I of ΓRBP and a given integer k (so that they are not considered a part of the input).

Generalized robust bin packing with trash (G-ΓRBP-T) Input:

I = (q, t, ), where q ∈ [n], t ∈ [(Γ -1)k] 0 , and ∈ [k + 1].

Output:

A feasible solution s is a partition of q, n into k -+ 3 sets, given as a triple (L, {b j : j ∈ , k }, T ), such that • for any j ∈ , k , f (b j ) ≤ 1 (the k -+ 1 regular bins must satisfy the fill constraints of ΓRBP-T) • |T | ≤ t (we only allow t items in the trash) • min(b ) = q (meaning that the deviating item of b is q)

Minimize: c(s) = a(L)
The objective of G-ΓRBP-T is to pack a part (defined by q, k ) of the instance I for ΓRBP-T given a fixed budget of resources (the number of bins and the size of the trash) while minimizing the sum only of nominal sizes of items in the leftover itemset L. The last constraint (the deviating item of b is q), which may appear somewhat artificial, will allow determining optimal solutions of ΓRBP-T from optimal solutions of G-ΓRBP-T, by carefully enumerating possibilities of largest deviating item to be packed each bin in an intelligent way; considering only those possibilities corresponding to solutions that are close to being feasible (a notion to be defined more precisely in the following) and whose objective value is at most that of an optimal solution. This enumeration scheme will be shown to be efficiently solvable by a DP. For convenience, the objective value of infeasible solutions s, including for example

s = (∅, ∅, ∅), is defined as c(s) = ∞.
The capacity constraints of the bins make G-ΓRBP-T hard to solve in general. Hence, following the spirit of Next-Fit-Decreasing introduced previously, we introduce below almost feasible solutions, which can exceed the capacity of each bin by one item. Definition 5.4 (almost feasible solution). We say that a bin b exceeds by at most one item iff f (b) > 1 and f (b \ {i}) ≤ 1 where i = max(b). Given an input I = (q, t, ) of G-ΓRBP-T, we say that a solution is almost feasible iff all the G-ΓRBP-T constraints are satisfied, except that for any j ∈ , k , we allow that b j exceeds by at most one item instead of f (b j ) ≤ 1. Definition 5.5 (an optimal almost-feasible solution). Given an input I = (q, t, ) of G-ΓRBP-T, we say that a solution s = (L, {b 1 , . . . , b k }, T ) is an optimal almost-feasible solution iff s is almost feasible with c(s) ≤ OPT(I ).

The relation between G-ΓRBP-T and ΓRBP-T is characterized in the following two lemmas. Let OPT(I ) be the optimal solution cost of G-ΓRBP-T instance I (notice that an optimal solution must also be feasible for the given problem).

Lemma 5.6. For any input I of ΓRBP and k such that (I, k, (Γ -1)k) is a yes input of ΓRBP-T, there exist a positive integer q ∈ [(Γ -1)k] such that OPT(q, (Γ -1)k -(q -1), 1) = 0.

Proof. Consider a ΓRBP-T yes-instance (I, k, (Γ -1)k) and corresponding partition solution {b 1 , . . . , b k }, T . Let q = min k j=1 b j (the item with smallest index that is packed in a bin). Let t = (Γ -1)k -(q -1). By definition of q and (4.1), items in [q -1] must be in T , and thus it means that the triple (∅, {b 1 , . . . , b k }, T \ [q -1]) is a feasible solution of G-ΓRBP-T (q, t, 1) with cost 0. Lemma 5.7. Let us fix I an input of ΓRBP and k an integer. For any q ∈ [(Γ -1)k], t = (Γ -1)k -(q -1), given an almost feasible 0 cost solution for G-ΓRBP-T instance I = (q, t, 1), a solution for ΓRBP-T instance (I, k, Γk) can be determined in polynomial time.

Proof. Let (∅, {b 1 , . . . , b k }, T ) be an optimal solution of G-ΓRBP-T instance I = (q, t, 1) with cost 0. For j ∈

[k], let b j = b j \ max(b j ). Let T = T ∪ [q -1] ∪ k j=1 max(b j
). We now have f (b j ) ≤ 1 for any j ∈ [k] (as b j exceeds by at most one item), and |T | ≤ Γk, so {b 1 , . . . , b k , T } is a yes-instance for ΓRBP-T (I, k, Γk) thus concluding the proof.

The above lemmas imply that any algorithm outputting an (approximately) optimal almost-feasible solution for G-ΓRBP-T can be used to devise approximation algorithm for ΓRBP with small values, which is illustrated by Algorithm 2. Proof. Let A be a G-ΓRBP-T algorithm that is guaranteed to output an optimal almost-feasible solution given instance (q, t, ), and whose output will be denotes by c(A(q, t, )). Let k * = OPT(I) be the optimal value of ΓRBP. Lemma 5.2 implies that (I, k * , (Γ -1)k * ) is a yes-instance of ΓRBP-T, so Lemma 5.6 implies that there exists q ∈ [(Γ -1)k * ] and t = (Γ -1)k * -(q -1) for which OPT(q, t, 1) = 0. Thus, the assumption implies that A outputs an almost feasible solution s that satisfies c(s) ≤ 0, so that c(s * ) ≤ k * . Applying Lemma 5.7, we can compute a solution of (I, k * , Γk * ) for ΓRBP-T, which corresponds to a solution using 3k * bins for ΓRBP following Lemma 5.3.

1 s * ← ([n], ∅, ∅), lb ← 1, ub ← n 2 while lb = ub do 3 k ← lb+ub
The next section describes a dynamic programming algorithm that outputs an optimal almost-feasible solution for G-ΓRBP-T.

5.3.

A DP algorithm for G-ΓRBP-T. We now define a DP scheme that given instance I = (q, t, ) provides an almost feasible solution s with c(s) ≤ OPT(I ). We start by providing a gentle description before formally defining the algorithm and proving its correctness. Let s * be an optimal solution for I , with bins ordered in non-increasing deviation size. The DP algorithm starts by enumerating g = (q , t ) ∈ q + 1, n × [t] 0 . One of the enumerated g must necessarily correspond to the optimal solution s * in the following sense:

• q = min(b * +1 )

• t is the number of items trashed from X in s * , where X = q, q -1 Because of the ordering of the bins in the optimal solution s * , the items of X must be packed in bins b * , L * or T * . As the DP algorithm mimics the optimal solution, it packs X in b , L and T . Specifically, the DP algorithm:

• Packs q to b (as required by the corresponding constraint of G-ΓRBP-T).

• Packs the remaining t largest nominal value items of X in the trash.

• Packs the remaining items of X into b until f (b ) > 1 or X = ∅.

• Packs the remaining items of X into L until X = ∅. We discuss next where the others items (of q , n ) are packed. Notice that in s * , bin b * may contain items of q , n , and thus the DP algorithm may also have to pack items of q , n into b . To allow for that possibility, we postpone the decision of which items of q , n to pack into b . Specifically, let ∆ be the size of the empty space in b after packing X as described above, and let L X = L ∩ X . After the previous steps, the DP algorithm makes a recursive call to get a solution s that packs q , n into regular bins, a trash, and a leftover itemset b0 . So far solution s has not benefited from the empty space ∆ . However, we can unpack items from b0 to b while ensuring that these items do not deviate in b (as all these items have index greater than q).

Our DP scheme is defined by Algorithm 3. An iteration of this algorithm is further illustrated in Figure 1.

1 s ← (∅, ∅, ∅) // Where c((∅, ∅, ∅)) = ∞ 2 if = k then 3 X ← q, n 4 T ← {min(t, n -q) largest items in terms of nominal value in X } 5 Pack X \ T in b k until f (b k ) > 1 or X \ T = ∅ 6 Pack the remaining items in L return: s = (L, {b k }, T ) 7 for g = (q , t ) ∈ q + 1, n × [min(t, q -q)] 0 do 8 Pack q in b 9 X ← q + 1, q -1 10 T ← {t largest items in terms of nominal value in X } 11 Pack X \ T in b until f (b ) > 1 or X \ T = ∅ 12
Pack the remaining items in X \ T in L 13 s = ( b0 , { b +1 , . . . , bk }, T ) ← DP(q , t -t , + 1) . We define similarly L X = L ∩ X . Let c 0 = a(L X ). Notice that c 0 corresponds to the total non deviating size of items packed in L after Step 12 of Algorithm 3. The following lemma establishes the correctness of the algorithm. Lemma 5.9. For any I = (q, t, ) input of G-ΓRBP-T, DP(I ) determines an optimal almost feasible solution.

Proof. The proof is by induction on . = k.. Consider a solution s = (L, {b k }, T ) output by the algorithm and optimal solution s *

= (b * 0 , {b * k }, T * ). In step 5, either b = {q} ∪ X \ T , in which case L = ∅ and ā(L) ≤ ā(b * 0 ), or f (b k ) > 1 ≥ f (b * k ). Then, as q = min(b k ) = min(b * k ) (since both solutions must satisfy this constraint) it implies that Γâ 1 (b * k ) = Γâ 1 (b k ), and so ā(b k ) > ā(b * k ).
Then, by the definition of T (in step 4 of the algorithm), also ā(T ) ≥ ā(T * ), so it must be that ā(L) < ā(b * 0 ). Inductive step.. Suppose now that, for some ˆ ∈ [k -1], the lemma holds for = ˆ + 1 and that s * = (b * 0 , {b * ˆ , . . . , b * k }, T * ) is optimal for G-ΓRBP-T (q, t, ˆ ). Let q * = min(b ˆ +1 ), let X * = q, q * -1 and observe that X * = ∅. Next, observe that (b * 0 \ X * , {b * ˆ +1 , . . . , b * k }, T * \ X * ) is optimal for G-ΓRBP-T (q * , t, ˆ + 1) where t = |T * \ X * | ≤ t. Otherwise, there would exist some feasible ( Ľ, { bˆ +1 , . . . , bk }, Ť ) that has ā( Ľ) < ā(b * 0 \ X * ) and then ( Ľ ∪ (b * 0 ∩ X * ), {{q} ∪ (b ˆ ∩ X * ), bˆ +1 , . . . , bk }, Ť ∪ (T * ∩ X * )) would be feasible for (q, t, ˆ ) with ā( b ∪ (b * 0 ∩ X * )) < ā(b * 0 ), thereby contradicting the optimality of s * .

By the inductive hypothesis ā( L) ≤ ā(b * 0 \ X * ) for some almost feasible ( L, { bˆ +1 , . . . , bk }, T ) that is output by the algorithm for (q * , t, ˆ + 1), where q * > q and t = |T * \ X * | ≤ t. Therefore,

(5.1) ā(b * 0 ) = ā(b * 0 \ X * ) + ā(b * 0 ∩ X * ) ≥ ā( L) + ā(b * 0 ∩ X * ).
Let us now consider the iteration g = (q * , t -t) of the main loop (starting in step 7) and L and b ˆ the corresponding bins of s g . By steps 8, 11 and 14 of the algorithm,

L \ L = (X * \ T ) \ b ˆ . In line 11 of the algorithm if f (b ˆ ) ≤ 1, then L ⊆ L and L \ L = ∅. Otherwise f (b ˆ ) > 1 ≥ f (b * ˆ ).
Recalling that both solutions must satisfy the given constraint, q = min(b

ˆ ) = min(b * ˆ ), it implies that Γâ 1 (b ˆ ) = Γâ 1 (b * ˆ ). Thus, ā(b ˆ ) > ā(b * ˆ ). Also, by definition of T , ā(T ) ≥ ā(T * ∩ X * ). Thus, ā(L \ L) ≤ ā(b * 0 ∩ X * ). Combined with (5.1) it finally implies that ā(b * 0 ) ≥ ā( L) + ā(b * 0 ∩ X * ) ≥ ā(L \ L) + ā( L) = ā(L).
Lemma 5.10. Suppose that a i , âi ≤ 1/Γ for all i ∈ [n]. Then, Algorithm 2, with Algorithm 3 as a subroutine for determining almost-feasible optimal solutions for G-ΓRBP-T, is a 3-factor approximation algorithm for ΓRBP. Further, the running time complexity bound of this algorithm is in O(n 6 log(n)).

Proof. Proposition 5.8 and Lemma 5.9 applied to Algorithm 3 imply that there exists a 3-approximation for ΓRBP with small values. The different number of inputs of the DP is O(n 3 ), and the running time for a fixed input is O(n 3 ). Then, since the binary search of Algorithm 2 takes O(log n) iterations it follows that the total running time is in O(n 6 log n). Now, Lemma 5.10 together with Proposition 5.1 imply the following theorem using the ρ bp = 3 2 -approximation first-fit decreasing (FFD) algorithm, whose running time complexity bound is clearly dominated by that of the algorithm for ΓRBP with small values, for standard bin packing; see for example [START_REF] Simchi-Levi | New worst-case results for the bin-packing problem[END_REF].

Theorem 5.11. There exists a 4.5-approximation algorithm for ΓRBPwith an O(n 6 log(n)) runtime complexity bound.

In particular, combining FFD to to pack large items i ∈ [n], for which with āi > 1/Γ or âi > 1/Γ and Algorithm 2 using Algorithm 3 as a subroutine to pack the other (small) items, is an algorithm that satisfies the claim of Theorem 5.11. Also note that while not a focus of the current paper, in the asymptotic setting (as n tends to be large) using instead an asymptotic FPTAS for packing the large items, for > 0, an asymptotic approximation of 4 + could be guaranteed in running time that is polynomial in n and 1/ . Finally, although this result establishes a constant factor approximation for our problem when Γ is a part of the input, it can be observed that in the special case that Γ ≤ 2, the next-fit-decreasing approximation established by Theorem 4.6 may be preferred as a practical and fast O(n log n) algorithm with a 2Γ ≤ 4 approximation guarantee.

Conclusion.

This paper considered the bin-packing with item size uncertainty, following the robust optimization approach with two widely used variants of budget uncertainty sets, U Γ and U Ω . We have shown that problem ΓRBP belongs to a different complexity class than its deterministic counterpart, in particular that no asymptotic approximation scheme exists for that problem, unless P = N P. We have further provided constant factor approximation algorithms for problems ΩRBP and ΓRBP. While the algorithm devised for ΩRBP is based on a rather natural extension of next-fit-decreasing, the latter provided only non-constant ratios for ΓRBP. Instead we are able to devise a constant-ratio for ΓRBP through a more involved algorithm, separately handling "small" and "large" items, and applying a dynamic programming-based algorithm to the small items, or rather to a relaxed version of the problem where not all of the items need to be packed in bins.

Future work may address the remaining gap between approximations ratio lower and upper bounds for ΩRBP and ΓRBP with different values of Γ (in particular the cases considered in this paper of Γ = 1, 2 and general Γ > 2). While is is shown that ΓRBP may not have a polynomial asymptotic approximation scheme, it remains open whether the same may be true for ΩRBP. A different line of research could seek to extend our results to natural extensions of the uncertainty sets considered herein. For instance, an important variant of these sets would consider smoothing constraints, correlating the amount of deviation of subsets of dependent items. In particular, results obtained for the U Ω model seem to naturally extend to this smoothed variant. Applications may include, for example, medical procedures that last longer when the staff performing them is not well rested. For such medical-procedure scheduling applications, also of interest are two stage models in which some of the items can be repacked at a cost, once the uncertainty is revealed. Clearly, if the cost of repacking is higher than unity, then the single stage approximation results of this paper would apply as a special case of this general setting.

initialization: j = 1 1 2 k j=1 {b 1 j , b 2 j } Algorithm 1 :

 1211 Pack items (with smaller index first) in b j until f U (b j ) > 1 or n ∈ b j . If n / ∈ b j then j ← j + 1 and repeat Step 1. Otherwise, k ← j proceed to Step 2. Pack the last item of each bin in a new bin: for any j, let i = max(b j ), b 1 j = b j \ {i}, and b 2 j = {i} return : Next-Fit-Decreasing(I) 1.3. Structure of the paper.

Lemma 3 . 1 .

 31 Suppose that the items are ordered according to (3.1). Then k ≤ k * . Proof. Consider an optimal solution b * 1 , . . . , b * k * and the subset of optimal bins given by G

  and k ≤ k * as claimed. The lemma combined with Step 2 of Next-Fit-Decreasing immediately imply the following theorem. Theorem 3.2. If the items are ordered according to (3.1) then Next-Fit-Decreasing is a 2-approximation algorithm for ΩRBP.

Theorem 4 . 1 .

 41 Suppose that the items are ordered according to (4.1). Then Next-Fit-Decreasing is a 2Γ-approximation algorithm for ΓRBP. Recall that k is the number of bins used in Step 1 and let s = (b 1 , . . . , b k ) be the bins output at the end of Step 1. Let s * = {b * 1 , . . . , b * k * } be an optimal solution. Define i * j = max(Γ(b * j )), and also let i j = max(Γ(b j )) for each j ∈ [k ]. The key element in proving Theorem 4.1 is the following counterpart of Lemma 3.1, a result which immediately implies an approximation-factor guarantee of 2Γ. Lemma 4.2. Suppose that the items are ordered according to (4.1). Then k ≤ Γk * .

Lemma 4 . 4 .

 44 |M | ≤ k * . Proof. For convenience in the following without loss of generality let M = {b 1 , . . . , b k } where k = |M | ≤ k , and let {b * 1 , . . . , b * k * * } ⊆ s * be the set of optimal solution bins each containing Γ deviations, with bins ordered in non-increasing smallest deviation size, so âi * k ≥ âi * k+1 for each k ∈ [k * * -1], and where k * * ≤ k * . Assume for the sake of deriving a contradiction that k

  [k * * ] âΓ (b j ) ≤ j∈[k * * ] âΓ (b * j ). To do so we show by induction on k = 1, . . . , k * * -1 for a fixed instance and corresponding algorithm bins and optimal solution bins, b 1 , . . . , b k and b * 1 , . . . , b k * * , respectively. For k = 1, and all i ∈ b 1 by the fact that these items either do not deviate in s * or are in the set of smallest deviations, {i * 1 , . . . , i * k * * }, it follows that âi ≤ âi * 1 , so âΓ (b 1 ) ≤ âΓ (b * 1 ). Now assume j∈[k] âΓ (b j ) ≤ j∈[k] âΓ (b * j ) in order to prove that j∈[k+1] âΓ (b j ) ≤ j∈[k+1] âΓ (b * j ). By the induction hypothesis j∈[k]

The last inequality followed from k - 1 ≥

 1 k * * and âΓ (b j ) ≤ 1 for j ∈ [k ], implying that k -1 j=k * * +1 âΓ (b j ) ≤ k -k * * . Proof of Lemma 4.2. Lemmas 4.3 and 4.4 immediately imply that k ≤ Γk * , thus proving the claim of Lemma 4.2.

Theorem 4 . 5 .

 45 If the items are ordered according to (4.1) or (3.1), then the approximation ratio of Next-Fit-Decreasing for ΓRBP is at least 2Γ 3 . Proof. Let us define an instance where the ordering (4.1) can lead to Step 1 of Next-Fit-Decreasing using k = Γ bins while OPT = 3. Every row of the Γ × Γ matrix below corresponds to the set of items in a bin (after the Step 1) of Next-Fit-Decreasing algorithm (4.2)

2 Γ(1+c) , = 1 Γ(Γ 2 +Γ- 1 )(

 2121 consider the instance given by the following Γ × Γ matrix:( , c ) ( , c ) . . . ( , c ) , c ) ( , c ) . . . ( , c )It can be verified that Next-Fit-Decreasing opens a bin for each row, since (1 + c) + (Γ -1)(1 + c) > 1.

Corollary 4 . 6 .

 46 If the items are ordered according to (4.1) and Γ = 1 or Γ = 2, then Next-Fit-Decreasing is a 2-factor or 4-factor approximation algorithm, respectively, for ΓRBP.

Lemma 5 . 2 .

 52 For any input I of ΓRBP where k * = OPT(I), (I, k * , (Γ -1)k * ) is a yes instance of ΓRBP-T.

  Y j = b * j \ b j to the trash. Note that Y j is either the set of Γ -1 largest deviating items of b * j , or otherwise b * j = Γ(b * j ) and |Γ(b * j )| < Γ. So, (b 1 , . . . , b k * , T ) is a yes instance for the ΓRBP-T problem since evidently f (b j ) = a(b j ) + ãΓ (b j ) ≤ a(b * j ) + âΓ (b * j ) ≤ 1 and also |T | ≤ (Γ -1)k * .

2 4 I 6 if c(s) = 0 then 7 ub ← k 8 s 9 else 10 lb ← k + 1 11Algorithm 2 :Proposition 5 . 8 .

 246789101258 ← (1, Γ(k -1), 1) 5 Let s = (L, {b 1 , . . . , b k }, T ), be an optimal almost-feasible solution returned by the G-ΓRBP-T algorithm given input I * ← s (accordingly T * ← T ) Pack items in T * into bins b * k+1 , . . . , b * k according to Observation 1. return: b * 1 , . . . , b * k Algorithm for ΓRBP with small values. Algorithm 2 is a 3-approximation for ΓRBP with small values.

14

  Unpack b0 into b until f (b ) > 1 or all items of b0 are unpacked, then unpack the potentially remaining items of b0 into L.

15 s

 15 g ← (L, {b , b +1 . . . , bk }, T ∪ T ) 16 if c(s g ) < c(s) then s ← s g return: s Algorithm 3: DP(q, t, )Let us introduce notations to describe the packing of the DP. Let b X = b ∩X be the items of X packed in b by the DP, and let ∆ = 1-f (b X ) (∆ could be negative).

Fig. 1 .

 1 Fig.1. DP algorithm handling guess (q , t ), starting from item q.

In general, if we have a polynomial time additive approximation algorithm using OP T +f (OP T ) bins and polynomial time ρ-approximation algorithm for ΓRBP with small values then our algorithm uses OP T (ρ + 1) + f (OP T ) bins for ΓRBP in polynomial time.

This research has benefited from the support of the ANR project ROBUST [ANR-16-CE40-0018].


Student Member, IEEE Axel Laborieux

Marc Bocquet

Student Member, IEEE Hirtzlin Tifenn

Member, IEEE Jacques-Olivier Klein

Etienne Nowak

Member, IEEE Elisa Vianello

Member, IEEE Jean-Michel Portal

Member, IEEE Damien Querlioz
email: damien.querlioz@c2n.upsaclay.fr

Tifenn Hirtzlin

Implementation of Ternary Weights with Resistive RAM Using a Single Sense Operation per Synapse

Keywords: Neural Networks, Resistive Memory, Quantized Neural Networks, Low Voltage Operation, Sense Amplifier

The design of systems implementing low precision neural networks with emerging memories such as resistive random access memory (RRAM) is a significant lead for reducing the energy consumption of artificial intelligence. To achieve maximum energy efficiency in such systems, logic and memory should be integrated as tightly as possible. In this work, we focus on the case of ternary neural networks, where synaptic weights assume ternary values. We propose a two-transistor/two-resistor memory architecture employing a precharge sense amplifier, where the weight value can be extracted in a single sense operation. Based on experimental measurements on a hybrid 130 nm CMOS/RRAM chip featuring this sense amplifier, we show that this technique is particularly appropriate at low supply voltage, and that it is resilient to process, voltage, and temperature variations. We characterize the bit error rate in our scheme. We show based on neural network simulation on the CIFAR-10 image recognition task that the use of ternary neural networks significantly increases neural network performance, with regards to binary ones, which are often preferred for inference hardware. We finally evidence that the neural network is immune to the type of bit errors observed in our scheme, which can therefore be used without error correction.

I. INTRODUCTION

Artificial Intelligence has made tremendous progress in recent years due to the development of deep neural networks. Its deployment at the edge, however, is currently limited by the high power consumption of the associated algorithms [START_REF] Xu | Scaling for edge inference of deep neural networks[END_REF]. Low precision neural networks are currently emerging as a solution, as they allow the development of low power consumption hardware specialized in deep learning inference [START_REF] Hubara | Quantized neural networks: Training neural networks with low precision weights and activations[END_REF]. The most extreme case of low precision neural networks, the Binarized Neural Network (BNN), also called XNOR-NET, is receiving particular attention as it is especially efficient for hardware implementation: both synaptic weights and neuronal activations assume only binary values [START_REF] Courbariaux | Binarized neural networks: Training deep neural networks with weights and activations constrained to+ 1 or-1[END_REF], [START_REF] Rastegari | Xnor-net: Imagenet classification using binary convolutional neural networks[END_REF]. Remarkably, this type of neural network can achieve high accuracy on vision tasks [START_REF] Lin | Towards accurate binary convolutional neural network[END_REF]. One particularly investigated lead is to fabricate hardware BNNs with emerging memories such as resistive RAM or memristors [START_REF] Bocquet | In-memory and error-immune differential rram implementation of binarized deep neural networks[END_REF]- [START_REF] Lee | Adaptive precision cnn accelerator using radix-x parallel connected memristor crossbars[END_REF]. The low memory requirements of BNNs, as well as their reliance on simple arithmetic operations, make them indeed particularly adapted for "in-memory" or "near-memory" computing approaches, which achieve superior energy-efficiency by avoiding the von Neumann bottleneck entirely.

Ternary neural networks [START_REF] Alemdar | Ternary neural networks for resource-efficient ai applications[END_REF] (TNN, also called Gated XNOR-NET, or GXNOR-NET [START_REF] Deng | Gxnor-net: Training deep neural networks with ternary weights and activations without fullprecision memory under a unified discretization framework[END_REF]), which add the value 0 to synaptic weights and activations, are also considered for hardware implementations [START_REF] Ando | Brein memory: A 13-layer 4.2 k neuron/0.8 m synapse binary/ternary reconfigurable in-memory deep neural network accelerator in 65 nm cmos[END_REF]- [START_REF] Pan | Skyrmion-induced memristive magnetic tunnel junction for ternary neural network[END_REF]. They are comparatively receiving less attention than binarized neural networks, however. In this work, we highlight that implementing TNNs does not necessarily imply considerable overhead with regards to BNNs. We introduce a two-transistor/two-resistor memory architecture for TNN implementation. The array uses a precharge sense amplifier for reading weights, and the ternary weight value can be extracted in a single sense operation, by exploiting the fact that latency of the sense amplifier depends on the resistive states of the memory devices. This work extends a hardware developed for the energy-efficient implementation of BNNs [START_REF] Bocquet | In-memory and error-immune differential rram implementation of binarized deep neural networks[END_REF], where the synaptic weights are implemented in a differential fashion. We, therefore, show that it can be extended to TNNs without overhead on the memory array.

The contribution of this work is as follows. After presenting the background of the work (section II):

• We demonstrate experimentally, on a fabricated 130 nm RRAM/CMOS hybrid chip, a strategy for implementing ternary weights using a precharge sense amplifier, which is particularly appropriate when the sense amplifier is operated at low supply voltage (section III).

• We analyze the bit errors of this scheme experimentally and their dependence on the RRAM programming conditions (section V). • We verify the robustness of the approach to process, voltage, and temperature variations (section IV). • We carry simulations that show the superiority of TNNs over BNNs on the canonical CIFAR-10 vision task, and evidence the error resilience of hardware TNNs (section VI). • We discuss the results, and compare our approach with the idea of storing three resistance levels per device. Partial and preliminary results of this work have been presented at a conference [START_REF] Laborieux | Low power in-memory implementation of ternary neural networks with resistive ram-based synapse[END_REF]. This journal version adds the experimental characterization of bit errors in our architecture, supported by a comprehensive analysis of the impact on process, voltage, and temperature variations, and their impact at the neural network level, together with a detailed analysis of the use of ternary networks over binarized ones.

II. BACKGROUND

The main equation in conventional neural networks is the computation of the neuronal activation A j = f ( i W ji X i ) , where A j , the synaptic weights W ji , and input neuronal activations X i assume real values, and f is a non-linear activation function. Binarized neural networks (BNNs) are a considerable simplification of conventional neural networks, in which all neuronal activations (A j , X i ) and synaptic weights W ji can only take binary values meaning +1 and -1. Neuronal activation then becomes:

A j = sign i XN OR (W ji , X i ) -T j , (1)
where sign is the sign function, T j is a threshold associated with the neuron, and the XN OR operation is defined in Table I. Training BNNs is a relatively sophisticated operation, during which each synapse needs to be associated with a real value in addition to its binary value (see Appendix).

Once training is finished, these real values can be discarded, and the neural network is entirely binarized. Due to their reduced memory requirements, and reliance on simple arithmetic operations, BNNs are especially appropriate for in-or near-memory implementations. In particular, multiple groups investigate the implementation of BNN inference with resistive memory tightly integrated at the core of CMOS [START_REF] Bocquet | In-memory and error-immune differential rram implementation of binarized deep neural networks[END_REF]- [START_REF] Lee | Adaptive precision cnn accelerator using radix-x parallel connected memristor crossbars[END_REF].

Usually, resistive memory stores the synaptic weights W ji . However, this comes with a significant challenge: resistive memory is prone to bit errors, and in digital applications, is typically used with strong error-correcting codes (ECC). ECC, which requires large decoding circuits [START_REF] Gregori | On-chip error correcting techniques for new-generation flash memories[END_REF], goes against the principles of in-or near-memory computing. For this reason, [START_REF] Bocquet | In-memory and error-immune differential rram implementation of binarized deep neural networks[END_REF] proposes a two-transistor/two-resistor (2T2R) structure, which reduces resistive memory bit errors, without the need for ECC decoding circuit, by storing synaptic weights in a differential fashion. This architecture allows the extremely efficient implementation of BNNs, and using the resistive memory devices in very favorable programming conditions (low energy, high endurance). It should be noted that systems using this architecture function with row-by-row read operations, and do not use the in-memory computing technique of using the Kirchhoff current law to perform the sum operation of neural networks, while reading all devices at the same time [START_REF] Prezioso | Training and operation of an integrated neuromorphic network based on metal-oxide memristors[END_REF], [START_REF] Ambrogio | Equivalent-accuracy accelerated neural-network training using analogue memory[END_REF]. This choice limits the parallelism of such architectures, while at the same time avoiding the need of analog-to-digital conversion and analog circuits such as operational amplifiers, as discussed in detail in [START_REF] Hirtzlin | Digital biologically plausible implementation of binarized neural networks with differential hafnium oxide resistive memory arrays[END_REF].

In this work, we show that the same architecture can be used for a generalization of BNNs -ternary neural networks (TNNs) 1 , where neuronal activations and synaptic weights A j , X i , and W ji can now assume three values: +1, -1, and 0. Equation (1) now becomes:

A j = φ i GXN OR (W ji , X i ) -T j . ( 2
)
GXN OR is the "gated" XNOR operation that realizes the product between numbers with values +1, -1 and 0 (Table I).

φ is an activation function that outputs +1 if its input is greater than a threshold ∆, -1 if the input is lesser than -∆ and 0 otherwise. We show experimentally and by circuit simulation in sec. III how the 2T2R BNN architecture can be extended to TNNs with practically no overhead, in sec. V its bit errors, and in sec. VI the corresponding benefits in terms of neural network accuracy. , if the two devices are programmed in an (a) LRS / HRS (5kΩ/350kΩ) or (b) HRS/HRS (320kΩ/350kΩ) configuration.

TABLE I TRUTH TABLES OF THE XNOR AND GXNOR GATES

W ji X i XN OR -1 -1 1 -1 1 -1 1 -1 -1 1 1 1 W ji X i GXN OR -1 -1 1 -1 1 -1 1 -1 -1 1 1 1 0 X 0 X 0
In this work, we use the architecture of [START_REF] Bocquet | In-memory and error-immune differential rram implementation of binarized deep neural networks[END_REF], where synaptic weights are stored in a differential fashion. Each bit is implemented using two devices programmed either as low resistance state (LRS) / high resistance state (HRS) to mean weight +1 or HRS/LRS to mean weight -1. Fig. 1 presents the test chip used for the experiments. This chip cointegrates 130nm CMOS and resistive memory in the back-end-of-line, between levels four and five of metal. The resistive memory cells are based on 10nm thick hafnium oxide (Fig. 1(a)). All devices are integrated with a series NMOS transistor. After an initial forming step (consisting in the application of a voltage ramp from zero volts to 3.3V at a rate of 1000V/s, and with a current limited to a compliance of 200µA), the devices can switch between high resistance state (HRS) and low resistance state (LRS), through the dissolution or creation of conductive filaments of oxygen vacancies. Programming into the HRS is obtained by the application of a negative RESET voltage pulse (typically between 1.5V and 2.5V during 1µs). Programming into the LRS is obtained by the application of a positive SET pulse (also typically between 1.5V and 2.5V during 1µs), with current limited to a compliance current through the choice of the voltage applied on the transistor gate through the word line. This test chip is designed with highly conservative sizing, allowing the application of a wide range of voltages and electrical currents to the RRAM cells. The area of each bit cell is 6.6 × 6.9µm 2 . More details on the RRAM technology are provided in [START_REF] Hirtzlin | Digital biologically plausible implementation of binarized neural networks with differential hafnium oxide resistive memory arrays[END_REF].

Our experiments are based on a 2, 048 devices array incorporating all sense and periphery circuitry, illustrated in Fig. 1(b-c). The ternary synaptic weights are read using onchip precharge sense amplifiers (PCSA), presented in Fig. 2, and initially proposed in [START_REF] Zhao | High speed, high stability and low power sensing amplifier for mtj/cmos hybrid logic circuits[END_REF] for reading spin-transfer magnetoresistive random access memory. Fig. 3(a) shows an electrical simulation of this circuit to explain its working principle, using the Mentor Graphics Eldo simulator. These first simulations are presented in the commercial 130nm ultra-low leakage technology, used in our test chip, with a low supply voltage of 1.2V [START_REF] Dreslinski | Near-threshold computing: Reclaiming moore's law through energy efficient integrated circuits[END_REF], with thick oxide transistors (the nominal voltage in this process for thick oxide transistor is 5V). Since the technology targets ultra-low leakage applications the threshold voltages are significantly high (around 0.6V), thus a supply voltage of 1.2V significantly reduces the overdrive of the transistors (V GS -V T H ).

In the first phase (SEN=0), the outputs Q and Qb are precharged to the supply voltage V DD . In the second phase (SEN=V DD ), each branch starts to discharge to the ground. The branch that has the resistive memory (BL or BLb) with the lowest electrical resistance discharges faster and causes its associated inverter to drive the output of the other inverter to the supply voltage. At the end of the process, the two outputs are therefore complementary and can be used to tell which resistive memory has the highest resistance and therefore the synaptic weight. We observed that the convergence speed of a PCSA depends heavily on the resistance state of the two resistive memories. This effect is particularly magnified when the PCSA is used with a reduced overdrive, as presented here: the operation of the sense amplifier is slowed down, with regards to nominal voltage operation, and the convergence speed differences between resistance values become more apparent. Fig. 3(b) shows a simulation where the two devices, BL and BLb, were programmed in the HRS. We see that the two outputs converge to complementary values in more than 200ns, whereas less than 50ns were necessary in Fig. 3(a), where the These first simulations suggest a technique for implementing ternary weights using the memory array of our test chip. Similarly to when this array is used to implement BNN, we propose to program the devices in the LRS/HRS configuration to mean the synaptic weight 1, and HRS/LRS to mean the synaptic weight -1. Additionally, we use the HRS/HRS configuration to mean synaptic weight 0, while the LRS/LRS configuration is avoided. The sense operation is performed during a duration of 50ns. If at the end of this period, outputs Q and Qb have differentiated, causing the output of the XOR gate to be 1, output Q determines the synaptic weight (1 or -1). Otherwise, the output of the XOR gate is 0, and the weight is determined to be 0.

This type of coding is reminiscent to the one used by the 2T2R ternary content-addressable memory (TCAM) cell of [START_REF] Yang | Ternary content-addressable memory with mos 2 transistors for massively parallel data search[END_REF], where the LRS/HRS combination is used for coding 0, the HRS/LRS combination for coding 1, and the HRS/HRS combination for coding "don't care" (or X).

Experimental measurements on our test chip confirm that the PCSA can be used in this fashion. We first focus on one synapse of the memory array. We program one of the two devices (BLb) to a resistance of 100kΩ. We then program its complementary device BL to several resistance values, and for each of them perform 100 read operations of duration 50ns, using on-chip PCSAs.

These PCSAs are fabricated using thick-oxide transistors, designed for a nominal supply voltage of 5V , and here used with a supply voltage of 1.2V , close to their threshold voltage (0.6V ), to reduce their overdrive, and thus to exacerbate the PCSA delay variations. In the test chip, they are sized conservatively with a total area of 290µm 2 . The use of thick oxide transistors in this test chip allows us to investigate the behavior of the devices at high voltages, without the concern of damaging the CMOS periphery circuits. Fig. 4 plots the probability that the sense amplifier has converged during the read time. In 50ns, the read operation is only converged if the resistance of the BL device is significantly lower than 100kΩ.

To evaluate this behavior in a broader range of programming conditions, we repeated the experiment on 109 devices and their complementary devices of the memory array programmed, every 14 times, with various resistance values in the resistive memory, and performed a read operation in 50ns with an on-chip PCSA. The memory array of our test chip features one separate PCSA per column. Therefore, 32 different PCSAs are used in our results. Fig. 5(a) shows, for each couple of resistance values R BL and R BLb if the read operation was converged with Q = V DD (blue), meaning a weight of 1, converged with Q = 0 (red), meaning a weight of -1, or not converged (grey) meaning a weight of 0.

The results confirm that LRS/HRS or HRS/LRS configurations may be used to mean weights 1 and -1, and HRS/HRS for weight 0. When both devices are in HRS (resistance higher than 100kΩ, the PCSA never converges within 50ns (weight of 0). When one device is in LRS (resistance lower than 10kΩ, the PCSA always converges within 50ns (weight of ±1). The separation between the 1 (or -1) and 0 regions is not strict, and for intermediate resistance values, we see that the read operation may or may not converge in 50ns. Fig. 5(b) summarizes the different operation regimes of the PCSA.

(±1 → 0) (0 → ±1) Fig. 7(a)
< 10 -6 < 1% 6.5% Fig. 7(b) < 10 -6 < 1% 18.5%

IV. IMPACT OF PROCESS, VOLTAGE, AND TEMPERATURE VARIATIONS

We now verify the robustness of the proposed scheme to process, voltage, and temperature variation. For this purpose, we performed extensive circuit simulations of the operation of the sense amplifier, reproducing the conditions of the experiments of Fig. 5, using the same resistance values for the RRAM devices, and including process, voltage, and temperature variations. The results of the simulations are processed and plotted using the same format as the experimental results of Fig. 5, to ease comparison.

These simulations are obtained using the Monte Carlo simulator provided by the Mentor Graphics Eldo tool with parameters validated on silicon, provided by the design kit of our commercial CMOS process. Each point in the graphs of Fig. 6 therefore features different transistor parameters. We included global and local process variations, as well as transistor mismatch, in order to capture the whole range of transistor variabilities observed in silicon. In order to assess the impact of voltage and temperature variations, these simulations are presented in three conditions: slow transistors (0 • C temperature, and 1.1V supply voltage, Fig. 6(a)), experimental conditions (27 • C temperature, and 1.2V supply voltage, Fig. 6(b)), and fast transistors (60 • C temperature, and 1.3V supply voltage, Fig. 6(c)). The RRAM devices are modeled by resistors. Their process variations are naturally included through the use of different resistance values in Fig. 6. The impact of voltage variation on RRAM is naturally included through Ohm's law, and the impact of temperature variation, which is smaller than on transistors, is neglected.

In all three conditions, the simulation results appear very similar to the experiments. Three clear regions are observed: non-convergence of the sense amplifier within 50ns for devices in HRS/HRS, and convergence within this time to a +1 or -1 value for devices in LRS/HRS and HRS/LRS, respectively. However, the frontier between these regimes is much sharper in the simulations than in the experiments. As the different data points in Fig. 6 differ by process and mismatch variations, this suggests that process variation does not cause the stochasticity observed in the experiments of Fig. 5, and that they have little impact in our scheme.

We also see that the frontier between the different sense regimes in all three operating conditions remains firmly within the 10 -100kΩ range, suggesting that even high variations of voltage (±0.1V ) and temperature (±30 • C) do not endanger the functionality of our scheme. Logically, in the case of fast transistors, the frontier is shifted toward higher resistances, whereas in the case of slow transistors, it is shifted toward lower resistances. Independent simulations allowed verifying that this change is mostly due to the voltage variations: the temperature variations have an almost negligible impact on the proposed scheme.

We also observed that the impact of voltage variations increased importantly when reducing the supply voltage. For example, with a supply voltage of 0.7V instead of the 1.2V value considered here, variations of the supply voltage of ±0.1V can impact the mean switching delay of the PCSA, by a factor two. The thick oxide transistors used in this work have a nominal voltage of 5V , and a typical threshold voltage of approximately 0.6V . Therefore, although our scheme is especially appropriate for supply voltages far below the nominal voltage, it is not necessarily appropriate for voltages in the subthreshold regime, or very close to the threshold voltage.

V. PROGRAMMABILITY OF TERNARY WEIGHTS

To ensure reliable functioning of the ternary sense operation, we have seen that devices in LRS should be programmed to electrical resistance below 10kΩ, and devices in HRS to resistances above 100kΩ (Fig. 5(b)). The electrical resistance of resistive memory devices depends considerably on their programming conditions [START_REF] Hirtzlin | Digital biologically plausible implementation of binarized neural networks with differential hafnium oxide resistive memory arrays[END_REF], [START_REF] Grossi | Fundamental variability limits of filament-based rram[END_REF]. Fig. 7 shows the distributions of LRS and HRS resistances using two programming conditions, over the 2, 048 devices of the array, differentiating devices connected to bit lines and to bit lines bar. We see that in all cases, the LRS features a tight distribution. The SET process is indeed controlled by a compliance current that naturally stops the filament growth at a targeted resistance value [START_REF] Bocquet | Robust compact model for bipolar oxidebased resistive switching memories[END_REF]. An appropriate choice of the compliance current can ensure LRS below 10kΩ in most situations.

On the other hand, the HRS shows a broad statistical distribution. In the RESET process, the filament indeed breaks in a random process, making it extremely hard to control the final state [START_REF] Bocquet | Robust compact model for bipolar oxidebased resistive switching memories[END_REF], [START_REF] Ly | Role of synaptic variability in resistive memorybased spiking neural networks with unsupervised learning[END_REF]. The use of stronger programming conditions leads to higher values of the HRS.

This asymmetry between the variability of LRS and HRS means that in our scheme, the different ternary weight values feature different error rates naturally. The ternary error rates in the two programming conditions of Fig. 7(a) are listed in Table II. Errors of Type 1, where weight values of 1 and -1 are inverted are the least frequent. Errors of Type 2, where a weight value of 1 or -1 is replaced by a weight value of 0 are infrequent as well. On the other hand, due to the large variability of the HRS, weight values 0 have a significant probability to be measured as 1 or -1 (Type 3 errors): 6.5% in the conditions of Fig. 7(a), and 18.5% in the conditions of Fig. 7(b). Some resistive memory technologies with large memory windows, such as specifically optimized conductive bridge memories [START_REF] Vianello | Resistive memories for ultra-low-power embedded computing design[END_REF], would feature lower Type 3 error rates. Similarly, program-and-verify strategies [START_REF] Lee | Multi-level switching of triple-layered taox rram with excellent reliability for storage class memory[END_REF]- [START_REF] Xu | Understanding the trade-offs in multi-level cell reram memory design[END_REF] may reduce this error rate. Nevertheless, the higher error rate for zeros than for 1 and -1 weights is an inherent feature of our architecture. Therefore, in the next section, we assess the impact of these errors on the accuracy of neural networks.

VI. NETWORK-LEVEL IMPLICATIONS

We first investigate the accuracy gain when using ternarized instead of binarized networks. We trained BNN and TNN versions of networks with Visual Geometry Group (VGG) type architectures [START_REF] Simonyan | Very deep convolutional networks for large-scale image recognition[END_REF] on the CIFAR-10 task of image recognition, consisting in classifying 1,024 pixels color images among ten classes (airplane, automobile, bird, cat, deer, dog, frog, horse, ship, and truck) [START_REF] Krizhevsky | Learning multiple layers of features from tiny images[END_REF]. Simulations are performed using PyTorch 1.1.0 [START_REF] Paszke | Automatic differentiation in pytorch[END_REF] on a cluster of eight Nvidia GeForce RTX 2080 GPUs.

The architecture of our networks consists of six convolutional layers with kernel size three. The number of filters at the first layer is called N and is multiplied by two every two layers. Maximum-value pooling with kernel size two is used every two layers and batch-normalization [START_REF] Ioffe | Batch normalization: Accelerating deep network training by reducing internal covariate shift[END_REF] every layer. The classifier consists of one hidden layer of 512 units. For the TNN, the activation function has a threshold ∆ = 5 • 10 -2 (as defined in section II). The training methods for both the BNN and the TNN are described in the Appendix. The training is performed using the AdamW optimizer [START_REF] Kingma | Adam: A method for stochastic optimization[END_REF], [START_REF] Loshchilov | Fixing weight decay regularization in adam[END_REF], with minibatch size 128. The initial learning rate is set to 0.01, and the learning rate schedule from [START_REF] Loshchilov | Fixing weight decay regularization in adam[END_REF], [START_REF]Sgdr: Stochastic gradient descent with warm restarts[END_REF] (Cosine annealing with two restarts, for respectively 100, 200, 400 epochs) is used, resulting in a total of 700 epochs. Training data is augmented using random horizontal flip, and random choice between cropping after padding and random small rotations.

No error is added during the training procedure, as our device is meant to be used for inference. The synaptic weights encoded by device pairs would be set after the model has been trained on a computer. Fig. 8 shows the maximum test accuracy resulting from these training simulations, for different sizes of the model. The error bars represent one standard deviation of the training accuracies. TNNs always outperform BNNs with the same model size (and, therefore, the same number of synapses). The most substantial difference is seen for smaller model size, but a significant gap remains even for large models. Besides, the difference in the number of parameters required to reach a given accuracy for TNNs and BNNs increases with higher accuracies. There is, therefore, a definite advantage to use TNNs instead of BNNs.

Fig. 8 compared fully ternarized (weights and activations) with regards to fully binarized (weights and activations) ones. Table III lists the impact of weight ternarization for different types of activations (binary, ternary, and real activation). All results are reported on a model of size N = 128, trained on CIFAR-10, and are averaged over five training procedures. We observe that for BNNs and TNNs with quantized acti- vations, the accuracy gains provided by ternary weights over binary weights are 0.84 and 0.86 points and are statistically significant over the standard deviations. This accuracy gain is more important than the gain provided by ternary activations over binary activations, which is about 0.3 points. This bigger impact of weight ternarization over ternary activation may come from the ternary kernels having a better expressing power over binary kernels, which are often redundant in practical settings [START_REF] Courbariaux | Binarized neural networks: Training deep neural networks with weights and activations constrained to+ 1 or-1[END_REF]. The gain of ternary weights drops to 0.26 points if real activation is allowed (using rectified linear unit, or ReLU, as activation function, see appendix), and is not statistically significant considering the standard deviations.

Quantized activations are vastly more favorable in the context of hardware implementations, and in this situation, there is thus a statistically significant benefit provided by ternary weights over binary weights. We finally investigate the impact of bit errors in BNNs and TNNs to see if the advantage provided by using TNNs in our approach remains constant when errors are taken into account. Consistently with the results reported in section V, three types of errors are investigated: Type 1 errors are sign switches, e.g., +1 mistaken for -1, Type 2 errors are only defined for TNNs and correspond to ±1 mistaken for 0, and Type 3 errors are 0 mistaken for ±1, as illustrated in the inset schematic of Fig. 9(a).

Fig. 9(a) shows the impact of these errors on the test accuracy for different values of the error rate at inference time. These simulation results are presented on CIFAR-10 with a model size of N = 128. Errors are randomly and artificially introduced in the weights of the neural network . Bit errors are included at the layer level and sampled at each mini-batch of the test set. Type 1 errors switch the sign of a synaptic weight with a probability equal to the rate of type 1 errors. Type 2 errors set a non-zero synaptic weight to 0 with a probability equal to the type 2 error rate. Type 3 errors set a synaptic weight of 0 to ±1 with a probability equal to the type 3 error rate, the choice of the sign (+1 or -1) is made with 0.5 probability. Fig. 9 is obtained by averaging the test accuracy obtained for five passes through the test set for increasing bit error rate.

Type 1 errors have the most impact on neural network accuracy. As seen in Fig. 9(b), the impact of these errors is similar to the impact of weight errors in a BNN. On the other hand, Type 3 errors have the least impact, with bit error rates as high as 20% degrading surprisingly little the accuracy. This result is fortunate, as we have seen in section V that Type 3 errors are the most frequent in our architecture.

We also performed simulations considering all three types of error at the same time, with error rates reported in Table II corresponding to the programming conditions of Fig. 7(a) and 7(b). For Type 1 and Type 2 errors, we considered the upper limits listed in Table II. For the conditions of Fig. 7(a) (Type 3 error rate of 6.5%), the test accuracy was degraded from 92.2% to 92.05 ± 0.14%, and to 92.02 ± 0.17% for the conditions of 7(b) (Type 3 error rate of 18.5%), where the average and standard deviation is performed over 100 passes through the test set. We found that the slight degradation on CIFAR-10 test accuracy was mostly due to the Type 2 errors, although Type 3 errors are much more frequent.

The fact that mistaking a 0 weight for a ±1 weight (Type 3 error) has much less impact than mistaking a ±1 weight for a 0 weight (Type 2 error) can seem surprising. However, it is known, theoretically and practically, that in BNNs, some weights have little importance to the accuracy of the neural networks [START_REF] Laborieux | Synaptic metaplasticity in binarized neural networks[END_REF]. They typically correspond to synapses that feature a 0 weight in a TNN, whereas synapses with ±1 weights in a TNN correspond to "important" synapses of a BNN. It is thus understandable that errors on such synapses have more impact on the final accuracy of the neural network.

VII. COMPARISON WITH THREE-LEVEL PROGRAMMING

An alternative approach to implementing ternary weights with resistive memory can be to program the individual devices into there separate levels. This idea is feasible, as the resistance level of the LRS can to a large extent be controlled through the choice of the compliance current during the SET operation in many resistive memory technologies [START_REF] Hirtzlin | Digital biologically plausible implementation of binarized neural networks with differential hafnium oxide resistive memory arrays[END_REF], [START_REF] Bocquet | Robust compact model for bipolar oxidebased resistive switching memories[END_REF].

The obvious advantage of this approach is that it requires a single device per synapse. This idea also brings several challenges. First, the sense operation has to be more complex. The most natural technique is to perform two sense operations, comparing the resistance of a device under test to two different thresholds. Second, this technique is much more prone to bit errors than our technique, as states are not programmed in a differential fashion [START_REF] Hirtzlin | Digital biologically plausible implementation of binarized neural networks with differential hafnium oxide resistive memory arrays[END_REF]. Additionally, this approach does not feature the natural resilience to Type 1 and Type 2 errors, and Type 2 and Type 3 errors will typically feature similar rates. Finally, unlike ours, this approach is prone to resistive drift, inherent to some resistive memory technologies [START_REF] Li | Resistance drift in phase change memory[END_REF].

These comments suggest that the choice of a technique for storing ternary weights should be dictated by technology. Our technique is especially appropriate for resistive memories not supporting single-device multilevel storage, with high error rates, or resistance drift. The three-levels per devices approach would be the most appropriate with devices with well controlled analog storage properties.

VIII. CONCLUSION

In this work, we revisited a differential memory architecture for BNNs. We showed experimentally on a hybrid CMOS/RRAM chip that, its sense amplifier can differentiate not only the LRS/HRS and HRS/LRS states, but also the HRS/HRS states in a single sense operation. This feature allows the architecture to store ternary weights, and to provide a building block for ternary neural networks. We showed by neural network simulation on the CIFAR-10 task the benefits of using ternary instead of binary networks, and the high resilience of TNNs to weights errors, as the type of errors observed experimentally in our scheme is also the type of errors to which TNNs are the most immune. This resilience allows the use of our architecture without relying on any formal error correction. Our approach also appears resilient to process, voltage, and temperature variation if the supply voltage remains reasonably higher than the threshold voltage of the transistors.

As this behavior of the sense amplifier is exacerbated at supply voltages below the nominal voltage, our approach especially targets extremely energy-conscious applications such as uses within wireless sensors or medical applications. This work opens the way for increasing the edge intelligence in such contexts, and also highlights that the low voltage operation of circuits may sometimes provide opportunities for new functionalities.

APPENDIX: TRAINING ALGORITHM OF BINARIZED AND TERNARY NEURAL NETWORKS

During the training of BNNs and TNNs, each quantized (binary or ternary) weight is associated with a real hidden weight. This approach to training quantized neural network was introduced in [START_REF] Courbariaux | Binarized neural networks: Training deep neural networks with weights and activations constrained to+ 1 or-1[END_REF] and is presented in Algorithm 1.

The quantized weights are used for computing neuron values (equations ( 1) and ( 2)), as well as the gradients values in the backward pass. However, training steps are achieved by updating the real hidden weights. The quantized weight is then determined by applying to the real value the quantizing function Quantize, which is φ for ternary or sign for binary as defined in section II. The quantization of activations is done by applying the same function Quantize, except for real activation, which is done by applying a rectified linear unit (ReLU(x) = max(0, x)).

Quantized activation functions (φ or sign) have zero derivatives almost everywhere, which is an issue for backpropagating the error gradients through the network. A way around this issue is the use of a straight-through estimator [START_REF] Bengio | Estimating or propagating gradients through stochastic neurons for conditional computation[END_REF], which consists in taking the derivative of another function instead of the almost everywhere zero derivatives. Throughout this work, we take the derivative of Hardtanh, which is 1 between -1 and 1 and 0 elsewhere, both for binary and ternary activations.

The simulation code used in this work is available publicly in the Github repository: https://github.com/Laborieux-Axel/ Quantized VGG Algorithm 1 Training procedure for binary and ternary neural networks. W h are the hidden weights, θ BN = (γ l , β l ) are Batch Normalization parameters, U W and U θ are the parameter updates prescribed by the Adam algorithm [START_REF] Kingma | Adam: A method for stochastic optimization[END_REF], (X, y) is a batch of labelled training data, and η is the learning rate. "cache" denotes all the intermediate layers computations needed to be stored for the backward pass. Quantize is either φ or sign as defined in section II. " • " denotes the elementwise product of two tensors with compatible shapes. Input: W h , θ BN = (γ l , β l ), U W , U θ , (X, y), η. Output: W h , θ BN , U W , U θ .

1: W Q ← Quantize(W h )
Computing quantized weights 2: A 0 ← X Input is not quantized 3: for l = 1 to L do For loop over the layers 4:

z l ← W Q l A l Matrix multiplication 5:
A l ← γ l • z l -E(z l ) √

Var(z l )+

+ β l

Batch Normalization [START_REF] Ioffe | Batch normalization: Accelerating deep network training by reducing internal covariate shift[END_REF] 6:

if l < L then If not the last layer 7:

A l ← Quantize(A l ) Activation is quantized

Fig. 1 .

 1 Fig. 1. (a) Electron microscopy image of a hafnium oxide resistive memory cell (RRAM) integrated in the backend-of-line of a 130nm CMOS process. (b) Photograph and (c) simplified schematic of the hybrid CMOS/RRAM test chip characterized in this work. The white rectangle in (b) materializes a single PCSA.

0 IIIFig. 2 .Fig. 3 .

 023 Fig. 2. Schematic of the precharge sense amplifier fabricated in the test chip.

Fig. 4 .

 4 Fig. 4. Two devices have been programmed in four distinct programming conditions, presented in (a), and measured using an on-chip sense amplifier. (b) Proportion of read operations that have converged in 50ns, over 100 trials.

Fig. 5 .

 5 Fig. 5. For 109 device pairs programmed with multiple R BL /R BLb configuration, value of the synaptic weight measured by the on-chip sense amplifier using the strategy described in body text and 50ns reading time.

Fig. 6 .Fig. 7 .

 67 Fig. 6. Three Monte Carlo SPICE-based simulation of the experiments of Fig. 5, in three situations: (a) slow transistors (0 • C temperature, 1.1V supply voltage), (b) experimental conditions (27 • C temperature, 1.2V supply voltage), (c) fast transistors (60 • C temperature, 1.3V supply voltage). The simulations include local and global process variations, as well as transistor mismatch, in a way that each point in the Figure is obtained using different transistor parameters. All results are plotted in the same manner and with the same conventions as Fig. 5.

Fig. 8 .

 8 Fig. 8. Simulation of the maximum test accuracy reached during one training procedure, averaged over five trials, for BNNs and TNNs with various model sizes on the CIFAR-10 dataset. Error bar is one standard deviation.

Fig. 9 .

 9 Fig. 9. Simulation of the impact of Bit Error Rate on the test accuracy at inference time for model size N = 128 TNN in (a) and BNN in (b). Type 1 errors are sign switches (e.g. +1 mistaken for -1), Type 2 errors are ±1 mistaken for 0, and Type 3 errors are 0 mistaken for ±1, as described in the inset schematics. Errors are sampled at each mini batch and the test accuracy is averaged over five passes through the test set. Error bars are one standard deviation. The bit error rate is given as an absolute rate.

  end for 10: ŷ ← A L 11: C ← Cost(ŷ, y) Compute mean loss over the batch 12:(∂ W C, ∂ θ C) ← Backward(C, ŷ, W Q , θ BN , cache) Cost gradients 13: (U W , U θ ) ← Adam(∂ W C, ∂ θ C, U W , U θ ) 14: W h ← W h -ηU W 15: θ BN ← θ BN -ηU θ 16: return W h , θ BN , U W , U θ

TABLE III COMPARISON

 III OF THE GAIN IN TEST ACCURACY FOR A N = 128 MODEL SIZE ON CIFAR-10 OBTAINED BY WEIGHT TERNARIZATION INSTEAD OF BINARIZATION FOR THREE TYPES OF ACTIVATION QUANTIZATION.

			Activations
		Binary	Ternary	Full Precision
	Weights
	Binary	91.19 ± 0.08	91.51 ± 0.09	93.87 ± 0.19
	Ternary	92.03 ± 0.12	92.35 ± 0.05	94.13 ± 0.10
	Gain of ternarization	0.84	0.86	0.26

In the literature, the name "Ternary Neural Networks" is sometimes also used to refer to neural networks where the synaptic weights are ternarized, but the neuronal activations remain real or integer[START_REF] Mellempudi | Ternary neural networks with fine-grained quantization[END_REF],[START_REF] Nurvitadhi | Can fpgas beat gpus in accelerating next-generation deep neural networks?[END_REF].

ACKNOWLEDGMENTS

The authors would like to thank M. Ernoult and L. Herrera Diez for fruitful discussions.

This work was supported by the ERC Grant NANOINFER (715872) and the ANR grant NEURONIC (ANR-18-CE24-0009).