Shilpa Pendyala
email: spendya2@mail.usf.edu

Srinivas Katkoori
email: katkoori@mail.usf.edu

Interval Arithmetic and Self Similarity Based Subthreshold Leakage Optimization in RTL Datapaths

Keywords: Sub-threshold Leakage Optimization, Minimum Leakage Input Vector, Interval Arithmetic, Self Similarity, RTL Datapath Optimization

We propose top-down and bottom-up interval propagation techniques for identifying low leakage input vectors at primary inputs of an RTL datapath. Empirically, we observed self-similarity in the leakage distribution of adder/multiplier modules i.e., leakage distribution at the sub-space level is similar to that at the entire input space. We exploit this property to quickly search low leakage vectors. The proposed module library leakage characterization is scalable and is demonstrated on adders/multipliers. Given an RTL datapath, interval propagation is carried out with the low leakage intervals of the module instances with primary inputs. The reduced interval set is further processed with simulated annealing, to arrive at the best low leakage vector set at the primary inputs. Experimental results for various DSP filters simulated in 16nm CMOS technology with top-down and bottom-up approaches yield leakage savings of 93.6% and 89.2% respectively with no area, timing, or control overheads.

Introduction and Motivation

Excessive subthreshold leakage power consumption is a serious concern in deepsubmicron technology nodes [START_REF]SIA: International technology roadmap for semiconductors (itrs)[END_REF]. Input vector control technique is widely used for subthreshold leakage optimization due to its low latency overhead. However, determination of minimum leakage vector (MLV) for large circuits is a difficult problem by itself as it requires excessive simulation time [START_REF] Abdollahi | Leakage current reduction in cmos vlsi circuits by input vector control. Very Large Scale Integration[END_REF]. An MLV is an input vector that puts the circuit in lowest leakage state possible i.e., it is an optimal vector. We refer to sub-optimal MLVs as low leakage vectors (LLVs) that put the circuit with leakage close to the optimal leakage.

For a given datapath, applying MLV to each RTL module incurs significant area and control overhead. The additional area stems from input multiplexers needed to apply the minimum leakage vector. Control overhead is incurred from the select lines of these multiplexers. As data values propagate through RTL modules, we have proposed [START_REF] Pendyala | Interval arithmetic based input vector control for rtl subthreshold leakage minimization[END_REF] to identify a set of input vectors such that they not only put the module instances at primary inputs (PIs) into low leakage, but also result in low leakage input vectors at internal module instances. In [START_REF] Pendyala | Self similarity and interval arithmetic based leakage optimization in rtl datapaths[END_REF] (conference version of this chapter), we leverage self-similarity of modulelevel leakage distributions. In this book chapter, we further extend the idea by experimenting with bottom up interval propagation. We experimented with five DSP filters and obtained average leakage savings of 93% with top down approach and 89.2% with bottom up approach. In both cases, we did not require any internal control points for any design, thus the proposed approaches incur no overhead in terms of area, control, or delay.

The overall optimization flow is as follows: First, we characterize the module library to gather low leakage vector intervals with two-phase Monte-Carlo simulation. Next, we propagate interval sets of module instances at primary inputs through the datapath. We have two choices: top-down or bottom-up. In top-down interval sets at PIs are propagated to primary outputs (POs) and in case of bottom-up, interval sets at POs are propagated back to PIs. The interval propagation yields a reduced set of low leakage intervals at the PIs. A simulated annealing algorithm is devised to find the best input vector set from the reduced interval set. Next, we briefly describe the self-similarity based scalable module leakage characterization for bit-sliced adder/multipliers.

As the search space grows exponentially with the module bit-width, finding low leakage vectors by exhaustive simulation is infeasible. In order to make the problem tractable, we exploit the self-similarity property of leakage distribution of a given functional unit. Briefly, given a stochastic distribution, it is said to be self-similar [START_REF] Mandelbrot | Fractal Geometry of Nature[END_REF][START_REF] Barnsley | Fractals Everywhere[END_REF], if any arbitrary sub-distribution (built by choosing contiguous samples) is similar to the original distribution. Hurst parameter is a commonly used metric to determine self-similarity of a distribution. For two module types, adder and multiplier, we empirically observe that their leakage distributions are self-similar. Given a confidence level (α) and an error tolerance (β), Halter and Najm [START_REF] Halter | A gate-level leakage power reduction method for ultra-lowpower cmos circuits[END_REF] derived a formula to determine the minimum number of random vectors needed to find a LLV that is no worse than β% of vectors with α% confidence, provided the leakage distribution is normal. While the goal of any MLV heuristic is to find one vector, our goal is to find as many LLVs as possible so that we can expand them into LLV data intervals. For adder and multiplier, we also empirically observe that their leakage distributions are normal. This normal leakage distribution along with self-similarity of a module enables us to partition the input space into small sub-spaces and then randomly sample each subspace with fixed number of vectors.

The chapter organization as follows: Sec. 2 presents background, related work, and terminology. Sec. 3 describes the module library characterization. Sec. 4 presents the proposed top-down and bottom-up interval propagation approaches followed by simulated annealing. Sec. 5 reports experimental results. Finally, Sec. 6 draws conclusions.

Background, Related Work, and Terminology

We first review input vector control techniques from the literature. We then provide a brief overview of fractal theory and self-similarity. Lastly, we present an overview of interval arithmetic.

Input Vector Control Techniques

As the leakage depends only on the current input vector, during the idle mode, we can apply the minimum leakage vector (MLV). Thus, MLV needs to be determined a priori and incorporated into the circuit. This technique is known as the Input Vector Control (IVC). For an n-input module, as the input space grows exponentially (2 n ), MLV determination heuristics have been proposed [START_REF] Abdollahi | Leakage current reduction in cmos vlsi circuits by input vector control. Very Large Scale Integration[END_REF][START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF]. For an IVC technique, the area penalty occurs due to additional hardware needed to incorporate MLV into the circuit. Delay penalty is incurred if this additional hardware is in the critical path of the circuit.

To the best of our knowledge, all the proposed IVC techniques are at the logic level. Abdollahi, Fallah, and Pedram [START_REF] Abdollahi | Leakage current reduction in cmos vlsi circuits by input vector control. Very Large Scale Integration[END_REF] propose gate-level leakage reduction with two techniques. The first technique is an input vector control wherein SAT based formulation is employed to find the minimum leakage vector. The second technique involves adding nMOS and pMOS transistors to the gate in order to increase the controllability of the internal signals. The additional transistors increase the stacking effect leading to leakage current reduction. The authors report over 70% leakage reduction with up to 15% delay penalty. Gao and Hayes [START_REF] Gao | Exact and heuristic approaches to input vector control for leakage power reduction[END_REF] present integer linear programming (ILP) and mixed integer linear programming (MILP) approaches, wherein MILP performs better than ILP and is thirteen times faster. Average leakage current is about 25% larger than minimum leakage current. IVC technique does not work effectively for circuits with large logic depth. Yuan and Qu [START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF] have proposed a technique to replace the gates of worst leakage state with other gates in active mode. A divide-andconquer approach is presented that integrates gate replacement, an optimal MLV searching algorithm for tree circuits, and a genetic algorithm to connect the tree circuits. Compared with the leakage achieved by optimal MLV in small circuits, the gate replacement heuristic and the divide-and-conquer approach can reduce on average 13% and 17% leakage, respectively.

Fractals and Self Similarity

Fractals are shapes made of parts similar to the whole [START_REF] Mandelbrot | Fractal Geometry of Nature[END_REF]. Property of scaling is exhibited by fractals, which means the degree of irregularity in them tends to be identical at all scales [START_REF] Mandelbrot | Fractal Geometry of Nature[END_REF][START_REF] Barnsley | Fractals Everywhere[END_REF]. Many examples of fractal behaviors are observed in nature. Figure 1 shows a romanesco broccoli in which we can observe the self similar property of its shape. In VLSI structures, H-tree based clock signal network, employed for zero-skew, exhibits fractal behavior (Fig. 2).

The main classifications of fractals are time or space, self-similar or selfaffine, and deterministic or stochastic. Space fractals are structures exhibiting the fractal property in the space domain. While time fractals are those exhibiting in time domain. Self-similar fractals have symmetry over entire scale, which is identical to recursion, i.e., a pattern inside another pattern. The details spread across finer and finer scales with certain constant measurements. Hurst parameter (0 < H < 1), is a measure of the correlation or long range dependence in the data which leads to the fractal behavior of the data set [START_REF] Radjassamy | Faster power estimation of cmos designs using vector compaction -a fractal approach[END_REF]. A random process has H value equal to 0.5. If H is less than 0.5 then the process exhibits anti-persistence. If H value is between 0.5 and 1, then the process has long term persistence [START_REF] Mandelbrot | Fractal Geometry of Nature[END_REF]. A stochastic process can be said to exhibit fractal behavior [START_REF] Leland | Statistical analysis and stochastic modeling of self-similar data traffic[END_REF] if the H value is between 0.5 and 1. Methods to estimate Hurst parameter are described in [START_REF] Radjassamy | Faster power estimation of cmos designs using vector compaction -a fractal approach[END_REF]. In this work, we employed the R/S plot method. Photo Credit [START_REF] Sullivan | Wikimedia commons[END_REF].

Self similarity has been leveraged in estimation and optimization problems in diverse fields. We give three examples. Premarathne et al [START_REF] Premarathne | Network traffic self similarity measurements using classifier based hurst parameter estimation[END_REF] employ selfsimilarity to detect anamalous events in network traffic by showing that traffic's self-similarity property is temporarily disturbed in the event of an attack. Radjassamy and Carothers [START_REF] Radjassamy | Faster power estimation of cmos designs using vector compaction -a fractal approach[END_REF] proposed a vector compaction technique to generate a compact vector set representative of original vector set such that it mimics power-determining behavior of the latter. Such vector compaction can speed up power estimation of circuits. Qian and Rasheed [START_REF] Qian | Hurst exponent and financial market predictability[END_REF] proposed Hurst parameter based financial market analysis wherein series with high H are predicted more accurately than those with H close to 0.5. For more examples, interested reader is referred to [START_REF] Falconer | Fractal Geometry: Mathematical Foundations and Applications[END_REF].

Interval Arithmetic

Interval arithmetic (IA) [START_REF] Moore | Methods and applications of interval analysis[END_REF] is concerned with arithmetic operations such as addition and subtraction on intervals. The intervals can be either discrete or continuous. IA has been extensively applied in error bound analysis arising in numerical analysis.

In this work, we are concerned with integer arithmetic therefore we restrict our discussion to integer intervals. An interval I = [a, b] represents all integers a ≤ i ≤ b. Further, the above interval is a closed interval as it includes both extremal values. We can have an open interval, such as I = (a, b) where a < i < b. The width of an interval is the difference between the extremal values |b -a|. If the interval width is zero, then the interval is referred to as a degenerate interval (for example [START_REF]SIA: International technology roadmap for semiconductors (itrs)[END_REF][START_REF]SIA: International technology roadmap for semiconductors (itrs)[END_REF]). We can represent a given integer, say a, as a degenerate interval [a, a].

Given two intervals U = [a, b] and V = [c, d],the following equations hold:

U + V = [a + c, b + d] (1)
U -V = [a -d, b -c] (2)
U * V = [min(a * c, a * d, b * c, b * d), max(a * c, a * d, b * c, b * d)] (3)
U ÷ V = [min(a ÷ c, a ÷ d, b ÷ c, b ÷ d), max(a ÷ c, a ÷ d, b ÷ c, b ÷ d)] (4)

Notation and Problem Formulation

-A DFG, G(V, E), is a directed graph such that v i ∈ V represents an operation and e = (v i , v j ) ∈ E represents a data transfer from operation v i to v j .

-A low leakage interval set, L(t, w), is the set of all low leakage intervals of a given module of type t and w. -A low leakage output interval set, LO(t, w), is the set of corresponding outputs of all low leakage input intervals of a given module of type t and w. -P(V, t, w) is the leakage power function which calculates the total leakage of the filter. -Problem Formulation: Given the following inputs: (1) a data flow graph G(V, E); (2) set of low leakage interval sets, t,w L(t, w), and (3) t,w LO(t, w), for all distinct operations of type t and width w, we need to identify best low leakage vector on primary inputs and a set of control points C such that the objective functions, vi∈V P(V, T (v i ), W(v i )) and C, are minimized, where C is the set of control points.

Self Similarity Based Monte Carlo Characterization for Low Leakage Intervals

In this Section, we first study the leakage profiles of adder and multiplier modules. Then we introduce a Monte Carlo leakage characterization technique based on self-similarity to extract low leakage interval set of a functional unit.

Leakage Profile and Scope for Optimization

Figure 3 shows the leakage current distribution of an 8 bit ripple carry adder based on simulation with all possible (exhaustive) vectors and 1000 vectors. Similarly, Fig. 4 shows the leakage current distributions of an 8 bit parallel multiplier for exhaustive and 1000 vectors respectively. The data has been generated by Synopsys Nanosim simulations of the CMOS layouts in 16nm technology node with PTM spice models [START_REF] Cao | Asu predictive technology model website[END_REF][START_REF] Zhao | New generation of predictive technology model for sub-45nm design exploration[END_REF][START_REF] Zhao | New generation of predictive technology model for sub-45 nm early design exploration[END_REF]. Similar normal leakage distribution plots were obtained for 16 bit adder and 16 bit multiplier with 1000 random vectors as shown in Figs. 5 and 6 respectively. We observed that even with increased bit width, the leakage distribution is normal for both adder and multiplier. Several prior works in literature for eg., [START_REF] Evmorfopoulos | A monte carlo approach for maximum power estimation based on extreme value theory[END_REF] and [START_REF] Qiu | Maximum power estimation using the limiting distributions of extreme order statistics[END_REF] report similar nomal power distributions. Based on this empirical evidence, we assume that the leakage distribution of adder and multiplier modules with any bit width is normal.

Consider the exhaustive simulation of 8 bit adder and multiplier. The leakage current ranges are [0.084µA, 4.3µA] and [1.4µA, 56µA] respectively. Thus, the approximate max-to-min leakage current ratio for 8b adder and 8b multiplier are 51 and 40 respectively. For 10% tolerance (ε=0.1), the number of distinct LLVs for the adder is 119. Thus, the percentage of input space that puts the adder in a low-leakage state is (119/(2 8 × 2 8 ))×100 = 0.18%. These vectors can be merged into 80 low leakage intervals. Similarly, for multiplier, number of low leakage vectors is 490 and the size of the interval set is 329. The percentage of input space that puts the multiplier in a low-leakage state for ε=0.1, is (490/(2 8 × 2 8 ))×100 = 0.74%. Based on these numbers, we can see only a small percentage of input space can result in significant leakage reduction. Our next challenge is to locate all these low leakage intervals in the entire input space.

Self Similarity of Leakage Distributions in n bit Adders and Multipliers

We empirically observed that the sub-threshold leakage distributions of adders and multipliers exhibit self similarity property. For example, in Fig. 7 for an 8 bit adder we show the Hurst values for various sub spaces obtained by partitioning the input space in three different ways. Recall that a distribution is self similar if Hurst value is between 0.5 and 1. For each partitioning scheme, we see that the leakage distribution of each sub space is self similar to the parent distribution. Similar results have been obtained in the case of 8 bit multiplier, 16 bit adder, and 16 bit multiplier. Due to limited space, we do not include these plots. Based on these empirical results, we assume in general that the leakage distribution of a n bit adder or multiplier is self-similar. The concept of self similarity enables us to develop a scalable methodology to identify low leakage intervals for a given module.

Monte Carlo Based Low Leakage Interval Search

We propose a two-stage Monte Carlo (MC) approach to deal with large input space. Typically, an MC based approach has four steps: (a) input space determination; (b) input sampling based on a probability distribution; (c) computation of property of interest; and (d) result aggregation.

Given a confidence level (α) and error tolerance (β) and assuming a normal leakage distribution, Halter and Najm [START_REF] Halter | A gate-level leakage power reduction method for ultra-lowpower cmos circuits[END_REF] derived a formula (Eq. 5) to determine the minimum number of random vectors needed to guarantee with (100×α)% confidence that the number of vectors with leakage lower than the lowest leakage found is (100×β)%. For example, if α=0.99 and β=0.01, then 460 random vectors are sufficient to give us 99% confidence that only 1% of vectors will have better leakage than the observed.

n = log e (1 -α)/ log e (1 -β) (5)
In Secs. 3.1 and 3.2 we empirically observed that the leakage distributions of n bit adders/multipliers are normal and self-similar. Consequently, we can now partition a module instance's input space into sub-spaces and then sample each sub-space with a fixed number of vectors as determined by Eq. 5 for user specified confidence and error tolerance levels.

Given the SPICE-level model of an n bit module instance, we perform two successive MC runs. The property of interest is the leakage power.

Stage I -Coarse grained MC run: The input space under consideration is the entire space, i.e., 2 2n input vectors, which is partitioned into equal sized sub-spaces (as illustrated in Fig. 7). Let us assume α=0.99 and β=0.01. Then, we uniformly sample each sub-space to identify 460 random vectors and then simulate with the vectors. For result aggregation in this stage, from each subspace, we collect 5% of the vectors that yield low leakage for further consideration in Stage II.

Stage II -Fine grained MC run: From the leakage profiles of the functional units, we also observed that low leakage values are clustered. Hence, the sampling in this stage is biased in the neighborhood of low leakage vectors identified in previous stage. The result aggregation involves merging input vectors to create set of low leakage intervals.

Run Time Complexity. We would like to estimate the run time complexity of the characterization for a module of size n.

Stage I: The total run time of Stage I is S × K × T (n), where S is the number of sub-spaces resulting from partitioning the entire input space, K is the minimum number of vectors required for user-given confidence (α) and error tolerance (β) levels, and T (n) is the simulation time for one vector. Note that K is a constant for fixed α and β values. The user can keep the number of subspaces fixed i.e., S is a constant. T (n) is proportional to the gate complexity. In case of ripple carry adder, T (n) = O(n) as gate complexity grows linearly with bit width, while for parallel multiplier, T (n) = O(n 2 ). Therefore, the complexity of Stage I is O(n) and O(n 2 ) for adder and multiplier respectively. If the user chooses to linearly scale the number of sub-spaces with the module complexity (i.e., S = O(n)), then the run time complexity increases to O(n 2 ) and O(n 3 ) for adder and multiplier respectively.

Stage II: We perform two steps: (1) local search around the vectors found in Stage I; and (2) then merge the low leakage vectors into low leakage intervals. The run-time complexity of first step is same as that of stage I, as we sample fixed number of vectors in the neighborhood of each vector from Stage I. The worst-case run time complexity of step 2 is same as that of a two-key sorting algorithm, O(nlogn), since, we sort all input vectors and then merge immediate neighbors into intervals. Therefore, the complexity of Stage II in case of an adder is O

(n) + O(nlogn) = O(nlogn), while for multiplier it is O(n 2 ) + O(nlogn) = O(n 2 ).
Since the two stages are performed sequentially, the overall runtime complexity of MC based leakage interval characterization procedure is O(nlogn) and O(n 2 ) for n bit adder and multiplier respectively.

Proposed Approach

For low leakage interval propagation, we propose two variants: top down and bottom up approaches. The motivation to propose the two variants is to compare the amount of leakage savings possible by starting at PIs and at POs.

-In top down approach, we carry out the propagation in two iterations. In first iteration, we start with a raw set of low leakage intervals at PIs and propagate them to the POs i.e., forward propagation. The raw low leakage intervals at PIs are obtained through characterization presented in detail in Sec. 3. In second iteration, these output intervals are propagated backwards reducing the input interval set at each intermediate node and thus ending up with a sparser set of low leakage intervals at PIs. This sparser interval set is further processed with simulated annealing to identify the best LLV. -In bottom up approach, we start with raw set of low leakage output intervals at POs and propagate them all the way to the PIs where we end up with a minimized interval set. This minimized set is again processed with simulated annealing to arrive at the best LLV.

We present two motivating examples illustrating top down and bottom up approaches.

Motivating Examples

Example 1: Top Down Approach In Fig. 8, we show an example DFG with two adders (A1, A2) and one multiplier (M1). Let us say the low leakage vector sets are: L(+, 8) = {( [START_REF] Abdollahi | Leakage current reduction in cmos vlsi circuits by input vector control. Very Large Scale Integration[END_REF][START_REF] Pendyala | Self similarity and interval arithmetic based leakage optimization in rtl datapaths[END_REF], [START_REF] Barnsley | Fractals Everywhere[END_REF][START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF]), ( [START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF][START_REF] Sullivan | Wikimedia commons[END_REF], [START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF][START_REF] Sullivan | Wikimedia commons[END_REF]), ( [START_REF] Premarathne | Network traffic self similarity measurements using classifier based hurst parameter estimation[END_REF][START_REF] Zhao | New generation of predictive technology model for sub-45 nm early design exploration[END_REF], [START_REF] Premarathne | Network traffic self similarity measurements using classifier based hurst parameter estimation[END_REF]24])} and L( * , 8) = {( [START_REF] Pendyala | Interval arithmetic based input vector control for rtl subthreshold leakage minimization[END_REF][START_REF] Pendyala | Self similarity and interval arithmetic based leakage optimization in rtl datapaths[END_REF], [START_REF] Mandelbrot | Fractal Geometry of Nature[END_REF][START_REF] Barnsley | Fractals Everywhere[END_REF]), ([9, 10], [START_REF] Mandelbrot | Fractal Geometry of Nature[END_REF][START_REF] Barnsley | Fractals Everywhere[END_REF]), ( [START_REF] Eppstein | Wikimedia commons[END_REF]24], [START_REF] Mandelbrot | Fractal Geometry of Nature[END_REF][START_REF] Barnsley | Fractals Everywhere[END_REF])}.

We start applying low leakage vector interval sets for A1 and A2 i.e., { [START_REF] Abdollahi | Leakage current reduction in cmos vlsi circuits by input vector control. Very Large Scale Integration[END_REF][START_REF] Pendyala | Self similarity and interval arithmetic based leakage optimization in rtl datapaths[END_REF], [START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF][START_REF] Sullivan | Wikimedia commons[END_REF], [START_REF] Premarathne | Network traffic self similarity measurements using classifier based hurst parameter estimation[END_REF][START_REF] Zhao | New generation of predictive technology model for sub-45 nm early design exploration[END_REF]} on the first input and { [START_REF] Barnsley | Fractals Everywhere[END_REF][START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF], [START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF][START_REF] Sullivan | Wikimedia commons[END_REF], [START_REF] Premarathne | Network traffic self similarity measurements using classifier based hurst parameter estimation[END_REF]24]} on the second input. Using Eq. ( 1) of interval arithmetic, we compute A1's and A2's output range to be { [START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF][START_REF] Sullivan | Wikimedia commons[END_REF], [START_REF] Falconer | Fractal Geometry: Mathematical Foundations and Applications[END_REF]24], [28, 44]}. As we want to apply low leakage vector to M1, we need to reduce the interval sets generated by A1 and A2 to those identified as LLVs for M1. Thus, the interval set of M1 for input 1 is { [START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF][START_REF] Sullivan | Wikimedia commons[END_REF], [START_REF] Falconer | Fractal Geometry: Mathematical Foundations and Applications[END_REF]24], [28, 44]} ∩ { [START_REF] Pendyala | Interval arithmetic based input vector control for rtl subthreshold leakage minimization[END_REF][START_REF] Pendyala | Self similarity and interval arithmetic based leakage optimization in rtl datapaths[END_REF], [START_REF] Gao | Exact and heuristic approaches to input vector control for leakage power reduction[END_REF][START_REF] Radjassamy | Faster power estimation of cmos designs using vector compaction -a fractal approach[END_REF], [START_REF] Eppstein | Wikimedia commons[END_REF]24]} = {[9, 10], [START_REF] Falconer | Fractal Geometry: Mathematical Foundations and Applications[END_REF]24]} and for input 2 is { [START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF][START_REF] Sullivan | Wikimedia commons[END_REF], [START_REF] Falconer | Fractal Geometry: Mathematical Foundations and Applications[END_REF]24], [28, 44]} ∩ { [START_REF] Mandelbrot | Fractal Geometry of Nature[END_REF][START_REF] Barnsley | Fractals Everywhere[END_REF]}= ∅. In Fig. 8, for sake of clarity, we show the intervals corresponding to multiplier in bold font, while those for adders in regular font. At the second input of M1 as we obtained an empty set, we will introduce a control point to put M1 in low leakage mode. The control point consists of a multiplexer that can be used to force a low leakage vector in idle state. Generally speaking, if the interval intersection results in an empty set, we will insert a control point and start with an entire low leakage vector set from that point.

To determine the LLV at PIs, a backward propagation of minimized interval sets is implemented. Fig. 9 illustrates this step for this example. The intervals available at input 1 of M1 are fed as outputs to A1 and these intervals are propagated to the inputs of A1. This gives a minimized interval set on which simulated annealing algorithm is applied to find the best LLV.

Example 2: Bottom Up Approach Now consider the same DFG for bottom up approach, however with output low leakage vector sets: L(+, 8) = {[2, 8]} and L( * , 8) = { [START_REF] Falconer | Fractal Geometry: Mathematical Foundations and Applications[END_REF]24]}. For a given output low leakage set { [START_REF] Falconer | Fractal Geometry: Mathematical Foundations and Applications[END_REF]24]} to M1, the corresponding input leakage set is {( [START_REF] Abdollahi | Leakage current reduction in cmos vlsi circuits by input vector control. Very Large Scale Integration[END_REF][START_REF] Pendyala | Self similarity and interval arithmetic based leakage optimization in rtl datapaths[END_REF], [START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF][START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF])}. Similarly, for an adder's output set {[2, 4]}, input sets are {( [START_REF]SIA: International technology roadmap for semiconductors (itrs)[END_REF][START_REF] Abdollahi | Leakage current reduction in cmos vlsi circuits by input vector control. Very Large Scale Integration[END_REF], [START_REF] Abdollahi | Leakage current reduction in cmos vlsi circuits by input vector control. Very Large Scale Integration[END_REF][START_REF] Abdollahi | Leakage current reduction in cmos vlsi circuits by input vector control. Very Large Scale Integration[END_REF])} and { [START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF][START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF]}, the corresponding input sets are {([2, 2], [START_REF] Pendyala | Self similarity and interval arithmetic based leakage optimization in rtl datapaths[END_REF][START_REF] Pendyala | Self similarity and interval arithmetic based leakage optimization in rtl datapaths[END_REF])}. We apply output low leakage vector interval sets on primary outputs (Fig. 10). We start applying output low leakage vector set for M1 i.e., { [START_REF] Falconer | Fractal Geometry: Mathematical Foundations and Applications[END_REF]24]} at output. It propagates corresponding input vectors { [START_REF] Abdollahi | Leakage current reduction in cmos vlsi circuits by input vector control. Very Large Scale Integration[END_REF][START_REF] Pendyala | Self similarity and interval arithmetic based leakage optimization in rtl datapaths[END_REF]} to input 1 and { [START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF][START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF]} to input 2 and M1 is set to low leakage mode. The input intervals at M1 are intersected with L(+, 8) and L( * , 8). For the resulting output intervals { [START_REF] Abdollahi | Leakage current reduction in cmos vlsi circuits by input vector control. Very Large Scale Integration[END_REF][START_REF] Pendyala | Self similarity and interval arithmetic based leakage optimization in rtl datapaths[END_REF]} and { [START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF][START_REF] Yuan | A combined gate replacement and input vector control approach for leakage current reduction. Very Large Scale Integration[END_REF]} the corresponding input intervals of A1 and A2 are {( [START_REF]SIA: International technology roadmap for semiconductors (itrs)[END_REF][START_REF] Abdollahi | Leakage current reduction in cmos vlsi circuits by input vector control. Very Large Scale Integration[END_REF], [START_REF] Abdollahi | Leakage current reduction in cmos vlsi circuits by input vector control. Very Large Scale Integration[END_REF][START_REF] Abdollahi | Leakage current reduction in cmos vlsi circuits by input vector control. Very Large Scale Integration[END_REF])} and {([2, 2], [START_REF] Pendyala | Self similarity and interval arithmetic based leakage optimization in rtl datapaths[END_REF][START_REF] Pendyala | Self similarity and interval arithmetic based leakage optimization in rtl datapaths[END_REF])}, respectively. The best LLV is found by processing the reduced set with simulated annealing algorithm described in Sec. 4.2.

A2

{[1, 2]} {[2, 2]} Interval Propagation Direction A1 M1 {[8, 8]} {[2, 2]} {[4, 4]} {[16, 24]} {[8, 8]} {[2, 8]} ∩ = {[8, 8]} {[2, 4]} {[2, 8]} ∩ = {[2, 4]} {[2, 4]} {[2, 4]} {[8, 8]}

Low Leakage Vector Determination

Figure 11 shows the pseudo-code of the proposed heuristic for low leakage vector determination. It accepts an input DFG (directed acyclic graph) and low leakage vector sets for distinct types and operations obtained from the characterization procedure as described in Sec. 3. Breadth First Search is used to cover all the nodes of the graph. Flag is used to determine the direction of propagation (top down or bottom up). First, the graph is topologically sorted (line 5) to yield a sorted list L. A set C that collects the control points, is initialized (line 6). The for loop in lines 8-27 visits each node in the order specified by L. If a node is a PI node (i.e., both inputs to the node are primary), then the intervals on both inputs are intialized to the appropriate low leakage vector sets (line 13). Both inputs are added to the set C (line 14). On line 16, we call a function Interval Propagate() that accepts an ordered interval pair and the operation type of the node (i.e., T (v i )). Interval Propagate() implements the interval arithmetic equations as mentioned in Sec. 2.3 and returns an appropriate output interval I c . In line 21, we invoke Interval Intersection() that checks if the computed interval is contained in low leakage vector set of the successor v j . If the check succeeds, then we move onto to the next node in the list. If the check fails, then a new control point is inserted by resetting the inputs of node v j to its low leakage vector set (lines 22 -26) and adding the inputs of v j to control point set. On line 28, Back Propagate() further reduces interval set at primary outputs, LOreduced(t, w), by performing similar propagation in backward direction to primary inputs as shown in Fig. 12. Back Propagate() function (line 30) performs backward propagation on the complete low leakage output interval set, LO(t, w), to obtain reduced interval set at primary inputs. The algorithm is similar to find LLV(), therefore we do not elaborate in detail. On line 32 of find LLV(), we invoke Simulated Annealing Leakage() (Fig. 13) to find the best LLV from the reduced LLV set. The initial temperature value T emp is set to 100 and cooling coefficient γ to 0.99 on lines 5 and 6 respectively. The number of iterations at each temperature is equal to (500 x T emp) (line 11). An LLV is initially chosen at random (line 7) from the reduced interval set. In each iteration, a new LLV is chosen (line 12) from the neighborhood of the current LLV. The neighborhood window size is equal to the current temperature (T emp). Est Leakage() function (line 26) calculates the leakage of the DFG for any given input vector. output llv is the best LLV found by the algorithm. In the find LLV() algorithm, we assume only one successor for each node. This assumption is made to simplify the presentation of the algorithm. It is straightforward to extend the algorithm to multiple successors.

Run

Experimental Results

We report the experiment results obtained by applying the best LLV from top down and bottom up approaches on five datapath-intensive benchmarks, namely, IIR, FIR, Elliptic, Lattice, and Differential Equation Solver. Only functional unit sharing (resource optimization technique during high level synthesis) is assumed. As described in Section 3 the library is characterized a priori and the low leakage vectors sets saved. Top down and bottom up techniques process these low leakage vector sets to obtain a reduced set. The best LLV is further obtained by applying simulated annealing algorithm on the reduced set. The leakage power values are measured at the layout level using Synopsys Nanosim. We employ the Predictive Technology Models for 16nm technology node generated by the online model generation tool available on the ASU PTM website [START_REF] Cao | Asu predictive technology model website[END_REF]. The simulations were carried out on a SunOS workstation (16 CPUs, 96GB RAM).

To obtain the experimental results, we initially vary the value of ε and determine the optimum (lowest) value for maximum leakage savings (i.e., leakage increase in any module up to ε% is tolerated). Figs. 14 and15 show the variation of leakage savings with tolerance in different designs for top down and bottom up approaches, respectively. From Fig. 14, it can be observed that the optimum tolerance value for top down approach is 10% (i.e., leakage increase in any module up to 10% is tolerated). The low leakage vector is captured within 10% tolerance itself. Even if the tolerance is increased above 10%, it is observed that leakage savings do not increase any further. From Figure 15, the optimum tolerance value for bottom up approach is observed to be 15%.

Tables 1 and2 report the results for top down and bottom up approaches respectively. Column 3 presents the leakage value obtained from interval propagation with simulated annealing. We can see a significant improvement compared to the random case in Column 2. Column 4 reports leakage savings. The top down approach achieved average leakage savings of 93.6% for 10% tolerance with no area overhead. On the other hand, bottom up approach achieved 89.2% average leakage savings for 15% tolerance with no area overhead.

We conducted an experiment to compare the proposed interval arithimetic followed by pure SA based approach. Table 3 reports the simulation speed up obtained when we use simulated annealing on reduced set of intervals from top down propagation as opposed to pure simulated annealing. Columns 2 and 4 present the leakage values. Execution time is reported in columns 3 and 5 in Table 3. It shows that simulated annealing with top down propagation is much faster than simulated annealing alone to obtain a similar solution quality. The difference in solution quality is also presented in column 6 of Table 3. For EWF benchmark, a solution that is much better than pure simulated annealing solution is obtained with interval propagation and simulated annealing combined. For the rest of the benchmarks, the difference between solutions is small. Column 7 reports speed up resulting from IA with as much as 124X in case of FIR filter. The results of this experiment demonstrate that Interval Arithmetic greatly helps in finding a good low leakage vector quickly. Table 4 compares the solution quality obtained by pure simulated annealing(SA) for 10 minutes, 1 hour, 2 hours, and 3 hours with that by top down approach with simulated annealing (an average of 10 minutes). Columns 2-5 report the leakage obtained by pure SA. Column 6 reports the leakage obtained from top down approach. It is observed that the leakage found by interval propagation with SA in 10 minutes is better than the solution found by pure SA in 3 hours. These results reiterate the efficacy of interval propagation based approach.

Conclusion

We have formulated the low leakage vector identification technique based on interval propagation and successfully demonstrated that significant subthreshold leakage savings can be obtained with no area overhead (i.e., no internal control points). We also proposed a self similarity based module characterization procedure that is scalable with module complexity. Both top down and bottom up approaches are equally effective (although in case of the benchmarks tested in this work, top down performs slightly better than bottom up approach).

Fig. 1 .

 1 Fig. 1. Romanesco Broccoli -An example of self-similar behavior occurring in nature. Photo Credit [12].

Fig. 2 .

 2 Fig.2. H tree -An example of a VLSI structure exhibiting self-similar behavior. Photo Credit[START_REF] Eppstein | Wikimedia commons[END_REF].

Fig. 3 .Fig. 4 .

 34 Fig. 3. Leakage current distribution for 8 bit adder.

Fig. 5 .

 5 Fig. 5. Leakage current distribution for 16 bit adder (1000 random vectors).

Fig. 6 .

 6 Fig. 6. Leakage current distribution for 16 bit parallel multiplier (1000 random vectors).

Fig. 7 .

 7 Fig. 7. Eight bit adder -Measured H values of leakage distributions at sub-space level with: (a) vertical partitioning into 4 sub-spaces; (b) horizontal and vertical partitioning into four sub-spaces; and (c) sixteen sub-spaces. Note that in each partition H value satisfies the self-similarity condition (0.5 < H < 1).

Fig. 9 .

 9 Fig. 9. Example 1: Top down approach -backward propagation of interval sets.

Fig. 10 .

 10 Fig. 10. Example 2: Bottom up approach -backward propagation of output intervals.

Fig. 11 .

 11 Fig. 11. Algorithm to determine Low Leakage Vector.

Fig. 12 .

 12 Fig. 12. Interval back propagation algorithm.

Fig. 13 .

 13 Fig. 13. Simulated annealing algorithm to find the best LLV.

  Time Complexity. In the first step of find LLV(), topological sort has a time complexity of O(|V | + |E|). The maximum number of edges in a DAG is (|V |)(|V | -1)/2. Hence the time complexity of topological sort is O(|V | 2 ). The for loop (lines 8-27) runs |V | times. If the node is PI, initial interval set is initialized at the node. This initialization takes a constant time. We also have the Interval Propagate() and Interval Intersection() which take time proportional to I 2 , where I is the number of low leakage intervals. The number of intervals I is a constant obtained from characterization as described in Subsec. 3.3. Back Propagate() function on lines 28 and 31 has the same complexity as that of the code on lines 8-27. Finally, simulated annealing algorithm takes a constant time based on the quality of solution desired. Hence, the run time complexity of the algorithm is O(|V | 2 ).

Fig. 14 .

 14 Fig. 14. Leakage vs. Tolerance -Top down Approach.

Fig. 15 .

 15 Fig. 15. Leakage vs. Tolerance -Bottom up Approach.

Table 1 .

 1 Power Savings -Top down Approach.

	Design	Leakage (µA)	Savings (%)
		Random Top down
	Diffeq ( 2+, 5*)	197.9	16.42	91.70
	EWF (26+, 8*)	654.0	64.13	90.19
	FIR ( 4+, 5*)	220.9	9.95	95.50
	IIR ( 4+, 5*)	266.8	11.94	95.50
	Lattice ( 8+, 5*)	226.4	11.01	95.10

Table 2 .

 2 Power Savings -Bottom up Approach.

	Design	Leakage (µA)	Savings (%)
		Random Bottom up
	Diffeq (2+, 5*)	197.9	21.6	89.0
	EWF (26+, 8*)	654.0	72.8	88.9
	FIR (4+, 5*)	220.9	25.8	88.3
	IIR (4+, 5*)	266.8	25.3	90.5
	Lattice (8+, 5*)	226.4	23.9	89.4

Table 3 .

 3 Speed up with Interval Propagation.

	Design	Interval + SA	SA		Leakage	Speed
		Leakage Time Leakage Time Improvement	up
		(µA)	(min)	(µA)	(min)	(%)
	Diffeq	20.60	12	18.75	170	+8.98	14x
	EWF	55.24	13	71.62	812	-29.65	62x
	FIR	16.50	2	16.30	247	+1.21	124x
	IIR	16.55	18	14.46	244	+12.63	14x
	Lattice	23.00	6	21.83	440	+5.09	73x

Table 4 .

 4 Simulated Annealing (SA) Only vs. Interval Propagation + SA.

	Design		Leakage with SA		Interval + SA
		10min	1hr	2hr	3hr	Leakage Time
		(µA) (µA) (µA) (µA)	(µA)	(min.)
	Diffeq	32.08 27.65 25.24 22.20	20.60	12
	EWF	99.12 83.50 87.18 85.95	55.24	13
	FIR	34.74 31.80 33.10 26.06	16.50	2
	IIR	45.53 27.56 29.90 28.16	16.55	18
	Lattice 41.24 36.64 32.71 33.21	23.00	6


Dominik Auras
email: auras@ice.rwth-aachen.de

Uwe Deidersen

Rainer Leupers

Gerd Ascheid

A Parallel MCMC-based MIMO Detector: VLSI Design And Algorithm

Stochastic detection for multi-antenna (MIMO) systems promises communications performance close to max-log detection for certain SNR regimes, especially when the system iterates between detector and channel decoder following the Turbo Principle. In this work, we propose a parallel VLSI architecture for soft-input soft-output Markov chain Monte Carlo based stochastic MIMO detection. It features runtime adaptibility to varying channel conditions, effectively allowing us to adjust the invested effort. Besides the details of our area-throughput efficient design, like the low-level algorithm and micro-architecture design, we also provide an extensive data set from our experiments regarding the detector's communications performance and relate it to our VLSI implementation results. The provided data analysis highlights the architecture's run-time adaptibility and demonstrates how we can trade off throughput for improved communications performance.

Introduction

With the wide-spread adoption of MIMO (multi-antenna) technology in current and future wireless communication systems, such as those based on the IEEE 802.11n standard [START_REF]Wireless LAN medium access control (MAC) and physical layer (PHY) specifications amendment 5: Enhancements for higher throughput. IEEE Std[END_REF], academia and industry are searching for MIMO detectors with reasonable implementation complexity and algorithmic performance. Especially for systems using bit-interleaved coded modulation with iterative decoding (BICM-ID) [START_REF] Guillén I Fàbregas | Bit-interleaved coded modulation[END_REF], a major challenge for VLSI implementation is the required soft-input soft-output (SISO) detector, since optimal detection has an exponential complexity.

Iterative MIMO decoding can yield impressive algorithmic performance gains in terms of significantly reduced signal-to-noise ratio (SNR) requirements to achieve a certain fixed error rate [START_REF] Hochwald | Achieving near-capacity on a multiple-antenna channel[END_REF]. This SNR gain has several possible uses amongst others: we can extend the transmission range, we can serve more users (i.e. tolerating more interference), we can lower the transmission power to save energy (and at the same time reduce interference to other users), or transmit at a higher throughput in the same bandwidth.

Possible detectors can be roughly put into two categories: linear detectors, e.g. MMSE-filter based [START_REF] Auras | A novel reduced-complexity soft-input softoutput mmse mimo detector: Algorithm and efficient vlsi architecture[END_REF][START_REF] Studer | ASIC implementation of soft-input soft-output MIMO detection using MMSE parallel interference cancellation. Solid-State Circuits[END_REF][START_REF] Auras | A novel class of linear mimo detectors with boosted communications performance: Algorithm and vlsi architecture[END_REF], and non-linear detectors e.g. [START_REF] Hochwald | Achieving near-capacity on a multiple-antenna channel[END_REF][START_REF] Witte | Efficiency and flexibility trade-offs for soft-input soft-output spheredecoding architectures[END_REF][START_REF] Chen | VLSI implementation of a high-throughput iterative fixed-complexity sphere decoder[END_REF][START_REF] Sun | Trellis-search based soft-input soft-output MIMO detector: Algorithm and VLSI architecture[END_REF]. Basically, linear detectors try to suppress noise using linear filtering, then decode the estimate.

In contrast to this, non-linear detectors perform a search, e.g. a randomly guided one, in the space of possibly transmitted data vectors. Stochastic detection based on Markov chain Monte Carlo (MCMC) methods [START_REF] Senst | A rao-blackwellized markov chain monte carlo algorithm for efficient MIMO detection[END_REF] belongs to this class. It enables small configurable detectors that can cover a large design space. Furthermore, when iterating between detector and channel decoder, MCMC detection shows a communications performance close to max-log detection for certain SNR regimes [START_REF] Senst | A rao-blackwellized markov chain monte carlo algorithm for efficient MIMO detection[END_REF].

To date, only some research effort has been directed towards this field. There exist only a handful of publications on MCMC detector architectures at the moment [START_REF] Laraway | Implementation of a markov chain monte carlo based multiuser/MIMO detector[END_REF][START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF][START_REF] Auras | VLSI design of a parallel MCMC-based MIMO detector with multiplier-free gibbs samplers[END_REF][START_REF] Yuan | A hardware-efficient VLSI architecture for hybrid sphere-MCMC detection[END_REF]. None of them correlates communications performance with VLSI implementations results.

Related Work An MCMC-based SISO MIMO detector ASIC design supporting independent parallel Gibbs Samplers is presented in [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF]. Amongst other things, [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF] introduces an initialization scheme for the completely recursive, and thus simplified, computation of the detector states, and shows how to reuse the circuitry to draw independent first samples. However, a multiplier in the timing critical path yields a limited throughput and a relatively large area consumption.

In [START_REF] Laraway | Implementation of a markov chain monte carlo based multiuser/MIMO detector[END_REF], the authors propose an MCMC-based SISO MIMO detector architecture mapped on an FPGA. It features one multiplier-free Gibbs Sampler pipelined at the symbol vector level. The architecture uses a simple recursive metric computation, but requires one dot-product per cycle. The first sample of every chain needs to be generated externally.

The hybrid soft-output only MCMC detector architecture [START_REF] Yuan | A hardware-efficient VLSI architecture for hybrid sphere-MCMC detection[END_REF] combined with a hard-output fixed-complexity sphere detector (FSD) features parallel multiplier-free Gibbs Samplers that start with the best candidates found by the FSD. However, the design requires the QR-decomposition of the channel matrix, and the results are only given in terms of operation counts.

Contribution We present a complete redesign of the MCMC-based MIMO detector architecture presented in [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF], with multiplier-free Gibbs Samplers and further architectural improvements that result in a significant area reduction and timing improvement. Post-layout area and clock period reduce by about 50% and 40% respectively. In extension to our previous publication [START_REF] Auras | VLSI design of a parallel MCMC-based MIMO detector with multiplier-free gibbs samplers[END_REF], we additionaly provide our detector's communications performance results and present an analysis showing how to trade off throughput for improved communications performance at run-time.

Outline First, we introduce the general concept of MCMC-based MIMO detection (Sec. 3), describe the implemented algorithm (Sec. 4), then we propose the redesigned architecture (Sec. 5). Subsequently, we explicitly highlight the differences to the reference design [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF] in Sec. 6. Our implementation results are presented in Sec. 7. The analysis of the communications performance results is explained in Sec. 8.

Encoder

System Model

We consider a spatial-multiplexing N t × N r MIMO system with BICM-ID, as depicted in Fig. ] = N 0 I Nr . In the remainder, the time index n is dropped for convenience. Using iterative MIMO decoding following the Turbo Principle [START_REF] Hagenauer | The turbo principle in mobile communications[END_REF], detector and channel decoder exchange extrinsic information λ e = λ pλ a in terms of log-likelihood ratios (LLRs), where λ p are the detector's posterior LLRs and λ a are the prior LLRs fed back from the decoder.

MCMC-based MIMO Detection

The Markov chain Monte Carlo based MIMO detector class that we consider performs a randomly guided search in the space c ∈ {0, 1} KNt . It starts with a random candidate, then walks around randomly. On its way, it evaluates and saves metric values of the current candidates, which are later used to approximate the posterior LLRs. The random process (Monte Carlo) from which it draws new candidates evolves recursively (Markov chain). By design the search converges towards candidates of high probability [START_REF] Senst | A rao-blackwellized markov chain monte carlo algorithm for efficient MIMO detection[END_REF].

We select independent first samples c (q,0) ∈ {0, 1} KNt , one per chain q = 1 . . . N q , either randomly from the prior distribution c (q,0) ∼ p(c) = f (λ a ) or given by an external hard-output detector c (q,0) = c ext (usually for at most one chain). Every sample s = 1 . . . N s is drawn in KN t steps. The algorithm sequentially replaces every bit with 0 and 1, computes the metric for those two candidates, then selects one of them as the next partial sample.

Let ϕ : {0, 1} NtK → X Nt be a rule that maps bit labels onto symbol vectors x ∈ X Nt . We define the metric

µ(c) = - 1 N 0 y -Hϕ(c) 2 -c T λ a (1)
for the candidate c ∈ {0, 1} KNt , which is related to the posterior probability P (c|y, H, λ a ). Furthermore, let

c bβ = (c 1 , . . . , c b-1 , β , c b+1 , . . . , c KNt ) (2)
be the vector c with the b-th bit replaced by β. The detector approximates the posterior LLRs as

λ p b ≈ max q,s µ(c (q,s) b0 ) -max q,s µ(c (q,s) b1 ) (3)
where we search for the two maxima for every bit over all chains and samples.

Low-Level Algorithm

The presented algorithm implements the max-log variant of the Rao-Blackwellized MCMC detection algorithm with uniform sampling described in [START_REF] Senst | A rao-blackwellized markov chain monte carlo algorithm for efficient MIMO detection[END_REF]. Its basic idea is to recursively compute the metric in Eq. ( 1) by tracking the changes while drawing bits [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF]. First, we introduce the basic concepts required for understanding the algorithm, then describe the algorithm in detail.

For the theoretic background, the reader is kindly referred to [START_REF] Senst | A rao-blackwellized markov chain monte carlo algorithm for efficient MIMO detection[END_REF][START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF].

Basic Concepts

Matched Filter The algorithm in [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF] replaces H with the Gram matrix R = H H H and the received symbol vector y with the matched filter output y mf = H H y in the metric. This does not influence the posterior LLR calculation, however it allows to use the symmetry R = R H .

Gibbs Sampler (GS)

We realize the Markov chains with Gibbs Sampling. To this end, the GS draw bits sequentially according to an approximation of the marginal distribution P (c b |c 1 , . . . , c b-1 , c b+1 , . . . , c KNt ). The state of the q-th GS at the s-th sample after drawing the b-th bit is denoted as

c (q,s) b = (c (q,s) 1
, . . . , c (q,s) b , c (q,s-1) b+1

, . . . , c

(q,s-1) KNt ) (4)
and thus contains bits from the previous sample c (q,s-1) and the current sample c (q,s) .

Common Starting Point All chains start with c (-1) , which maps onto x (-1) with x t = 1 + j, i.e. we have ϕ(c (-1) ) = x (-1) . This concept enables the initialization of parallel independent Gibbs Samplers [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF].

Symbol Deltas When the GS state changes, at most one bit is different. We introduce the notation

|∆| 2 b (c) = |ϕ n (c b1 )| 2 -|ϕ n (c b0 )| 2 ∆ b (c) = ϕ n (c b1 ) -ϕ n (c b0 ) (5)
where ϕ n is the mapping rule for the n-th antenna, and the b-th bit belongs to the n-th antenna.

Recursive Dot-Product The algorithm tracks the current value of

S = y mf -Rϕ(c (q,s) b ) ( 6
)
where R is the matrix R with the diagonal set to zero. Starting from S (-1) = y mf -Rx (-1) , it updates S recursively when c (q,s) b changes.

Recursive Metric Computation

We introduce an arbitrary offset such that µ(c (-1) ) = 0, which cancels out in Eq. ( 3). Let the distance update be

δ (q,s) b = Re{r nn }|∆| 2 b (c (q,s-1) ) -2Re{S * n ∆ b (c (q,s-1) )} (7)
where the b-th bit belongs to the n-th antenna, then the metric update is

∆µ = 1 N 0 δ (q,s) b + λ a b ( 8
)
which we either subtract from or add to the current metric µ(c (q,s) ), depending on the bit flip direction, if the b-th bit changes.

Log-domain Bit Probability

The term

γ = 1 ηN 0 δ (q,s) b + λ a b (9)
expresses the probability of the next bit being 1 in the log-domain, where the temperature parameter η mitigates lock-in effects in the high-SNR regime [START_REF] Senst | A rao-blackwellized markov chain monte carlo algorithm for efficient MIMO detection[END_REF].

For the conversion to the linear domain, we apply a piece-wise linear approximation to logistic(γ) = 1/(1+e -γ ) as in [START_REF] Laraway | Implementation of a markov chain monte carlo based multiuser/MIMO detector[END_REF][START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF]. To this end, the GS simply limits γ to the range [-4, 4) and compares -γ to a uniformly distributed pseudo-random number u ∼ U (-4, 4) in the same range.

Overall Algorithm Design

Front-end Processing

First, choose Γ = 2 α /(ηN 0 ) with α such that Γ ∈ [0.5, 1). We assume η = 2.

The FEP computes

R = Γ H H H S (-1) = Γ H H y -Rx (-1) (10)
as described in Sec. 4.1 but scaled by Γ .

Gibbs Sampler

Alg. 1 describes how the GS sequentially draws bits of the candidate sequence c (q,s) . GS and Metric Update share the term δ (q,s) b computed in line 6. Note the back-shifting with α to compensate the normalization of Γ . For the first sample (s = 0), only the prior LLRs are used, in order to draw c (q,0) ∼ λ a (line 7). The saturation in line 8 produces a threshold in the range [-4, 4) (cf. Sec. 4.1). The comparison to a uniformly distributed pseudorandom number in the same range (line 13) yields the new bit value. Afterwards, we need to update the S state (lines [START_REF] Yuan | A hardware-efficient VLSI architecture for hybrid sphere-MCMC detection[END_REF][START_REF] Hagenauer | The turbo principle in mobile communications[END_REF][START_REF] Patel | Vlsi implementation of a WiMAX/LTE compliant low-complexity high-throughput soft-output K-Best MIMO detector[END_REF].

Algorithm 1: Gibbs Sampler input: S (-1) , R, c ext , λ a , Chain Index q output: c (q,s) b , c (q,s-1) b , δ (q,s) b 1 c (q,-1) = c (-1) 2 S ← S (-1) 3 for s = 0 to Ns do 4 for b = 1 to NtK do 5 n ← (b -1)/K + 1 6 δ (q,s) b = Re{rnn}|∆| 2 b (c (q,s-1) ) -2Re{S * n ∆ b (c (q,s-1) )} 2 -α 7 γ ← λ a b + 0 if s = 0 δ (q,s) b otherwise 8 γ ← saturate(-4, 4, γ) 9 draw u ∼ U (-4, 4) 10 if s = 0 and q = 1 then /* first sample, first chain */ 11 c (q,s) b = c ext b 12 else 13 c (q,s) b = sign(u + γ) 14 ∆St ← rtn∆ b (c (q,s-1) ) ∀ t = 1 . . . Nt, t = n 15 if c (q,s-1) b = c (q,s) b then 16 St ← St + ∆St if c (q,s-1) b = 1 -∆St if c (q,s-1) b = 0 ∀ t = n

Metric Update

Alg. 2 recursively computes the current candidate's metric µ(c (q) b ), using the state µ (q) , and produces the two metrics for the current bit µ(c b0/1 ). As stated earlier, we arbitrarily set the metric for the common starting point to zero (line 1). Lines 4 to 9 show the underlying metric update. Of the two possible states, one is identical to the current state, and thus has the same metric value (line 4). The other one is updated according to the direction of the bit flip (lines 6 and 8). In line 9, we select one of the two as the new current metric. It remains unaltered if the bit does not change.

Algorithm 2: Metric Update input: c (q,s) b , c (q,s-1) b , δ (q,s) b , λ a , Chain Index q output: µ(c (q,s) b0 ), µ(c (q,s) b1 ) 1 µ (q) ← 0 2 for s = 0 to Ns do 3 for b = 1 to NtK do 4 µ(c (q,s) b0 ) = µ(c (q,s) b1 ) = µ (q) 5 if c (q,s-1) b = 0 then 6 µ(c (q,s) b1 ) = µ(c (q,s-1) b ) -(ηδ (q,s) b + λ a b ) 7 else 8 µ(c (q,s) b0 ) = µ(c (q,s-1) b ) + (ηδ (q,s) b + λ a b ) 9 µ(c (q,s) b ) = µ(c (q,s) b0 ) if c (q,s) b = 0 µ(c (q,s) b1 ) if c (q,s) b = 1 10 µ (q) ← µ(c (q,s) b )

LLR Computation

Alg. 3 searches for the maximum metrics among all chains, then compares these local maxima with the current global maxima. It excludes the s = 0 step, which is the transition from c (-1) to c (q,0) , from the search (line 3). The computation of the extrinsic LLRs in line 7 is included, as it can be easily implemented in hardware.

Algorithm 3: LLR Computation input: µ(c (q,s) b0 ), µ(c (q,s) b1 ), λ a output: λ e 1 µ max b0 ← -∞ ∀ b = 1 . . . NtK 2 µ max b1 ← -∞ ∀ b = 1 . . . NtK 3 for s = 1 to Ns do /* Note: ignore input for s = 0 */ 4 for b = 1 to NtK do /* For every bit index */ 5 µ max b0 ← max(µ max b0 , max q (µ(c (q,s) b0 ))) 6 µ max b1 ← max(µ max b1 , max q (µ(c (q,s) b1 ))) 7 λ e b = µ max b0 -µ max b1 -λ a b ∀ b = 1 . . . KNt

VLSI Architecture

Overview

The macro pipeline of FEP-Circuit and MCMC core, shown in Fig. 3, constitutes the proposed MCMC detector. Both components require multiple clock cycles per input vector, but double buffering between FEP and Core ensures that the computations can overlap. The MCMC core in turn contains four stages connected via registers. The stages exchange information in every clock cycle. They effectively run in a pipeline manner.

The FSM and the multiplexers (e.g. λ a b , and for the column of R) are part of the Mux stage. There are N p GS-Circuits implementing Alg. 1. For every GS-Circuit, there is one corresponding M-Circuit executing Alg. 2. The L-Circuit performs the LLR Computation in Alg. 3. Every GS/M-Circuit can run several chains sequentially. For example N q = 8 chains can be run on N p = 4 GS/M-Circuits by executing two chains sequentially per GS/M-Circuit. We can also turn off some GS/M-Circuits, e.g. run N q = 4 chains on N p = 8 GS/M-Circuits with four inactive circuits.

FEP-Circuit

The architecture, depicted in Fig. 4, contains in total five multipliers. Using four of these, the dot-product for the terms H H y and R = H H H requires N r cycles per complex entry. We need only the lower triangular of R due to R H = R. The architecture computes either one complex off-diagonal entry, or two real diagonal entries in parallel. The fifth multiplier alternatingly multiplies real and imaginary parts with Γ = 2 α N0η . In parallel, we multiply the entries of R with x (-1) t = 1 + j (cf. Sec. 4.1) using only adders and multiplexers, and accumulate the results to obtain S. , then finishes in the write-enable control for the S registers. The M-Circuit, shown in Fig. 7, implements Alg. 2 using a write-enabled register for the current metric, which is updated when we flip the current bit. The multiplication with η is implemented as a constant shift.

GS/M-Circuit

c (q,s) b ∆ b (c (q,s-1) ) control signals |∆| 2 b (c (q,s-1) ) δ (q,s) b
+ +/- µ ×η c (q,s) b c (q,s-1) b µ(c (q,s) b ) µ(c (q,s) b0 ) µ(c (q,s) b1 ) δ (q,s) b λ a b c (q,s-1) b Fig. 7. M-Circuit
Update-S-Circuit The Update-S-Circuit shown in Fig. 8 has (N t -1) complexvalued ∆-multipliers, i.e. 2(N t -1) times Fig. 6(b). Using the multiplexers 3 and 4 , we can update all N t elements of S, however only N t -1 change per clock cycle. The entries of R e.g. r 1n , r 2n are selected in the Mux stage. Similar to the GS-Circuit, the adder-subtractor control 1 considers ∆ < 0, if |∆| is imaginary, and additionally the old bit c

L-Circuit

The L-Circuit shown in Fig. 9 contains two register files (RFs) for the current maximum metrics with KN t entries each. We use tokens propagating alongside the data to indicate whether a value is valid. The Compare Select (CS) elements select the maximum of the valid inputs. The registers also store tokens per entry, which are reset to zero when the processing of a symbol vector starts. After the scalar subtractor, we saturate the extrinsic LLRs to limit their dynamic range. The saturation has a positive influence on the communications performance. The proposed architecture is a complete redesign of [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF]. This section explicitly highlights the architectural modifications. The original and new timing critical path are located in the GS-Circuit.

Multiplier-Free Gibbs Sampler: Similar to [START_REF] Laraway | Implementation of a markov chain monte carlo based multiuser/MIMO detector[END_REF], we move the multiplication with 1/(ηN 0 ) out of the GS into the FEP, by scaling R and S with Γ . This removes the multiplier from the detector's critical path, but increases the required word lengths.

Dynamic Scaling: The normalization of Γ ∈ [0.5, 1) allows to use smaller word lengths, mitigating the previously mentioned increase. Consequently, we need an arithmetic shifter in the GS-Circuit at the previous location of the multiplier, which reverts the normalization.

Pipelined Input Multiplexers: Our MCMC detector selects the column of R and the entry of λ a in the new Mux stage in front of the GS stage. While this removes those multiplexers from the detector's critical path, it adds an additional latency cycle.

Reduced Update-S-Circuit: We remove two ∆-multipliers (one per real and imaginary part) from the Update-S-Circuit, since in every cycle one of the entries of S does not change. This requires multiplexers for the resource sharing, which are however not in the critical path and are smaller than the removed ∆-multipliers.

Shared Maximum Metric Register File: The RFs are moved from the M-Circuit [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF] to the L-Circuit. This reduces the required RFs from N p to one. We also add a pipeline register after the L-Circuit to improve timing, which requires another extra latency cycle. Also, our M-Circuit in Fig. 5 has one addersubtractor instead of two adders, similar to [START_REF] Laraway | Implementation of a markov chain monte carlo based multiuser/MIMO detector[END_REF].

Adder-Subtractor Units: These new units right after the ∆-multipliers in the GS-and the Update-S-Circuit, replace the original adders and the conditional negation units. The control selects addition or subtraction depending on the sign of ∆, if ∆ is imaginary, the old bit c (q,s-1) b and if Im{r tn } = -Im{r nt }. Simplified Delta Multiplier: Our ∆-multipliers, used for γ and S, compute the absolute value |∆|. This removes one multiplexer stage from the critical path.

Postponed Conjugation: We are storing only the lower half of R. Due to the hermitian property of R, we have Im{r tn } = -Im{r nt }. The control of the subsequent adder-subtractor units considers the required negation, instead of an explicit conjugation [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF].

Results

With the word lengths given in Sec. 7.1 and the throughput equations in Sec. 7.2, we first compare our model to the reference architecture [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF] based on gate-level synthesis results, then we present post-layout results for different design-time variants of our architecture. Sec. 8 presents the algorithmic evaluations.

Simulation Setup

A 802.11n-like 4 × 4 MIMO system is considered assuming a spatially uncorrelated Rayleigh channel, perfect channel knowledge and a max-log BCJR decoder. For all results, we assumed a rate-5/6 tail-biting binary convolutional code with generator polynomials 0133 and 0171 and puncturing, a random interleaver and 64-QAM modulation (K = 6). The frame length of 2160 information bits equals the interleaver's length, which is one OFDM symbol for this setup. For every data point, we simulated at least 10 5 frames. The average signal-to-noise ratio (SNR) per receive antenna is defined as SNR = E[ Hx 2 ]/(N r N 0 ). The required word lengths for an SNR loss of ≤ 0.1dB compared to the floating-point model at a frame error rate (FER) of 10% are: The first chain (q = 0) is always initialized with the result of an hard-output zero-forcing MIMO detector. We assume N q = 8 chains with N s = 8 samples per chain (i.e. N gs = 64 in [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF]) for the next three sections, but vary those parameters in Sec. 8.

Architecture

Our parameterized architecture implementation currently supports up to 4 × 4 MIMO and 64-QAM. MIMO mode and QAM scheme can be configured at runtime within the supported range, which in turn can be configured at design-time.

Each GS/M-pair can process up to 16 chains sequentially, with up to 16 samples per chain. The FEP-Circuit requires

n fep = N r ((N t + 1)N t /2 + N t /2 ) + 3 (11)
cycles for its computation. This is slightly faster than the FEP-Circuit in [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF]. The MCMC core runs for

n gs = N q N p (N s + 1)KN t + 5 (12)
cycles. Compared to [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF], we need two extra latency cycles (cf. Sec. 6). The code bit throughput of the architecture is θ c = KNt ngs f clk assuming n gs ≥ n fep and sufficient input data.

Synthesis Results

We synthesized the design with Synopsys Design Compiler I-2013.12-SP2 in topographical mode using a 1.0V standard-performance standard cell library for the UMC 90nm SP-RVT LowK CMOS process. One gate-equivalent (GE) is the area of one 2-input drive-1 NAND gate. Fig. 10 compares the four instances N p = {1, 2, 4, 8} to [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF]. While the most efficient design in [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF] has an AT exec -product of 181.7 kGEµs, our proposed design achieves 50.0 kGEµs, which is 3.6 times more efficient. Tbl. 1 lists the synthesis results for our fastest design instance and the reference design [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF]. The FEP is larger (5 kGE), while the GS is smaller (-6.2 kGE per GS), since we moved the multiplier from the GS to the FEP. The additional area of the new arithmetic shifter is partially compensated for by the other improvements. The Update-S-Circuit becomes smaller (-2.5 kGE) since we save one complex ∆-multiplier and use |∆| now. The saving effect is larger than the additional area from the multiplexers required for the resource sharing. The M-Circuit exhibits only about 7.4% of the original area, since we moved the RFs to the L-Circuit, which consequently became larger (10 kGE). The remainder of the area (-12.9 kGE) is occupied amongst others by the R column multiplexers. The area is reduced because the multiplexers are no longer in the timing critical path.

In total, the redesigned architecture takes on only about 48% of the original area for N p = 8. The saving depends on the number of GS/M-Circuits. The critical path was shortened by about 40%, i.e. the maximum clock frequency increased from 312 MHz to 526 MHz.

Layout Results

A layout was obtained with Cadence SoC Encounter 9.1 for each configuration's fastest design instance in order to further study the proposed architecture's implementation complexity and to enable more precise comparison with future related work. All following area figures are taken from the layout results, depicted in Fig. 11. The consumed area slightly increased, while the achievable clock frequency decreased. It is interesting that the throughput mainly depends on the number of parallel GS/M-Circuits and the chain parameters, i.e.

θ c = KN t n gs f clk ≈ N p N q (N s + 1) f clk (13)
as can be seen in Fig. 11. The largest instance, for 64-QAM, N t = 4 and N p = 8, requires 149.5 kGE or 0.47 mm 2 and achieves a maximum clock frequency of 479 MHz, yielding a code bit throughput of 52 Mbit/s. The fastest instance in terms of throughput supports 4-QAM, N t = 2 and has N p = 8 GS/M-Circuits. It occupies in total an area of 70.7 kGE or 0.22 mm 2 and runs at 664 MHz, which results in a throughput of 66 Mbit/s.

To determine the smallest instance, which should be the lower corner of the covered design space, 1 in Fig. 11, we synthesized the detector with N t = 2, 4-QAM and one GS/M-Circuit for a target of 100 MHz. This ASIC consumes 19.2 kGE or 0.06 mm 2 , runs at 165 MHz and yields a 2.27 Mbit/s throughput. The FEP-Circuit and MCMC core require 10.9 kGE and 8.3 kGE respectively. Further word length optimizations could yield additional area reductions.

Tbl. 2 compares our work to a selection of reported MIMO detector implementations. We make three observations. First, in terms of hardware efficiency expressed in Mbit/s/kGE, the MCMC detector resides in about the same order of magnitude as the single-tree-search sphere decoder (STS-SD) [START_REF] Witte | Efficiency and flexibility trade-offs for soft-input soft-output spheredecoding architectures[END_REF], though our architecture is more than two times more efficient than our reference architecture [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF]. The MCMC detector exhibits a deterministic run-time, which eases the receiver system design, while the SD can in principle always achieve near-capacity performance at the cost of a strongly varying run-time. Secondly, b We assume 100 visited nodes per symbol vector [START_REF] Studer | ASIC implementation of soft-input soft-output MIMO detection using MMSE parallel interference cancellation. Solid-State Circuits[END_REF]. c Scaled to 90nm CMOS technology assuming t pd ∼ 1/s d Area for optional initial detection not included e Area for chip IO interface excluded f Area for required QRD not included the MCMC detectors (and the STS-SD) are about one order of magnitude less efficient than the linear [START_REF] Auras | A novel reduced-complexity soft-input softoutput mmse mimo detector: Algorithm and efficient vlsi architecture[END_REF][START_REF] Studer | ASIC implementation of soft-input soft-output MIMO detection using MMSE parallel interference cancellation. Solid-State Circuits[END_REF], iterative-linear [START_REF] Auras | A novel class of linear mimo detectors with boosted communications performance: Algorithm and vlsi architecture[END_REF] detectors, and most notably the fixed-complexity sphere decoder (FCSD) [START_REF] Chen | VLSI implementation of a high-throughput iterative fixed-complexity sphere decoder[END_REF], which achieves close-to-optimal communications performance at a deterministic run-time. In this perspective, the FCSD [START_REF] Chen | VLSI implementation of a high-throughput iterative fixed-complexity sphere decoder[END_REF] is the best choice. In case that a particularly small implementation is needed, the MCMC might have an advantage, depending on how well the FCSD scales. Lastly, there are three cases for the preprocessing circuitry. Some implementations include it [START_REF] Auras | A novel reduced-complexity soft-input softoutput mmse mimo detector: Algorithm and efficient vlsi architecture[END_REF][START_REF] Studer | ASIC implementation of soft-input soft-output MIMO detection using MMSE parallel interference cancellation. Solid-State Circuits[END_REF][START_REF] Auras | A novel class of linear mimo detectors with boosted communications performance: Algorithm and vlsi architecture[END_REF], it is optional for the MCMC detectors [START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF], and definitely required for the other reported work [START_REF] Witte | Efficiency and flexibility trade-offs for soft-input soft-output spheredecoding architectures[END_REF][START_REF] Chen | VLSI implementation of a high-throughput iterative fixed-complexity sphere decoder[END_REF][START_REF] Sun | Trellis-search based soft-input soft-output MIMO detector: Algorithm and VLSI architecture[END_REF][START_REF] Patel | Vlsi implementation of a WiMAX/LTE compliant low-complexity high-throughput soft-output K-Best MIMO detector[END_REF]. This of course makes the area-throughput efficiency comparison difficult.

Algorithmic Considerations

In this section, we put the code bit throughput θ c , as an implementation property of our architecture, in relation to our design's communications performance in terms of SNR required to achieve a 10% frame error rate. With this data, we can determine for example appropriate run-time parameters, or an appropriate run-time strategy to adapt them. Depending on the optimization criterion, the parameter choices might be different. Possible criteria are for example spectral efficiency or energy efficiency (as future work, we plan to perform energy estimations). The first part of this section gives a general overview, while the second part explains in more detail the iterative receiver figures.

In the remainder, we use the post-layout implementation results of the 64-QAM, N t = 4, N p = 8 instance that runs at 479 MHz. The simulation setup that we select resembles the highest-throughput mode of the 802.11n standard, which requires a high SNR. However, our experiments show that the MCMCbased detection performs best in a mid-range SNR regime, in combination with lower-order modulation schemes. Thus this can be considered as kind of a worstcase scenario for the MCMC detector.

We assume the same simulation setup as in Sec. 7.1. Additionally, we perform up to two detector-decoder iterations, i.e. per frame, we execute the MCMC detector and BCJR decoder twice. This gives us four run-time parameters: the number of chains N q1 and samples N s1 in the first iteration and respectively N q2 , N s2 for the second iteration. The short-hand notation GS 1 8x6 denotes N q1 = 8 and N s1 = 6, similarly we use GS 2 N q2 xN s2 . We simulated the parameter set N q1/2 ∈ {8, 16} and N s1/2 = {1, 2, . . . , 16}. Thus all N p = 8 GS/M-Circuits are always active. The total number of samples per iteration defined as N gs1/2 = N q1/2 • N s1/2 is our measure for the invested effort.

Fig. 12 shows four curves: two for the first iteration, and two for the second. The last part of this section explains how we determine the two second-iteration curves. They are pareto-optimal in terms of SNR versus throughput.

Clearly in Fig. 12, we can identify the existence of a run-time tradeoff between SNR and throughput. As could be expected, more effort (i.e. more samples, more chains) results in a better algorithmic performance (lower SNR). An SNR gain has several possible uses amongst others: we can extend the transmission range, we can serve more users (more interference), or we can also lower the transmission power to save energy (and reduce interference to other users).

In the non-iterative case (first iteration), we observe a vanishing gain beyond five samples, both for eight and 16 chains. At around 33.5 dB, it is better to use 16 instead of eight chains. Interestingly, this switches from GS 1 8x8 at 33.52 dB to GS 1 16x4 at 32.53 dB. The total number of samples for both configurations is 64, but we gain about 1 dB SNR while approximately maintaining the throughput. It is not completely identical due to the pipeline delays of the architecture.

Instead of using GS 1 8x6 after GS 1 8x5, a good decision would be to switch to the second iteration, therefore never using 16 chains in the non-iterative case. This yields a large SNR gain of about 2.7 dB at a similar throughput. At this transition point, we switch from GS 1 8x5 to GS 1 8x2-GS 2 8x2. The throughputs drops slightly from 77.15 Mbit/s to 74.65 Mbit/s. With N gs1 = 40 compared to N gs1 + N gs2 = 32, the MCMC detector's effort remains very similar.

MCMC-based detection benefits greatly from iterative MIMO decoding. Switching from one to two iterations yields SNR gains as large as 6 dB. While in the first iteration we achieve only about 31.7 dB, all SNR operating points of the second iteration are lower than 31 dB. A possible explanation is that the guidance from the channel decoder, in terms of prior LLRs, is the contributing factor for this. It helps the MCMC-based detection in two ways: we select the initial samples c (q,0) ∼ p(c) = f (λ a ), and the transition probability γ depends on λ a . This seems to let the chains converge faster (in less samples) to interesting regions.

It follows a closer look on the second iteration. There are four parameters, N q1/2 and N s1/2 . For a given SNR, we determine the parameter combination that yields the highest throughput. These pareto-optimals points are shown in Fig. 13. For the two second-iteration curves in Fig. 13, we fix the number of chains in the first iteration N q1 = 8 and N q1 = 16 respectively, then optimize over the remaining three parameters.

For our calculations, we assume that the channel decoder and the buffering between decoder and detector cause no additional delay. This is a somewhat ideal scenario, since it might give us a large area consumption e.g. of the buffers, but definitely provides us with an upper bound for the achievable throughput. Thus, the throughput is given as θ c,2 = KN t n gs,1 + n gs,2 f clk [START_REF] Yuan | A hardware-efficient VLSI architecture for hybrid sphere-MCMC detection[END_REF] with n gs,1/2 = N q1/2 Np N s1/2 + 1 KN t + 5 and fixing N p = 8 here. We observe that more effort is required in the first iteration. For example, around 27 dB, the two configurations GS 1 8x7 and GS 2 8x4 are in use, i.e. N gs1 = 56 total samples for the first iteration, and N gs2 = 32 for the second. At about 26.8 dB, we switch from eight to 16 chains in the first iteration. It appears that multiple short chains are favorable for the first iteration. Only at around 24 dB, the detector should switch from eight to 16 chains in the second iteration. It is also the point where the effort significantly rises (near 23.6 dB), especially for the second iteration. This could be an indication for switching to three iterations.

From a pure SNR-throughput perspective, we can say that two iterations are better than a single. As previously stated, we observe large SNR gains from iterating, and the best non-iterative operating point is off by about 0.7 dB compared to the worst second-iteration point. However, this of course ignores the hardware cost caused by the required buffering and the increased throughput requirement on the detector and decoder architectures. A realistic comparison depends on the overall objective, i.e. lowest energy, small area, best spectral efficiency, and on additional constraints, like minimum supported bandwidth. While this is out of scope here, we think that our data outlines the run-time adaptibility of the MCMC-based MIMO detection architecture. It also shows that it performs particularly well in iterative receivers, therefore it could be a reasonable candidate to consider in the design of such a system.

Conclusions & Outlook

We have presented synthesis and layout results of the proposed MCMC detector architecture. The area reduction of up to 52% and the shorter clock period by up to 40% indicate that the proposed architectural modifications to the reference design are effective. Our extensive data set for the communications performance further highlights the available tradeoff between signal-to-noise ratio and architecture throughput. With its run-time adaptibility covering a large design space, our detector is effectively able to cope with a lot of channel conditions at the appropriate effort. Though being a stochastic detector, its completely deterministic run-time eases scheduling at the system level, i.e. inside a complex iterative receiver.

Still, the architecture suffers from a relatively low but deterministic throughput, which stems from the MCMC detection method itself. The main advantage appears to be its simple scalability through N p and configurability through N t and K. This allows the architecture to cover a large design space. Practically, only the availability of sufficient data might limit the architectural parallelism.

As future work, we plan to correlate algorithmic performance with energy consumption, which might reveal another tradeoff capability of the proposed design.

Fig. 1 .

 1 Fig. 1. Assumed MIMO BICM-ID System Model. Detector and decoder iteratively exchange information to improve the final decoding result.

Fig. 2 Fig. 2 .

 22 Fig.2depicts the algorithm partitioned into four different parts: the Front-end Processing (FEP), that transforms the channel observations, the parallel Gibbs Samplers (GS) realizing the Markov chains, the Metric Update (M) tracking the current metric state, and the LLR Computation, which searches for the two maximum metric values per bit.

Fig. 3 .

 3 Fig. 3. Architecture design of the MCMC detector. The n-th column rn of R and λ a b

Fig. 5 Fig. 4 .

 54 Fig. 5 depicts the GS-Circuit. The |∆| 2 -multiplier, depicted in detail in Fig. 6(a), exploits the limited range of |∆| 2 ∈ {-3, -2, . . . , 3} × {8, 16} which assumes only 14 different values for 4-/16-/64-QAM. The factor ∆ is either

Fig. 5 .Fig. 6 .

 56 Fig. 5. GS-Circuit. The arithmetic shifter (ASH) reverts the normalization of Γ .

4 SnFig. 8 .

 48 Fig. 8. Update-S-Circuit. Example for Nt = 4 antennas. All units exist for the real and for the imaginary parts respectively (not drawn).

Fig. 9 .

 9 Fig. 9. L-Circuit

  [integer.fractional] y [7.8], H [3.8], λ a [5.4], 1/N 0 [6.11], R [6.10], S [9.9], δ [17.6], µ [19.5], γ [3.29], 2 α δ [14.6], α [4.0], λ e [8.4]. All are signed, per entry, and for real and imaginary part identical.

Fig. 10 .

 10 Fig. 10. Area vs. execution time based on the MCMC detector's synthesis results, comparing this work to [12], assuming Nt = 4, K = 6

1 2 Fig.

 2 Fig. Area vs. throughput based on the MCMC detector's layout results. For each design-time configuration, the ASIC with the fastest clock is shown. As an example, the 16-QAM 2 × 2 design supports one or two antennas and 4-or 16-QAM at run-time.

a

  We assume Ns = 64, Ngs = 8, Np = 8[START_REF] Deidersen | A parallel VLSI architecture for markov chain monte carlo based MIMO detection[END_REF].

Fig. 12 .

 12 Fig. 12. Code bit throughput over SNR required to achieve a 10% frame error rate

Fig. 13 .

 13 Fig.[START_REF] Auras | VLSI design of a parallel MCMC-based MIMO detector with multiplier-free gibbs samplers[END_REF]. Iterative receiver: two detector-decoder activations per symbol vector. For the throughput curve, dots and crosses denote if eight or 16 chains are used in the first of the two iterations. For the two total samples curves, the dots and crosses likewise denote if eight or 16 chains are used in the respective iteration. The figure shows only the pareto-optimal points in terms of SNR versus throughput determined from the data set with two iterations.

  1. A message b ∈ {0, 1} N b is encoded with rate r = N b /N

c and interleaved, yielding the code word c ∈ {0, 1} Nc . Let X ⊂ C be a modulation alphabet with K = log 2 |X | bits per symbol. The code word is partitioned into multiple subvectors c n ∈ {0, 1} KNt . They are subsequently mapped to symbol vectors x n ∈ X Nt that are transmitted independently. Assuming a frequency-flat fading channel characterized by H n ∈ C Nr×Nt , the received symbol vector at time n is y n = H n x n + w n where w n ∈ C Nr is a white Gaussian noise process with E[w n w H n

Table 1 .

 1 MCMC Detector Synthesis Results

	Component		This Work [12]
	FEP-Circuit		16.0	11.0 kGE
	GS-Circuit	8×	10.7	16.9 kGE
	M-Circuit	8×	0.9	12.2 kGE
	L-Circuit	Np = 8	13.3	3.3 kGE
	Miscellaneous		5.0	17.9 kGE
	Update-S-Circuit (cont. in GS)	5.2	7.7 kGE
	Total	Np = 8 127.1 265.0 kGE
	Clock frequency		526	312 MHz
	Cycles (Nq = Ns = Np = 8)	221	219
	Average throughput		57.4	34.2 Mbit/s
	Area efficiency		0.45	0.13 Mbit/s/kGE

Table 2 .

 2 Comparison to other reported MIMO detectors

	This work
	[12]
	[6]

Acknowledgements This work has been supported by the Ultra High-Speed Mobile Information and Communication (UMIC) Research Centre, RWTH Aachen University.