Petros Arun Singh
email: singhak@eurecom.fr

Joakim Jaldén
email: jalden@kth.se

Achieving a vanishing SNR-gap to exact lattice decoding at a subexponential complexity

The work identifies the first lattice decoding solution that achieves, in the general outage-limited MIMO setting and in the high-rate and high-SNR limit, both a vanishing gap to the error-performance of the (DMT optimal) exact solution of preprocessed lattice decoding, as well as a computational complexity that is subexponential in the number of codeword bits. The proposed solution employs lattice reduction (LR)-aided regularized (lattice) sphere decoding and proper timeout policies. These performance and complexity guarantees hold for most MIMO scenarios, all reasonable fading statistics, all channel dimensions and all full-rate lattice codes.

In sharp contrast to the above very manageable complexity, the complexity of other standard preprocessed lattice decoding solutions is revealed here to be extremely high. Specifically the work is first to quantify the complexity of these lattice (sphere) decoding solutions and to prove the surprising result that the complexity required to achieve a certain rate-reliability performance, is exponential in the lattice dimensionality and in the number of codeword bits, and it in fact matches, in common scenarios, the complexity of ML-based solutions. Through this sharp contrast, the work was able to, for the first time, rigorously demonstrate and quantify the pivotal role of lattice reduction as a special complexity reducing ingredient.

Finally the work analytically refines transceiver DMT analysis which generally fails to address potentially massive gaps between theory and practice. Instead the adopted vanishing gap condition guarantees that the decoder's error curve is arbitrarily close, given a sufficiently high SNR, to the optimal error curve of exact solutions, which is a much stronger condition than DMT optimality which only guarantees an error gap that is subpolynomial in SNR, and can thus be unbounded and generally unacceptable for practical implementations.

I. INTRODUCTION

The work applies to the general setting of outage-limited MIMO communications, where MIMO techniques offer significant advantages in terms of increased throughput and reliability, although at a cost of a potentially much higher computational complexity for decoding at the receivers. This high complexity brings to the fore the need for efficient decoders that tradeoff error-performance with complexity in a better manner than computationally expensive decoders like the strictly optimal maximum-likelihood (ML) decoder.

Specifically in terms of ML-based decoding, the use of the brute-force ML decoder, introduces a complexity that scales exponentially with the number of codeword bits. If on the other hand, a small gap to the exact ML performance is acceptable, then different branch-and-bound algorithms such as the sphere decoder (SD) have been known to accept reduced computational resources. Despite the reduced complexity of sphere decoding, recent work in [START_REF] Jaldén | Sphere decoding complexity exponent for decoding full rate codes over the quasi-static mimo channel[END_REF] has revealed that, to achieve a vanishing error-gap to optimal ML solutions, even such branch-and-bound algorithms generally require computational resources that, albeit significantly smaller than those required by a brute-force ML decoder, again grow exponentially in the rate and the dimensionality, and remain prohibitive for several MIMO scenarios.

This high complexity required by ML-based decoding solutions, serves as further motivation for exploring other families of decoding methods. A natural alternative is lattice decoding obtained by simply removing the constellation boundaries of the ML-based search, an action that loosely speaking exploits a certain symmetry which in turn may yield faster implementations. It is the case though that even with lattice decoding, the computational complexity can be prohibitive: finding the exact solution to the lattice decoding problem is generally an NP hard problem (cf. [START_REF] Micciancio | The hardness of the closest vector problem with preprocessing[END_REF]). At the same time though, the other extreme of very early terminations of lattice decoding, such as linear solutions, have been known to achieve computational efficiency at the expense though of a very sizable, and often unbounded, gap to the exact solution of the lattice decoding problem.

In this work we explore lattice decoding solutions that, in conjunction with terminating policies, strike the proper balance between this exponential complexity and exponential gap.

A. System model

We consider the general m×n point-to-point multiple-input multiple-output model given by

y = √ ρHx + w (1)
where x ∈ R m , y ∈ R n and w ∈ R n respectively denote the transmitted codewords, the received signal vectors, and the additive white Gaussian noise with unit variance, where the parameter ρ takes the role of the signal to noise ratio (SNR), and where the fading matrix H ∈ R n×m is assumed to be random, with elements drawn from arbitrary statistical distributions. We consider that one use of (1) corresponds to T uses of some underlying "physical" channel. We further assume the transmitted codewords x to be uniformly distributed over some codebook X ∈ R m , to be statistically independent of the channel H, and to satisfy the power constraint

E{ x 2 } ≤ T. (2)

B. Rate, reliability and complexity in outage-limited MIMO communications

In terms of error performance, we let P e denote the probability of codeword error, and we consider the rate,

R = 1 T log |X |, (3)
in bits per channel use (bpcu), where |X | denotes the cardinality of X .

Regarding complexity, we let N max describe the computational resources, in floating point operations (flops) per T channel uses, that the transceiver is endowed with, in the sense that after N max flops, the transceiver must simply terminate, potentially prematurely and before completion of its task. We note that naturally, N max is intimately intertwined with the desired P e and R, and that any attempt to significantly reduce N max may be at the expense of a substantial degradation in error-performance.

In the high SNR regime, a given encoder X r and decoder D r are said to achieve a multiplexing gain r (cf. [START_REF] Zheng | Diversity and multiplexing: A fundamental tradeoff in multiple-antenna channels[END_REF] In the same high SNR regime, the complexity is here chosen to take the form

c(r) := lim ρ→∞ N max log ρ , (5)
which is henceforth denoted as the complexity exponent.

Noting that R = r log ρ, we observe that c(r) > 0 implies a complexity that is exponential in the rate. Remark 1: A reasonable question at this point would pertain as to why the computational resources N max scale with ρ and are dependent on r, to which we note that the complexity of decoding is generally dependent on the density of the codebook, which in turn depends on ρ and R. Furthermore this dependence of the complexity exponent (and by extension of N max ) on r, reflects a potential ability to regulate the computational resources depending on the rate. Finally the fact that both P e and N max are represented as polynomial functions of ρ, simply stems from the fact that both P e and |X | naturally scale as polynomial functions of ρ. Specifically we quickly note that c(r) captures the entire complexity range 0 ≤ c(r) ≤ rT of all reasonable transceivers, with c(r) = 0 corresponding to the fastest possible transceiver (requiring a subexponential number of flops per T channel uses), and with c(r) = rT corresponding to the optimal but arguably slowest, full-search uninterrupted ML decoder 1 in the presence of a canonical code with multiplexing gain r, i.e., with |X r | = 2 RT = ρ rT . If this canonical code though is linear, searching the entire codebook can be avoided by algorithmic solutions like the sphere decoder (SD) which can provide substantial complexity reductions at a potential small loss in error performance. Such solutions take advantage of the linear nature of the code that is defined by a generator matrix G and a shaping region R ′ . Specifically for r ≥ 0, a (sequence of) full-rate linear (lattice) code(s) X r is given by X r = Λ r ∩ R ′ where Λ r ρ -rT κ Λ and Λ {Gs | s ∈ Z κ }, where Z κ denotes the κ = min{m, n} dimensional integer lattice, where R ′ is a compact convex subset of R κ that is independent of ρ, and where G ∈ R m×κ is full rank and independent of ρ. For the class of lattice codes considered here, the codewords take the form

x = ρ -rT κ Gs, s ∈ S κ r Z κ ∩ ρ rT κ R, (6)
where R ⊂ R κ is a natural bijection of the shaping region R ′ that preserves the code, and where R contains the all zero vector 0.

As noted before, despite the reduced complexity of sphere decoding of such lattice codes (as compared to brute-force ML decoding), recent work in [START_REF] Jaldén | Sphere decoding complexity exponent for decoding full rate codes over the quasi-static mimo channel[END_REF] has revealed that even such branch-and-bound algorithms generally require computational resources that grow exponentially in the number of codeword bits and the dimensionality. As an indicative example of this high complexity, we note that the work in [START_REF] Jaldén | Sphere decoding complexity exponent for decoding full rate codes over the quasi-static mimo channel[END_REF] showed that such SD algorithms, when applied for decoding a large family of high-performing codes including all known full-rate DMT optimal codes, over the n T × n R quasi-static MIMO channel with Rayleigh fading and n R ≥ n T , introduce a complexity exponent 2 of the form

c(r) = T n T r(n T -⌊r⌋ -1) + (n T ⌊r⌋ -r(n T -1)) + . (7)
In the above, ⌊r⌋ denotes the largest integer not greater than r. The exponent, which simplifies to c(r) = T nT r(n Tr) for integer values of r, reaches at r = n T /2 (for even values of n T ) an overall maximum value of n T T /4 which, for the aforementioned codes is equal to κ/8, corresponding to complexity in the order of 2

1 8 κ log ρ = ρ κ/8 = |X |.
At any fixed multiplexing gain, these required computational resources can be seen to be in the order of 2

RT ( n T -r n T
) flops which reveals a complexity that is exponential in the number of codeword bits, and a corresponding exponential slope of nT-r nT . 1 We here note that strictly speaking, Xr, Dr may potentially introduce a complexity exponent larger than rT . In such a case though, Xr, Dr may be substituted by a lookup table implementation of Xr and an unrestricted ML decoder. This encoder-decoder will jointly require resources that are a constant multiple of |Xr| . = ρ rT as it has to construct and visit all possible |Xr| codewords, at a computational cost of a bounded number of flops per codeword visit. It is noted that the number of flops per visited codeword is naturally independent of ρ.

2 Although premature at this point, we hasten to note for the expert reader that this complexity indeed holds irrespective of the radius updating policy, irrespective of the decoding ordering, and as we will see later on, holds even in the presence of MMSE preprocessing.

C. Transition to lattice decoding for reducing complexity

As mentioned, this high complexity of ML based (constrained) decoders, motivates consideration of other decoder families, with a natural alternative being the unconstrained (naive) lattice decoder which takes the general form

xL = arg min x∈Λr y - √ ρHx 2 . ( 8
)
Naturally when xL / ∈ X r , the decoder declares an error. The use of lattice decoding, and specifically of preprocessed lattice decoding in MIMO communications has received substantial attention from works like [START_REF] Agrell | Closest point search in lattices[END_REF], [START_REF] Murugan | A unified framework for tree search decoding: rediscovering the sequential decoder[END_REF] and [START_REF] Gamal | Lattice coding and decoding achieve the optimal diversity-multiplexing tradeoff of MIMO channels[END_REF], where the latter proved that lattice decoding in the presence of MMSE preprocessing achieves the optimal DMT for specific MIMO channels and statistics, and for DMT-optimal random codes. The use of lattice decoding as an alternative to computationally expensive ML based solutions, was recently further validated on the one hand by the aforementioned work in [START_REF] Jaldén | Sphere decoding complexity exponent for decoding full rate codes over the quasi-static mimo channel[END_REF], [START_REF] Jaldén | The complexity of sphere decoding perfect codes under a vanishing gap to ml performance[END_REF] which revealed the large computational disadvantages of ML based solutions, and on the other hand by the work in [START_REF] Jaldén | DMT optimality of LR-aided linear decoders for a general class of channels, lattice designs, and system models[END_REF] which further confirmed the performance advantages of lattice decoding by showing that regularized (MMSE-preprocessed) 3 lattice decoding achieves the optimal DMT performance, for almost all MIMO scenarios and fading statistics, and all non-random lattice codes, irrespective of the codes' ML performance.

It is the case though that the aforementioned extreme complexity of exact lattice decoding solutions, in conjunction with the potentially unbounded error-performance degradation (gap) of very early terminations (as opposed to exact implementations) of lattice decoding, bring to the fore the need for balanced approximations of lattice decoding solutions that better balance the very sizable complexity and gap. Specifically for any simplified variant D r of the baseline (exact) MMSEpreprocessed lattice decoder, this gap can, in the high SNR regime, be quantified as

g L (c) lim ρ→∞ P e P (x = x) (9)
where P (x = x) describes the probability of error of the exact MMSE-preprocessed lattice decoder, where P e denotes the probability of error of D r , and where c (i.e., c(r)) is the complexity exponent that describes the (asymptotic rate of increase of the) computational resources required to achieve this performance gap. Generally a smaller computational complexity exponent c implies a larger gap g L (c). The clear task has remained for some time to construct decoders that optimally traverse this tradeoff between g and c, i.e., that reduce the performance gap to the exact lattice decoding solution, with reasonable computational complexity. Equivalently for N max (g) denoting the computational resources in flops required to achieve a certain gap g to the baseline exact MMSE-preprocessed lattice decoder, the above task can be described, in the high SNR regime, as trying to minimize lim ρ→∞ log N max (g) log ρ .

This will be achieved later on.

D. Contributions

We first show that the computational complexity required by the MMSE-preprocessed (unconstrained) lattice sphere decoder, asymptotically matches the complexity of the (constrained) ML-based (MMSE-preprocessed or not) sphere decoders, and is commonly exponential in the dimensionality and the number of codeword bits. This is established for a large class of codes of arbitrary error-performance, a large class of fading statistics, and specifically for the quasi-static MIMO channel -for example the complexity required for DMT optimal lattice sphere decoding, in the presence of a large family of DMT optimal codes, takes the previously seen simple piecewise linear form in [START_REF] Jaldén | The complexity of sphere decoding perfect codes under a vanishing gap to ml performance[END_REF]. In a parenthetical note, and deviating slightly from the spirit of this paper, we also provide a universal upper bound on the complexity of regularized lattice sphere decoding, which holds irrespective of the lattice code applied and irrespective of the fading statistics. This upper bound again takes the form in [START_REF] Jaldén | The complexity of sphere decoding perfect codes under a vanishing gap to ml performance[END_REF], matching that in the case of constrained ML-based sphere decoding, thus revealing the surprising fact that there exists no statistical channel behavior that will allow the removal of the bounding region to cause unbounded increases in the complexity of the decoder 4 .

With provable evidence of the very high complexity of regularized lattice decoding, we turn to the powerful tool of lattice reduction and seek to understand its effects on computational complexity. While there has existed a general agreement in the community that lattice reduction does reduce complexity, cf. [START_REF] Damen | On maximum-likelihood detection and the search for the closest lattice point[END_REF], this has not yet been supported analytically in any relevant communication settings. In fact, and quite opposite to common wisdom, it was recently shown that for a fixed-radius 5 sphere decoding implementation of the naive lattice decoder [START_REF] Seethaler | Tail behavior of sphere-decoding complexity in random lattices[END_REF], LR does not improve the sphere decoder complexity tail exponent.

What our present work shows is that lattice reduction reduces an ML-like exponentially increasing complexity, to very manageable subexponential values. We specifically proceed to prove that the LR-aided regularized lattice decoder, implemented by a fixed-radius sphere decoder and timeout policies that occasionally abort decoding and declare an error, achieves

g L (ǫ) = 1, lim ρ→∞ log N max (g) log ρ = 0 ∀ǫ > 0, g ≥ 1,
i.e., achieves a vanishing gap to the exact implementation of regularized lattice decoding and does so with a complexity exponent that vanishes to zero, which in turn implies subexponential complexity in the sense that the complexity scales slower than any conceivable exponential function. It is finally noted that this vanishing gap approach serves the practical purpose of an analytical refinement over basic diversity analysis which generally fails to address potentially massive gaps between theory and practice.

E. Notation

We use . = to denote the exponential equality, i.e., we write 4)). In this paper we use • to denote the smallest integer not smaller than the argument, • to denote the largest integer not larger than the argument, (•) H to denote the conjugate transpose of (•), (•) + to denote max{0, (•)} and vec(•) to denote the operation whereby the columns of the argument (•) are stacked to form a vector.

f (ρ) . = ρ B to

II. MMSE-PREPROCESSED LATTICE SPHERE DECODING COMPLEXITY

We proceed to describe the preprocessed lattice decoder, its sphere decoding implementation, and for a practical setting of interest that includes the quasi-static MIMO channel and common codes, to establish the decoder's computational complexity.

A. Lattice sphere decoding

Combining (1) and ( 6) yields the equivalent model

y = M r s + w (10)
where

M r = ρ 1 2 -rT κ HG ∈ R n×κ (11)
is a function of the multiplexing gain 6 r.

Consequently the corresponding naive lattice decoder in (8) takes the form (see for example [START_REF] Jaldén | DMT optimality of LR-aided linear decoders for a general class of channels, lattice designs, and system models[END_REF], also [START_REF] Damen | On maximum-likelihood detection and the search for the closest lattice point[END_REF])

ŝL = arg min ŝ∈Z κ y -Mŝ 2 . ( 12
)
As a result though of neglecting the boundary region, the above decoder declares additional errors if ŝL / ∈ S κ r , resulting in possible performance costs. These costs motivated the use of MMSE preprocessing which essentially regularizes the decision metric to penalize vectors outside the boundary constraint S κ r (cf. [START_REF] Jaldén | DMT optimality of LR-aided linear decoders for a general class of channels, lattice designs, and system models[END_REF]). Specifically the MMSE-preprocessed lattice decoder is obtained by implementing an unconstrained search over the MMSE-preprocessed lattice, and takes the form

ŝr-ld = arg min ŝ∈Z κ Fy -Rŝ 2 , (13)
where F and R are respectively the MMSE forward and feedback filters such that F = R -H M H , where

R H R = M H M + α 2 r I, (14)
6 For simplicity of notation we will, in most cases, denote Mr with M.

where α r = ρ -rT κ and where R is an upper-triangular matrix (more details can be found in Appendix D). For r Fy, the model transitions from [START_REF] Damen | On maximum-likelihood detection and the search for the closest lattice point[END_REF] to

r = R -H M H Ms + R -H M H w = R -H (R H R -α 2 r I)s + R -H M H w = Rs -α r 2 R -H s + R -H M H w = Rs + w ′ (15)
where

w ′ = -α 2 r R -H s + R -H M H w (16)
is the equivalent noise that includes self-interference (first summand) and colored Gaussian noise. Consequently the corresponding regularized lattice decoder takes the form

ŝr-ld = arg min ŝ∈Z κ r -Rŝ 2 , (17)
which is then solved by the sphere decoder which recursively enumerates all lattice vectors ŝ ∈ Z κ within a given sphere of radius ξ > 0, i.e., which identifies as candidates the vectors ŝ that satisfy

r -Rŝ 2 ≤ ξ 2 . ( 18
)
The algorithm specifically uses the upper-triangular nature of R to recursively identify partial symbol vectors ŝk

, k = 1, • • • , κ, for which r k -R k ŝk 2 ≤ ξ 2 , (19)
where ŝk and r k respectively denote the last k components of ŝ and r, and where R k denotes the k×k lower-right submatrix of R. Clearly any set of vectors ŝ ∈ Z κ , with common last k components that fail to satisfy [START_REF] Windpassinger | Low-complexity nearmaximum-likelihood detection and precoding for MIMO systems using lattice reduction[END_REF], may be excluded from the set of candidate vectors that satisfy [START_REF] Yao | Lattice-reduction-aided detectors for MIMO communication systems[END_REF]. The enumeration of partial symbol vectors ŝk is equivalent to the traversal of a regular tree with κ layers -one layer per symbol component of the symbol vectors, such that layer k corresponds to the kth component of the transmitted symbol vector 7 s. There is a one-to-one correspondence between the nodes at layer k and the partial vectors ŝk . We say that a node is visited by the sphere decoder if and only if the corresponding partial vector ŝk satisfies [START_REF] Windpassinger | Low-complexity nearmaximum-likelihood detection and precoding for MIMO systems using lattice reduction[END_REF], i.e., there is a bijection between the visited nodes at layer k and the set

N k {ŝ k ∈ Z k | r k -R k ŝk 2 ≤ ξ 2 }. (20)

B. Complexity of MMSE-preprocessed lattice sphere decoding

Consequently the total number of visited nodes (in all layers of the tree) is given by

N SD = κ k=1 N k , (21)
where N k |N k | is the number of visited nodes at layer k of the search tree. The total number of visited nodes is commonly taken as a measure of the sphere decoder complexity. It is easy to show that in the scale of interest the SD complexity exponent c(r) would not change if instead of considering the number of visited nodes, we considered the number of flops spent by the decoder 8 .

Naturally the total number of visited nodes is a function of the search radius ξ. We here use a fixed radius, which may result in a non-zero probability that the transmitted symbol vector s is not in N κ . Consequently we must choose a radius that strikes the proper balance between decreasing the aforementioned probability and at the same time sufficiently decreasing the size of N κ . Towards this we note that for the transmitted symbol vector s, the metric in [START_REF] Zhao | Diversity and multiplexing tradeoff in general fading channels[END_REF] satisfies

r -Rs 2 = w ′ 2 ,
which means that if w ′ > ξ, then the transmitted symbol vector is excluded from the search, resulting in a decoding error. As Lemma 2 will later argue taking into consideration the self-interference and non-Gaussianity of w ′ , we can set ξ = √ z log ρ, for some z > d(r) such that

P w ′ 2 > ξ 2 < ρ -d(r) ,
which implies a vanishing probability of excluding the transmitted information vector from the search, and a vanishing degradation of error performance. We here note that the MMSE-preprocessed lattice sphere decoder differs from its ML-based equivalent in two aspects: the presence of MMSE preprocessing and the absence of a bounding region to constrain the search. These two aspects are generally perceived to have an opposite effect on the complexity. On the one hand, MMSE preprocessing, which we recall from [START_REF] Wübben | Lattice reduction: A survey with with applications to wireless communication[END_REF] to introduce unpruned sets

N k {ŝ k ∈ Z k | r k -R k ŝk 2 ≤ ξ 2 }, k = 1, • • • , κ,
is associated to reduced complexity in lattice-based SD solutions (cf. [START_REF] Seethaler | Tail behavior of sphere-decoding complexity in random lattices[END_REF]) due to the resulting penalization of faraway lattice points (cf. [START_REF] Jaldén | DMT optimality of LR-aided linear decoders for a general class of channels, lattice designs, and system models[END_REF]). On the other hand, the absence of boundary constraints can be associated to increased complexity as it introduces an unbounded number of candidate vectors. We proceed to show that in terms of the complexity exponent, under common MIMO scenarios and codes, these two aspects exactly cancel each other out, and that consequently MMSEpreprocessed lattice sphere decoding introduces a complexity exponent that matches that of ML-based sphere decoding (cf. [START_REF] Jaldén | Sphere decoding complexity exponent for decoding full rate codes over the quasi-static mimo channel[END_REF]), which it self is shown here to also match the complexity exponent of ML-based SD in the presence of MMSE preprocessing 9 . Before proceeding we note that this analysis is specific to sphere decoding, and that it does not account for any other ML based solutions that could, under some (arguably rare) circumstances, be more efficient. A classical example of such rare circumstances would be a MIMO scenario, or equivalently 8 To see this, we consider that the cost of visiting a node, is independent of ρ. Once at a visited node, this same bounded cost includes the cost of establishing which children-nodes not to visit in the next layer. 9 We clarify that ML-based SD in the presence of MMSE preprocessing, corresponds to unpruned sets N k ∩ S k r where S k r is the k-dimensional set resulting from the natural reduction of S κ r from [START_REF] Gamal | Lattice coding and decoding achieve the optimal diversity-multiplexing tradeoff of MIMO channels[END_REF].

a set of fade statistics, that always generate diagonal channel matrices. Another example would be having codes drawn from orthogonal designs which introduce very small decoding complexity, but which are provably shown to be highly suboptimal except for very few unique cases like the n T = 2, n R = 1 quasi-static case [START_REF] Radon | Lineare Scharen Orthogonaler Matrizen[END_REF]. In light of this, in this section only, we mainly focus on the widely considered n T × n R (n R ≥ n T ) i.i.d. and quasi-static MIMO setting and on the large but specific family of full-rate (κ = 2 min{n T , n R }T = 2n T T ) threaded codes (cf. [START_REF] Gamal | Universal space-time coding[END_REF]- [START_REF] Oggier | Perfect space-time block codes[END_REF]), which includes all known DMT optimal codes as well as uncoded transmission (V-BLAST).

We proceed with the main Theorem of the section, which applies under natural detection ordering (cf. [START_REF] Jaldén | Sphere decoding complexity exponent for decoding full rate codes over the quasi-static mimo channel[END_REF], [START_REF] Murugan | A unified framework for tree search decoding: rediscovering the sequential decoder[END_REF]), and under the assumption of i.i.d. regular fading statistics 10 .

Theorem 1: The complexity exponent for MMSEpreprocessed lattice sphere decoding any full-rate threaded code over the quasi-static MIMO channel with i.i.d. regular fading statistics, is equal to the complexity exponent of ML-based SD with or without MMSE preprocessing.

Proof: See Appendix A.
We clarify that even though all three decoders are DMT optimal, the above result incorporates more than just DMT optimal decoding, in the sense that any timeout policy will tradeoff d(r) with c(r) identically for ML-based and latticebased sphere decoding. In other words the three decoders share the same d(r) and c(r) capabilities, irrespective of the timeout policy.

Furthermore, considering different SD detection orderings (cf. [START_REF] Murugan | A unified framework for tree search decoding: rediscovering the sequential decoder[END_REF]), the following extends the range of codes for which the ML-based and lattice-based SD share a similar complexity. The proof follows from the proof of Theorem 1 in Appendix A, and from Theorem 4 in [START_REF] Jaldén | Sphere decoding complexity exponent for decoding full rate codes over the quasi-static mimo channel[END_REF].

Corollary 1a: Given any full-rate code of arbitrary DMT performance, there is always at least one non-random fixed permutation of the columns of G, for which the complexity exponent of the MMSE-preprocessed lattice sphere decoder matches that of the ML based sphere decoder.

The following focuses on a specific example of practical interest.

Corollary 1b: The complexity exponent for DMT optimal MMSE-preprocessed lattice sphere decoding of minimum delay (T = n T ) DMT optimal threaded codes over the quasistatic MIMO channel with i.i.d. regular fading statistics, takes the following form

c r-ld (r) = r(n T -⌊r⌋ -1) + (n T ⌊r⌋ -r(n T -1)) + , (22)
which simplifies to

c r-ld (r) = r(n T -r) (23)
for integer values of r.

Proof: See Appendix B. Further evidence that connects the complexity behavior of MMSE-preprocessed lattice-based SD, with that of its MLbased counterpart, now comes in the form of a non-trivial universal bound that is shared by the two methods. This is particularly relevant because unconstrained lattice decoding could conceivably require unbounded computational resources given the unbounded number of candidate lattice points. Specifically the following universal upper bound on the complexity of regularized lattice-based SD, matches the upper bound in [START_REF] Jaldén | Sphere decoding complexity exponent for decoding full rate codes over the quasi-static mimo channel[END_REF] for the ML case, and it holds irrespective of the full-rate lattice code applied and irrespective of the fading statistics. The generality with respect to the fading statistics is important because it guarantees that no set of fading statistics, even those that always generate infinitely dense lattices, can cause an unbounded increase in the complexity due to removal of the boundary constraints.

Corollary 1c: Irrespective of the fading statistics and of the full-rate lattice code applied, the complexity exponents of MMSE-preprocessed lattice SD and of ML-based SD, are upper bounded by

c(r) = T n T r(n T -⌊r⌋ -1) + (n T ⌊r⌋ -r(n T -1)) + (24)
which simplifies to

c(r) = T n T r(n T -r) (25)
for integer r.

Proof: See Appendix B. The above results revealed the very high, ML-like complexity of MMSE-preprocessed lattice decoding. Coming back to the main focus of the paper, and after reverting to the most general setting of MIMO scenarios, statistics and full-rate lattice codes, we proceed to show how proper utilization of lattice sphere decoding and LR techniques can indeed reduce the complexity exponent to zero, at an error-performance cost that vanishes in the high SNR limit.

III. LR-AIDED REGULARIZED LATTICE SPHERE DECODING COMPLEXITY

Lattice reduction techniques have been typically used in the MIMO setting to improve the error performance of suboptimal decoders (cf. [START_REF] Yao | Lattice-reduction-aided detectors for MIMO communication systems[END_REF], [START_REF] Windpassinger | Low-complexity nearmaximum-likelihood detection and precoding for MIMO systems using lattice reduction[END_REF], see also [START_REF] Wübben | Lattice reduction: A survey with with applications to wireless communication[END_REF], [START_REF] Ling | On the proximity factors of lattice reduction-aided decoding[END_REF]). In the current setting the LR algorithm, which is employed at the receiver after the action of MMSE preprocessing, modifies the search of the MMSE-preprocessed lattice decoder, from ŝrld = arg min ŝ∈Z κ r -Rŝ

by accepting as input the MMSE-preprocessed lattice generator matrix R, and producing as output the matrix T ∈ Z κ×κ which is unimodular meaning that it has integer coefficients and unit-norm determinant, and which is designed so that RT is (loosely speaking) more orthogonal than R. As a result of this unimodularity, we have that T -1 Z κ = Z κ , and consequently the new search in [START_REF] Horn | Matrix Analysis[END_REF] corresponds to yet another lattice decoder, referred to as the LR-aided MMSEpreprocessed lattice decoder, which operates over a generally better conditioned channel matrix RT.

Finally with sphere decoding in mind, the LR algorithm is followed by the QR decomposition 11 of the new latticereduced MMSE-preprocessed matrix RT, resulting in a new upper-triangular model r = Rs + w ′′ [START_REF] Dembo | Large deviations techniques and applications[END_REF] and in the new LR-aided MMSE-preprocessed lattice search, which accepts the application of the sphere decoder, and which takes the form

slr-rld = arg min ŝ∈Z κ r -Rŝ 2 , (28)
where Q R = RT corresponds to the QR-decomposition of RT, where R is upper-triangular, where r QH r, s = T -1 s, and where w ′′ = QH w ′ . At the very end,

ŝlr-rld = Ts lr-rld , (29)
allows for calculation of the estimate of the transmitted symbol vector s in [START_REF] Damen | On maximum-likelihood detection and the search for the closest lattice point[END_REF].

We note here that this (exact) solution of the LR-aided MMSE-preprocessed lattice decoder defined by ( 28), [START_REF] Tulino | Random matrix theory and wireless communications[END_REF], is identical to the exact solution of the MMSE-preprocessed lattice decoder given by [START_REF] Zhao | Diversity and multiplexing tradeoff in general fading channels[END_REF], because

min ŝ∈Z κ r -Rŝ 2 = min ŝ∈Z κ r -RTT -1 ŝ 2 (a) = min ŝ∈Z κ r -Q RT -1 ŝ 2 (b) = min ŝ∈Z κ r -RT -1 ŝ 2 = min ŝ∈T -1 Z κ r -Rŝ 2 (c) = min ŝ∈Z κ r -Rŝ 2 , (30)
where (a) follows from the fact that Q R = RT, (b) follows from the rotational invariance of the Euclidean norm, and (c) follows from the fact that T -1 Z κ = Z κ . While though the two lattice decoding solutions (with and without LR) provide identical error performance in the setting of exact implementations, we proceed to show that, in terms of complexity, lattice reduction techniques, and specifically a proper utilization of the LLL algorithm [START_REF] Lenstra | Factoring polynomials with rational coefficients[END_REF], can provide dramatic improvements.

A. Complexity of the LR-aided regularized lattice sphere decoder

We are here interested in establishing the complexity of the LR-aided regularized lattice sphere decoder. Given that the costs of implementing MMSE preprocessing and of implementing the linear transformation in [START_REF] Tulino | Random matrix theory and wireless communications[END_REF] are negligible in the scale of interest 12 , we limit our focus on establishing the cost of lattice reduction, and then the cost of the SD implementation of the search in [START_REF] Raghavan | When is limited feedback for transmit beamforming beneficial?[END_REF]. Starting with the SD complexity, as in [START_REF] Wübben | Lattice reduction: A survey with with applications to wireless communication[END_REF], we identify the corresponding unpruned set at layer k to be

N k {ŝ k ∈ Z k | rk -Rk ŝk 2 ≤ ξ 2 }, (31)
and in bounding the size of the above, we first focus on understanding the statistical behavior of the k × k lowerright submatrices Rk of matrix

R (k = 1, • • • , κ),
, k = 1, • • • , κ, satisfies P σ min ( Rk ) . < ρ -ǫT κ . ≤ ρ -d L (r-ǫ) , for all r ≥ ǫ > 0. (32)
To bound the cardinality N k of N k (cf. ( 31)), and eventually the total number N SD = κ k=1 N k of lattice points visited by the SD, we proceed along the lines of the work in [START_REF] Jaldén | Sphere decoding complexity exponent for decoding full rate codes over the quasi-static mimo channel[END_REF], making the proper modifications to account for MMSE preprocessing, for the removal of the bounding region, and for lattice reduction.

Towards this we see that, after removing the boundary constraint, Lemma 1 in [START_REF] Jaldén | Sphere decoding complexity exponent for decoding full rate codes over the quasi-static mimo channel[END_REF] tells us that

N k |N k | ≤ k i=1 √ k + 2ξ σ i ( Rk ) , where σ min ( Rk ) = σ 1 ( Rk ) ≤ • • • ≤ σ k ( Rk )
are the singular values of Rk . Consequently we have that

N k ≤ √ k + 2ξ σ min ( Rk ) k . ( 33
)
As a result, for any Rk such that

σ min ( Rk ) . ≥ ρ -ǫT κ , (34)
and given that ξ = √ z log ρ for some finite z, then 12 Even though the work here focuses on decoding, we can also quickly state the obvious fact that the cost of constructing the codewords is also negligible in the scale of interest because it again only involves a finite-dimensional linear transformation (cf. ( 6)).

N k . ≤ √ k + 2 √ z log ρ ρ -ǫT κ k . = ρ ǫT k κ , (35)
which guarantees that the total number of visited lattice points is upper bounded as

N SD = κ k=1 N k . ≤ κ k=1 ρ ǫT k κ . = ρ ǫT . (36)
Consequently, directly from Lemma 1, we have that

P N SD ≥ ρ ǫT ≤ ρ -d L (r-ǫ) . (37)
A similar approach deals with the complexity of the LLL algorithm, which is known (cf. [START_REF] Jaldén | Worst-and average-case complexity of LLL lattice reduction in MIMO wireless systems[END_REF]) to be generally unbounded. Specifically drawing from [8, Lemma 2], under the natural assumption of power-limited channels 13 (cf. [START_REF] Jaldén | DMT optimality of LR-aided linear decoders for a general class of channels, lattice designs, and system models[END_REF]), under the natural assumption that d L (rǫ) > d L (r) for all ǫ > 0, and for N LR denoting the number of flops spent by the LLL algorithm, one can readily conclude that

P (N LR ≥ γ log ρ) ≤ ρ -d L (r-ǫ) , (38)
for any γ > 1 2 (d L (r-ǫ)). Consequently the overall complexity

N . = N SD + N LR ,
in flops, for the LR-aided MMSE preprocessed lattice sphere decoder, satisfies the following

P N ≥ρ ǫT . = P {N SD ≥ρ ǫT } ∪ {N LR ≥ρ ǫT } . ≤ ρ -d L (r-ǫ) . ( 39
)
Now going back to [START_REF] Murugan | A unified framework for tree search decoding: rediscovering the sequential decoder[END_REF], and having in mind appropriate timeout policies that bound N max while at the same time specifically guarantee a vanishing error performance gap to the exact solution of regularized lattice decoding, we can see that the complexity exponent c(r) takes the equivalent form recently introduced (for the ML case) in [START_REF] Jaldén | Sphere decoding complexity exponent for decoding full rate codes over the quasi-static mimo channel[END_REF]

c(r) = inf{x | -lim ρ→∞ log P (N ≥ ρ x ) log ρ > d L (r)}. ( 40
)
To see this we quickly note that for N max = ρ x where x = c(r)δ for any δ > 0, it is the case that (cf. ( 9)) lim ρ→∞ P(N ≥ρ x ) P(x L =x) → ∞. Finally applying (39) we see that for any positive ǫ 1 < ǫ, it is the case that

c(r) = inf{ǫ | -lim ρ→∞ log P N ≥ ρ ǫT +ǫ1 log ρ > d L (r)} (41)
which vanishes arbitrarily close to zero, resulting in a zero complexity exponent. What remains is to consider the error-performance gap in the presence the LR-aided regularized lattice SD with a timeout policy that interrupts at N max = ρ x for any vanishingly small x > 0.

13 This is a moderate assumption that asks that E H 2 F . ≤ ρ. We note that this holds true for any telecommunications setting.

B. Gap to the exact solution of MMSE-preprocessed lattice decoding

We here prove that the LR-aided regularized lattice sphere decoder and the associated time-out policies that guarantee a vanishing complexity exponent, also guarantee a vanishing gap to the error performance of the exact lattice decoding implementation. This result is motivated by potentially exponential gaps in the performance of other DMT optimal decoders (cf. [START_REF] Jaldén | DMT optimality of LR-aided linear decoders for a general class of channels, lattice designs, and system models[END_REF]), where these gaps may grow exponentially up to 2 κ 2 (cf. [START_REF] Babai | On Lovász' lattice reduction and the nearest lattice point problem[END_REF]) or may potentially be unbounded [START_REF] Loyka | Finite-snr diversity-multiplexing tradeoff via asymptotic analysis of large mimo systems[END_REF].

Towards establishing this gap, we recall that the exact MMSE-preprocessed lattice decoder in (13) makes errors when ŝr-ld = s. On the other hand the LLL-reduced MMSEpreprocessed lattice sphere decoder with run-time constraints, in addition to making the same errors (ŝ r-lr-ld = s), also makes errors when the run-time limit of ρ x flops becomes active, i.e., when N ≥ ρ x , as well as when a small search radius causes N κ = ∅. Consequently the corresponding performance gap to the exact regularized decoder, takes the form

g L (x) = lim ρ→∞ P ({ŝ r-lr-ld = s} ∪ {N ≥ ρ x } ∪ {N κ = ∅}) P (ŝ r-ld = s) .
To bound the above gap, we apply the union bound and the fact that

P (N κ = ∅) ≤ P ( w ′′ > ξ)
to get that

Furthermore from [START_REF] Taherzadeh | LLL reduction achieves the receive diversity in MIMO decoding[END_REF] we observe that

P (ŝ r-lr-ld = s) = P (ŝ r-ld = s) , (43)
and from (39) we recall that

P N ≥ρ ǫT . ≤ ρ -d L (r-ǫ)
which implies that for any x > 0 it holds that

lim ρ→∞ P (N ≥ ρ x ) P (ŝ r-ld = s) = 0. ( 44
)
Finally the last term in (42) relates to the search radius ξ, and to the behavior of the noise w ′′ which was shown in ( 16), [START_REF] Dembo | Large deviations techniques and applications[END_REF] to take the form

w ′′ = QH -α 2 r R -H s + R -H M H w . ( 45
)
The following lemma, whose proof is found in Appendix D, accounts for the fact that w ′′ includes self-interference and colored noise, to bound the last term in (42).

Lemma 2: There exist a finite z > d L (r) for which a search radius ξ = √ z log ρ guarantees that

lim ρ→∞ P ( w ′′ > ξ) P (ŝ r-ld = s) = 0. (46)
Consequently combining (43), ( 44) and ( 46) gives that g L (x) = 1, ∀x > 0. The following directly holds.

Theorem 2: LR-aided MMSE-preprocessed lattice sphere decoding with a computational constraint activated at ρ x flops, allows for a vanishing gap to the exact solution of MMSEpreprocessed lattice decoding, for any x > 0. Equivalently the same LR-aided decoder guarantees that g L (ǫ) = 1 and lim

ρ→∞ log N max (g) log ρ = 0 ∀ǫ > 0, g ≥ 1,
for all fading statistics, all MIMO scenarios, and all full-rate lattice codes.

IV. CONCLUSIONS

The work identified the first lattice decoding solution that achieves, in the most general outage-limited MIMO setting and the high rate and high SNR limit, both a vanishing gap to the error-performance of the (DMT optimal) exact solution of preprocessed lattice decoding, as well as a computational complexity that is subexponential in the number of codeword bits. The proposed solution employs lattice reduction (LR)aided regularized lattice sphere decoding and proper timeout policies. As it turns out, lattice reduction is a special ingredient that allows for complexity reductions; a role that was rigorously demonstrated here for the first time, by proving that without lattice reduction, for most common codes, the complexity cost for asymptotically optimal regularized lattice sphere decoding is exponential in the number of codeword bits, and in many cases it in fact matches the complexity cost of ML sphere decoding.

In light of the fact that, prior to this work, a vanishing error performance gap was generally attributed only to near-full lattice searches that have exponential complexity, in conjunction with the fact that subexponential complexity was generally attributed to early-terminated (linear) solutions which have though a performance gap that can be up to exponential in dimension and/or rate, the work constitutes the first proof that subexponential complexity need not come at the cost of exponential reductions in lattice decoding error performance.

APPENDIX A PROOF FOR THEOREM 1 AND COROLLARY 1A

In the following we begin by providing an upper bound on the complexity exponent of MMSE-preprocessed (unconstrained) lattice sphere decoding, where this bound holds for the general quasi-static MIMO channel, for all fading statistics and for any full-rate lattice code. We will then proceed to provide a lower bound on the complexity exponent of the same decoder, where this bound, under the extra assumptions of regular i.i.d. fading statistics and of layered codes, will in fact match the above mentioned upper bound to prove the theorem and the associated corollaries. Before proceeding with the bounds, we describe the n T ×n R (n R ≥ n T ) quasi-static pointto-point MIMO channel, and its corresponding association to the general MIMO channel model in [START_REF] Damen | On maximum-likelihood detection and the search for the closest lattice point[END_REF] and metric in [START_REF] Zhao | Diversity and multiplexing tradeoff in general fading channels[END_REF].

The aforementioned quasi-static channel model takes the form

Y C = √ ρH C X C + W C , ( 47
)
where

X C ∈ C n T ×T , Y C ∈ C n R ×T and W C ∈ C n R ×T
represent the transmitted, received and noise signals over a period of T time slots, and where H C ∈ C nR×nT represents the matrix of fade coefficients. The real-valued representation of (47) can be written as

Y R = √ ρH R X R + W R , (48)
where

Y R = ℜ{Y C } -ℑ{Y C } ℑ{Y C } ℜ{Y C } , H R = ℜ{H C } -ℑ{H C } ℑ{H C } ℜ{H C } , X R = ℜ{X C } -ℑ{X C } ℑ{X C } ℜ{X C } and W R = ℜ{W C } -ℑ{W C } ℑ{W C } ℜ{W C }
, and subsequent vectorization gives the real-valued model

y = √ ρ(I T ⊗ H R )x + w ( 49
)
where y = vec(Y R ), x = vec(X R ), and w = vec(W R ). The system model in ( 49) is of the familiar form

y = √ ρHx + w (50)
as in ( 1) with m = 2n T T , n = 2n R T , and where

H = I T ⊗ H R . (51)
As before the vectorized codewords x, associated to the fullrate code, take the form

x = ρ -rT κ Gs, s ∈ Z κ ∩ ρ rT κ R, (52)
where κ = 2 min{n T , n R }T = 2n T T = m, which allows us to rewrite the model as

y = Ms + w, (53)
for

M = ρ 1 2 -rT κ HG = ρ 1 2 -rT κ (I T ⊗ H R )G. ( 54
)
Finally the corresponding coherent MMSE-preprocessed lattice decoder for the transmitted symbol vector s, can be expressed to be (cf. [START_REF] Zhao | Diversity and multiplexing tradeoff in general fading channels[END_REF])

ŝr-ld = arg min ŝ∈Z κ r -Rŝ 2 , ( 55
)
where r = Q H 1 y and R ∈ C κ×κ is the upper-triangular matrix, where furthermore both Q 1 and R result from the thin QR decomposition of the (n + κ) × κ dimensional preprocessed channel matrix

M reg M α r I = QR = Q 1 Q 2 R ( 56
)
and where as before α r = ρ -rT κ .

A. Upper bound on complexity of regularized lattice SD

In establishing the upper bound, we consider Lemma 1 in [START_REF] Jaldén | Sphere decoding complexity exponent for decoding full rate codes over the quasi-static mimo channel[END_REF], which we properly modify to account for MMSE preprocessing and for the removal of the constellation boundaries, and get that the number N k of nodes visited at layer k by the MMSE-preprocessed lattice sphere decoder, is upper bounded as

N k = |N k | ≤ k i=1 √ 2k + 2ξ σ i (R k ) , (57)
where

σ i (R k ), i = 1, • • • , k denote the singular values of R k in increasing order.
Towards lower bounding σ i (R k ), we note that

σ i (R k ) ≥ σ i (R) = σ i (M reg ) = α 2 r + σ i (M H M), ( 58
)
where the first inequality makes use of the interlacing property of singular values of sub-matrices [START_REF] Horn | Matrix Analysis[END_REF]. Furthermore for

µ j - log σ j (H H C H C ) log ρ , j = 1, • • • , n T (59)
and

µ 1 ≥ • • • ≥ µ n T , we see that σ j (H C ) = ρ -1 2 µj
, and from (54) that

σ i (M) ≥ ρ 1 2 -rT κ σ min (G)σ (i) (I T ⊗ H R )) . = ρ 1 2 -rT κ σ l 2T (i) (H C ) = ρ -rT κ + 1 2 (1-µ l 2T (i) ) , (60)
where l T (i) i T , and where the asymptotic equality is due to the fact that σ min (G) . = ρ 0 . Substituting from (60) in (58) we now have that

σ i (R k ) . ≥ ρ -rT κ + 1 2 (1-µ l 2T (i) ) + , i = 1, • • • , κ. (61)
Corresponding to (57) we see that

√ 2k + 2ξ σ i (R k ) . ≤ ρ ( rT κ -1 2 (1-µ l 2T (i) ) + ) + ,
for any i = 1, • • • , 2n T T , and from (57) we have that

N k (µ) . ≤ ρ k i=1 ( rT κ -1 2 (1-µ l 2T (i) ) + ) + , (62)
where

µ = (µ 1 , • • • , µ n T ). It follows that N SD (µ) = κ k=1 N k (µ) . ≤ κ k=1 ρ k i=1 ( rT κ -1 2 (1-µ l 2T (i) ) + ) + . = ρ κ i=1 ( rT κ -1 2 (1-µ l 2T (i) ) + ) + . = ρ T n T j=1 r n T -(1-µj ) + + , (63)
where the last asymptotic equality is due to the multiplicity of the singular values. Now consider the set

T (x)    µ | T n T j=1 r n T -(1 -µ j ) + + ≥ x    , (64)
and note that for any y < x, then (63) and µ / ∈ T (y) jointly imply that N SD < ρ x , which in turn implies that P (µ / ∈ T (y)) ≤ P (N SD < ρ x ) and consequently that

-lim ρ→∞ log P (N SD ≥ ρ x ) log ρ ≥ -lim ρ→∞ log P (µ ∈ T (y)) log ρ . (65)
In evaluating the right hand side of (65) we note that T (y) is a closed set and thus, applying the large deviation principle (cf.

[27]), we have that

-lim ρ→∞ log P (µ ∈ T (y)) log ρ ≥ inf µ∈T (y) I(µ) (66)
for some rate function I(µ). Consequently from (65) and (66), it follows that

-lim ρ→∞ log P (N SD ≥ ρ x ) log ρ ≥ inf µ∈T (y) I(µ). (67)
This lower bound specified in (67) holds for any y < x. Consequently to get the tightest possible bound, we need to find sup y<x inf µ∈T (y) I(µ). As inf µ∈T (y) I(µ) is nondecreasing and left-continuous in y, it follows that

sup y<x inf µ∈T (y)
I(µ) = inf µ∈T (x)
I(µ).

Consequently

lim

ρ→∞ log P (N SD ≥ ρ x ) log ρ ≥ inf µ∈T (x) I(µ), (68)
which in conjunction with (40) gives that

c r-ld (r) ≤ c r-ld (r) inf{x| inf µ∈T (x) I(µ) > d L (r)} = sup{x| inf µ∈T (x) I(µ) ≤ d L (r)} = max{x| inf µ∈T (x) I(µ) ≤ d L (r)} (69)
where the above follows from the aforementioned fact that lim

ρ→∞ log P(N SD ≥ρ x ) log ρ
(and by extension also inf µ∈T (x) I(µ))

is continuous and nondecreasing in x, and from the fact that T (x) is a closed set. Consequently c r-ld (r) takes the form c r-ld (r) max

µ x (70a) s.t. T n T j=1 r n T -(1 -µ j ) + + ≥ x, (70b)
I(µ) ≤ d L (r), (70c)
µ 1 ≥ • • • ≥ µ n T ≥ 0. (70d)
Furthermore since T (x) is a closed set, the maximum x in (70) must be such that (70b) is satisfied with equality, in which case c r-ld (r) can be obtained as the solution to a constrained maximization problem according to

c r-ld (r) max µ T n T j=1 r n T -(1 -µ j ) + + (71a) s.t. I(µ) ≤ d L (r), (71b) µ 1 ≥ • • • ≥ µ n T ≥ 0. (71c)
Equivalently for µ * = (µ * 1 , • • • , µ * n T ) being one of the maximizing vectors 14 , i.e., such that µ * ∈ T (x) and I(µ * ) = d L (r), then c r-ld (r) takes the form

c r-ld (r) = T n T j=1 r n T -(1 -µ * j ) + + . (72)
As we will now show, the above bound is also shared by the ML-based sphere decoder, with or without MMSE preprocessing, irrespective of the full-rate code and the fading statistics. Directly from [1, Theorem 2], and taking into consideration that MMSE-preprocessed lattice decoding is DMT optimal for any code [START_REF] Jaldén | DMT optimality of LR-aided linear decoders for a general class of channels, lattice designs, and system models[END_REF], we recall that the equivalent upper bound for the ML-based sphere decoder, without MMSE preprocessing, takes the form

c ml (r) max µ T n T j=1 min r n T -1 + µ j , r n T + (73a) s.t. I(µ) ≤ d L (r), (73b)
µ 1 ≥ • • • ≥ µ n T ≥ 0. (73c)
Comparing ( 71) and ( 73) we are able to conclude that both the objective functions (71a) and (73a) as well as both pairs of constraints are identical. To see this, we first note that for

0 ≤ µ j ≤ 1, then min r n T -1 + µ j , r n T + = r n T -1 + µ j + , r n T -(1 -µ j ) + + = r n T -1 + µ j + ,
and furthermore we note that for µ j > 1, then

min r n T -1 + µ j , r n T + = r n T -(1 -µ j ) + + = r n T ,
which proves that c ml (r) and c r-ld (r) are identical.

In considering the case of MMSE-preprocessed ML SD, it is easy to see that the summands in the objective function in (73a) will be modified to take the form

min r n T -(1 -µ j ) + , r n T +
which can be seen to match (71a) for all µ j ≥ 0, which in turn concludes the proof that the upper bound c r-ld (r) for MMSE-preprocessed lattice SD is also shared by the ML-based sphere decoder, with or without MMSE preprocessing, irrespective of the full-rate code, and for all fade statistics represented by monotonic rate functions.

B. Lower bound on complexity of regularized lattice SD

We will here, under the extra assumptions of regular i.i.d. fading statistics and of layered codes with natural decoding order, provide a lower bound that matches the upper bound in (72). The same bound and tightness will also apply to any full-rate code, under the assumption of a fixed, worst case decoding ordering.

The goal here is to show that at layer k = 2qT , for some q ∈ [1, n T ], the sphere decoder visits close to ρ c r-ld (r) nodes with a probability that is large compared to the probability of decoding error P (s L = s) . = ρ -d L (r) , which from the expression of the complexity exponent (40), will prove that c r-ld (r) = c r-ld (r).

Going back to (72), we let q be the largest integer for which

r n T -(1 -µ * q ) + > 0, (74)
in which case (72) takes the form

c r-ld (r) = T q j=1 r n T -(1 -µ * j ) + . (75)
We recall from (59) that µ j = -

log σj (H H C H C ) log ρ , j = 1, • • • , n T ,
and that µ * ∈ T (x) satisfies I(µ * ) = d L (r) and maximizes (71a). We also note that without loss of generality we can assume that q ≥ 1 as otherwise c r-ld (r) = 0 (cf. ( 72)). Consequently it is the case that µ * j > 0 for j = 1, • • • , q. Furthermore given the monotonicity of the rate function I(µ), and the fact that the objective function in (71) does not increase in µ j beyond µ j = 1, we may also assume without loss of generality that µ * j ≤ 1 for j = 1, • • • , n T . As in [START_REF] Jaldén | Sphere decoding complexity exponent for decoding full rate codes over the quasi-static mimo channel[END_REF] we proceed to define two events Ω 1 and Ω 2 which we will prove to be jointly sufficient so that, at layer k = 2qT , the sphere decoder visits close to ρ c r-ld (r) nodes. These are given by

Ω 1 {µ * j -2δ < µ j < µ * j -δ, j = 1, • • • , q 0 < µ j < δ, j = q + 1, • • • , n T }, (76)
for a given small δ > 0, and

Ω 2 {σ 1 (I T ⊗ V H p )G |p ≥ u}, (77)
for some given u > 0, where for p n Tq then G |p denotes the first 2pT columns of G, and where V p denotes the last 2p columns of V obtained by applying the singular value decomposition on H R , i.e., H R = UΣV H , where

Σ diag{σ 1 (H R ), • • • , σ 2n T (H R )} with σ 1 (H R ) ≤ • • • ≤ σ 2n T (H R )
and VV H = I. Hence, V H p corresponds to the 2p largest singular values of H R . Note also that by choosing δ sufficiently small, and using the fact that µ * i > 0 for i = 1, • • • , q, we may without loss of generality assume that Ω 1 implies that µ j > 0 for all j = 1, • • • , n T .

Modifying the approach in [1, Theorem 1] to account for MMSE preprocessing and unconstrained decoding, the lower bound on the number of nodes visited at layer k by the sphere decoder, is given by

N k ≥ k i=1 2ξ √ kσ i (R k ) - √ k + . (78)
In the following, and up until (84), we will work towards upper bounding σ i (R k ) so that we can then lower bound N k . Towards this let

M reg |p ρ 1 2 -rT κ HG |p α r I |p ∈ R 2(n R +n T )T ×2pT
contain the first 2pT columns of M reg from (56), and note that

(M reg |p ) H M reg |p = ρ 1-2rT κ G H |p H H HG |p + α 2 r I ,
and that from (51) we get

(M reg |p ) H M reg |p = ρ 1-2rT κ G H |p (I T ⊗ H H R H R )G |p + α 2 r I.
Since

H H R H R = V(diag{σ 1 (H H R H R ), • • • , σ 2n T (H H R H R )})V H = V(diag{σ 1 (H H R H R ), • • • , σ 2n T (H H R H R )} -σ (2q+1) (H H R H R )diag{0, • • • , 0 2q , 1, • • • , 1 2p })V H + σ (2q+1) (H H R H R )V(diag{0, • • • , 0 2q , 1, • • • , 1 2p })V H ,
we have that

H H R H R σ (2q+1) (H H R H R )V(diag{0, • • • , 0 2q , 1, • • • , 1 2p })V H = σ (2q+1) (H H R H R )V(diag{0, • • • , 0 2q , 1, • • • , 1 2p }) (diag{0, • • • , 0 2q , 1, • • • , 1 2p })V H = σ (2q+1) (H H R H R )V p V H p
where the last equality follows from the fact that V p contains the last 2p columns of V and where A B denotes that A-B is positive-semidefinite. Since σ i (H H H) ∈ R and since the Kronecker product induces singular value multiplicity, it follows that

(M reg |p ) H M reg |p ρ 1-2rT κ σ (2q+1) (H H R H R )G H |p (I T ⊗ V p V H p )G |p + α 2 r I.
With respect to the smallest singular value of (M reg |p ) H M reg |p we have

σ 1 ((M reg |p ) H M reg |p ) ≥ ρ 1-2rT κ σ (2q+1) (H H R H R ) • σ 1 G H |p (I T ⊗ V p V H p )G |p + α 2 r
and consequently, given that H R ∈ Ω 2 , we have that

σ 1 (M reg |p ) ≥ ρ -rT κ u 2 ρσ l2(2q+1) (H H C H C ) + 1 . = ρ -rT κ ρ 1 2 (1-µq+1) + ≥ ρ -rT κ + 1 2 (1-δ) + , (79)
where the first inequality follows from (77), the exponential equality follows from (59) and from the fact that u > 0 is fixed and independent of ρ, and the last inequality follows from (76). From (54) we have that

σ i (M reg ) ≤ ρ -rT κ (1 + ρ(σ κ (G)σ l 2T (i) (H C )) 2 ) . = ρ -rT κ + 1 2 (1-µ l 2T (i) ) + , i = 1, • • • , 2n T T, (80)
where the asymptotic equality follows from the fact that σ κ (G) is fixed and independent of ρ. Furthermore (76) gives that for i = 1, • • • , 2qT then

σ i (M reg ) . ≤ ρ -rT κ +δ+ 1 2 (1-µ * l 2T (i) ) + , (81)
where we have made use of the fact that µ

* j ≤ 1 for j = 1, • • • , n T .
Given that µ * j > 0 for j = 1, • • • , q, then for sufficiently small δ and for i = 1, • • • , 2qT , we have that

- rT κ + 1 2 (1 -δ) + ≥ - rT κ + δ + 1 2 (1 -µ * l 2T (i) ) + ,
which means that for sufficiently small δ, a comparison of (79) and (81) yields

σ i (M reg ) < σ 1 (M reg |p ), for i = 1, • • • , 2qT .
The above inequality allows us to apply Lemma 3 in [START_REF] Jaldén | Sphere decoding complexity exponent for decoding full rate codes over the quasi-static mimo channel[END_REF], which in turn gives that

σ i (R k ) ≤ σ κ (M reg ) σ 1 (M reg |p ) + 1 σ i (M reg ), (82)
for i = 1, • • • , 2qT . Setting i = κ in (80) upper bounds the maximum singular value of M reg as σ κ (M reg ) . ≤ ρ -rT κ + 1 2 (1-µn T ) + ≤ ρ 1 2 -rT κ , (83)
where the last inequality is due to the fact that µ j ≥ 0.

Consequently combining (83) and (79) gives that

σ κ (M reg ) σ 1 (M reg |p ) + 1 . ≤ ρ 1 2 δ ,
which together with (81) and (82) gives that

σ i (R k ) . ≤ ρ -rT κ + 3 2 δ+ 1 2 (1-µ * l 2T (i) ) + , i = 1, • • • , 2qT. (84)
Consequently, going back to (78), we have that

2ξ √ kσ i (R k ) - √ k + . ≥ ρ ( rT κ -3 2 δ-1 2 (1-µ * l 2T (i) ) + ) > 0 (85)
and furthermore for i = 1, • • • , 2qT , we have that rT κ -3 2 δ -1 2 (1µ * l 2T (i) ) + > 0 directly from definition of q and for sufficiently small δ. As a result, for k ≤ 2qT we have that

N k . ≥ k i=1 ρ ( rT κ -3 2 δ-1 2 (1-µ * l 2T (i) )) (86) = ρ k i=1 ( rT κ -1 2 (1-µ * l 2T (i) ) + )-3 2 kδ , (87)
and setting k = 2qT we have that

N 2qT . ≥ ρ ( 2qT i=1 ( rT κ -1 2 (1-µ * l 2T (i) ) + )-3qT δ) (88) = ρ (T q j=1 ( rT κ -(1-µ * j ) + )-3qT δ) (89) = ρ (c r-ld (r)-3qT δ) , (90)
where the last equality follows from (75). Consequently

N SD ≥ N 2qT . ≥ ρ c r-ld (r)-3qT δ ,
for small δ > 0. Given that δ can be chosen arbitrarily small, and given that events Ω 1 and Ω 2 occur, then the number of nodes visited by the SD at layer 2qT is arbitrarily close to the upper bound of ρ c r-ld (r) . Now to show that c r-ld (r) ≥ c r-ld (r) -3qT δ, we just have to prove thatlim ρ→∞ P N SD .

≥ ρ c r-ld (r)-3qT δ log ρ < d L (r). Toward this we note that as (76) and (77) imply that N SD .

≥ ρ c r-ld (r)-3qT δ , it follows that

P N SD . ≥ ρ c r-ld (r)-3qT δ ≥ P (Ω 1 ∩ Ω 2 ) = P (Ω 1 ) P (Ω 2 )
where the equality follows from the i.i. We now turn to [1, Lemma 2] and recall that for the layered codes assumed here, as well as for any full-rate design and some non-random fixed decoding ordering (corresponding to a permutation of the columns of G), there exists a unitary matrix

V ′ p such that rank (I T ⊗ (V ′ p ) H )G |p = 2pT i.e., that σ 1 (I T ⊗ (V ′ p ) H )G |p > 0.
However, by continuity of singular values [START_REF] Horn | Matrix Analysis[END_REF] it follows for sufficiently small u > 0 (cf.(77)) that P (Ω 2 ) > 0, which implies 15 that P (Ω 2 ) . = ρ 0 as Ω 2 is independent of ρ. This in turn implies that

P N SD . ≥ ρ c r-ld (r)-3qT δ . ≥ P (Ω 1 ) . (91)
With Ω 1 being an open set, we have that

-lim ρ→∞ P (Ω 1 ) log ρ ≤ inf µ∈Ω1 I(µ), = q j=1 (|n T -n R | + 2j -1)(µ * j -2δ), = d L (r) -2(|n T -n R | + q)qδ, < d L (r), (92)
where the above follows from the monotonicity of the rate function

I(µ) = n T j=1 (|n T -n R | + 2j -1)µ i + n R n T t 2 µ n T , evaluated at {µ * 1 -2δ • • • , µ * q -2δ, 0, • • • , 0} = arg inf µ∈Ω1 I(µ),
and 16 also follows from the fact that, by definition, I(µ * ) = d L (r).

Consequently from (91) we have that

-lim ρ→∞ P N SD . ≥ ρ c r-ld (r)-3qT δ log ρ < d L (r), (93)
and directly from the definition of the complexity exponent, we have that c r-ld (r) ≥ c r-ld (r) -3qT δ. As the bound holds 15 In light of the fact that event V ′ p has zero measure, what the continuity of eigenvalues guarantees is that we can construct a neighborhood of matrices around V ′ p which are full rank, and which have a non zero measure. We also note that the matrices V ′ p can be created recursively, starting from a single matrix V ′ n T . 16 Recall that parameter t was previously introduced as a parameter that regulates the near zero behavior of the random variable. and since the left hand side of (103) cannot be negative, and furthermore given that c is independent of ŝ, we conclude that c . ≤ ρ b . We will now proceed to lower bound y -M r ′ŝ 2 + α 2 r ′ ŝ 2 and then use (103) to lower bound r -R r ′ŝ .

Towards lower bounding y -M r ′ŝ 2 + α 2 r ′ ŝ 2 we draw from Theorem 1 in [START_REF] Jaldén | DMT optimality of LR-aided linear decoders for a general class of channels, lattice designs, and system models[END_REF] and we let B be the spherical region given by

B {d ∈ R κ | d 2 ≤ Γ 2 }
where the radius Γ > 0 is independent of ρ and is chosen so that This may clearly be done for arbitrary ζ > 0. We will in the following temporarily assume that ν r ′ +ζ ≥ 1 and prove that, together with w 2 ≤ ρ b , the two conditions are sufficient for λ( Rr ′ )

.

≥ ρ ζT κ to hold. In order to bound the metric for ŝ ∈ Z κ where ŝ = s, we note that ν r ′ +ζ ≥ 1 implies that ∀d ∈ ρ B ∩ Z κ , we have from (104) that

1 4 M r ′ (ŝ -s) 2 = 1 4 M r ′ d 2 ≥ ρ 2ζT κ . ( 105
)
As w 2 ≤ ρ b , it follows that We know from the exponential tail of the Gaussian distribution that P w 2 > ρ b . = ρ -∞ and from Lemma 1 in [START_REF] Jaldén | DMT optimality of LR-aided linear decoders for a general class of channels, lattice designs, and system models[END_REF] that P (ν r ′ +ζ < 1)

. ≤ ρ -d M L (r ′ +ζ) . Hence P σ min ( Rr ) for all r ≥ ǫ > 0.

The association with the singular values

σ 1 ( Rr,k ) ≤ • • • ≤ σ k ( Rr,k )
is made using the interlacing property of singular values of sub-matrices, which gives that

σ i ( Rr,k ) ≥ σ i ( Rr ), i ≤ k = 1, • • • , κ, (113)
and for k = 1, • • • , κ, that P σ min ( Rr,k )

.

< ρ

-ǫT κ .

≤ ρ -d M L (r-ǫ) .

Finally from the DMT optimality of the exact implementation of the regularized lattice decoder [START_REF] Gamal | Lattice coding and decoding achieve the optimal diversity-multiplexing tradeoff of MIMO channels[END_REF], [START_REF] Jaldén | DMT optimality of LR-aided linear decoders for a general class of channels, lattice designs, and system models[END_REF], we have that P σ min ( Rr,k )

.

< ρ

-ǫT κ .

≤ ρ -d L (r-ǫ) .

This proves Lemma 1.

APPENDIX D PROOF FOR LEMMA 2

For a search radius that grows as ξ = √ z log ρ . = ρ 0 , we first prove that P w ′′ 2 > ξ 2 . ≤ ρ -z ′ for z > z ′ > d L (r). Towards establishing the properties of the equivalent noise w ′′ (cf. ( 45)), we consider an equivalent representation of the MMSE-preprocessed lattice decoder and let (cf. [START_REF] Wübben | Nearmaximum-likelihood detection of MIMO systems using MMSE-based lattice reduction[END_REF])

QR = Q 1 Q 2 R = M α r I ∈ R (n+κ)×κ (114)
be the thin QR factorization of the modified channel matrix, where

Q 1 = R -1 M ∈ R n×κ , Q 2 = α r R -1 ∈ R κ×κ and
where R H R = M H M + α 2 r I. It then follows that for F = Q H 1 , the MMSE-preprocessed lattice decoder is equivalent to lattice decoding in the presence of channel R and noise

w ′ = -α 2 r R -H s + R -H M H w = -α r Q H 2 s + Q H 1 w. ( 115
)
Consequently we calculate where the last equality follows after choosing the search radius such that z > z ′ > d L (r). This proves Lemma 2.

P w ′ > ξ ≤ P -α r Q H 2 s + Q H 1 w > ξ

2 (

 2 cf. (17)), to the new slr-rld = arg min ŝ∈Z κ r -RTŝ 2 ,

g

  L (x) ≤ lim ρ→∞ P (ŝ r-lr-ld = s) P (ŝ r-ld = s) + lim ρ→∞ P (N ≥ ρ x ) P (ŝ r-ld = s) + lim ρ→∞ P ( w ′′ > ξ) P (ŝ r-ld = s) .

d 1 + 4 M r ′ d 2 ,

 142 d 2 ∈ R for any d 1 , d 2 ∈ B. The existence of the set B follows by the assumption that 0 is contained in the interior of R. Now let ν r ′ min d∈ρ r ′ T κ B∩Z κ :d =0 1 and for given γ > ζ > 0 choose b > 0 such that 2ζT κ > b > 0.

4 M r ′ +ζ d 2 κ

 42 (r ′ +ζ)T κ B ∩ Z κ , d = 0 it is the case that 1 (a) follows from the fact that M r = ρ , ∀d ∈ ρ (r ′ +ζ)T κ B ∩ Z κ , d = 0. (104)As R is bounded, and asζ > 0, it holds that R ⊂ 1 2 ρ ζT κ B for all ρ ≥ ρ 1 , for a sufficiently large ρ 1 . This implies that s ∈ 1 2 ρ (r ′ +ζ)T κ B for ρ ≥ ρ 1 since s ∈ ρ r ′ T κ R. For s, d ∈ 1 2 ρ (r ′ +ζ)T κ B ∩Z κ , there exists an ŝ ∈ ρ (r ′ +ζ)T κ B ∩Z κ , ŝ = s, such that ŝ = d+s. Hence for any ŝ ∈ ρ (r ′ +ζ)T κ

≤

  ρ b allows for (103) to give that r -R r ′ŝ 2 . and (102), we haveλ(R r ′ ) ≥ r -R r ′ŝw .where the exponential inequality follows from (110). Furthermore we know thatλ(R r ) = ρ -γT κ λ(R r ′ ) ǫ = γζ, r ≥ ǫ > 0,and from (100) and (112) it follows that σ min ( Rr ) now note that the above implies that for ν r ′ +ζ ≥ 1 andw 2 ≤ ρ b then σ min ( Rr ) ν r ′ +ζ < 1) ∪ ( w 2 > ρ b )≤ P (ν r ′ +ζ < 1) + P w 2 > ρ b .

.

  ≤ ρ -d M L (r-ǫ)

≤ 2 = P w 2 > z 1 κ 2

 22 κ w + κK > ξ) = P κ w > (z log ρ) P κ w > (z 1 log ρ)1 log ρ (d) = P w 2 > z 2 log ρ . = ρ -z2(116)where (a) follows from the MMSE-preprocessed equivalent channel representation (cf. (114)), and where the inequalities in (b), (c) and (d) follow for some fixed K that upper bounds sup s∈S κ r α r s , and for some arbitraryz 1 , z 2 satisfying z > z 1 > z 2 > 0 independent of ρ. Consequently P w ′′ > ξ = P QH w ′ > ξ . ≤ ρ -z ′for some 0 < z ′ < z 2 , and as a result lim ρ→∞ P w ′′ > ξ P (ŝ r-ld = s) = lim ρ→∞ ρ (d L (r)-z ′ ) = 0,

  where we recall that R is the upper triangular code-channel matrix, after MMSE preprocessing and LLL lattice reduction. Towards this, and for d L (rǫ) denoting the diversity gain of the exact implementation of the regularized lattice decoder at multiplexing gain rǫ, we have the following lemma on the smallest singular value of Rk . The proof appears in Appendix C.Lemma 1: The smallest singular value σ min ( Rk ) of submatrix Rk

  d. assumption on the entries in H C , which makes the singular values of H H C H C independent of the singular vectors of H H C H C [28], [29], and which in turn also implies independence of the singular values of H H C H C (event Ω 1 ) from the singular vectors of H H R H R (event Ω 2 ).

  1 4 M r ′ d 2 ≥ w 2 for large ρ, and that y -M r ′ŝ 2 = M r ′ (sŝ) + w Given (107) and (109), for any ŝ ∈ Z κ such that ŝ = s, it is the case that y -M r ′ŝ 2 +α 2 r ′ ŝ

				On the other hand if ŝ / ∈ ρ we have that α 2 r ′ ŝ 2 ≥ 1 4 Γ 2 ρ (r ′ +ζ)T κ 2ζT κ , and consequently that B, then by definition of B
				y -M r ′ŝ 2 + α 2 r ′ ŝ	2 ≥	1 4	Γ 2 ρ	2ζT κ .	(108)
				From (107) and (108) we then conclude that
				y -M r ′ŝ 2 + α 2 r ′ ŝ	2 . ≥ ρ	2ζT κ .	(109)
					2 .
	2 . ≥ ρ	2ζT κ .	(106)
	2 . ≥ ρ	2ζT κ .	(107)

Consequently

y -M r ′ŝ 2 + α 2 r ′ ŝ

We will interchangeably use MMSE-preprocessed decoder and regularized decoder, with the first term being more commonly used, and with the second implying a more general family of decoders (cf.[START_REF] Jaldén | DMT optimality of LR-aided linear decoders for a general class of channels, lattice designs, and system models[END_REF] where the equivalence between the two decoders is discussed.). Even though in the asymptotic setting of interest, the two accept the same results throughout the paper, some extra error-performance gains can be achieved by proper optimization of the regularized decoder (cf.[START_REF] Pan | A lagrangian dual relaxation approach to ML MIMO detection: reinterpreting regularized lattice decoding[END_REF]).

In other words, this complexity bound holds even if the channel statistics are such that the channel realizations cause the decoder to always have to solve the hardest possible lattice search problem.

The radius here is considered fixed in the sense that it does not vary with respect to the channel realization and rate.

The i.i.d. regular fading statistics satisfy the general set of conditions as described in[START_REF] Zhao | Diversity and multiplexing tradeoff in general fading channels[END_REF], where a) the near-zero behavior of the fading coefficients h is bounded in probability as c 1 |h| t ≤ p(h) ≤ c 2 |h| t for some positive and finite c 1 , c 2 and t, where b) the tail behavior of h is bounded in probability as p(h) ≤ c 2 e -b|h| β for some positive and finite c 2 , b and β, and where c) p(h) is upper bounded by a constant K.

A more proper statement would be that the QR decomposition is performed by the LR algorithm it self.

In general, (71) does not have a unique optimal point because (a) + is constant in a for a ≤ 0.

The research leading to these results has received funding from the European Research Council under the European Community's Seventh Framework Programme (FP7/2007-2013) / ERC grant agreement no. 228044, from the Swedish Foundation for Strategic Research (SSF) / grant ICA08-0046, from FP7/2007-2013 grant agreement no. 257616 (CONECT), and from the Mitsubishi Electric R&D Centre Europe project Home-eNodeBS. A. Singh and P.

for arbitrarily small δ > 0, it follows that c r-ld (r) = c r-ld (r). Directly from [START_REF] Jaldén | Sphere decoding complexity exponent for decoding full rate codes over the quasi-static mimo channel[END_REF]Theorem 4] which analyzes the ML-based complexity exponent c ml (r), together with the fact that the ML-based sphere decoder, with or without MMSE preprocessing, shares the same upper bound c r-ld (r) as the MMSEpreprocessed lattice decoder, gives that c ml (r) = c r-ld (r), which in turns implies that c r-ld (r) = c ml (r).

This establishes Theorem 1 and Corollary 1a.

APPENDIX B PROOF FOR COROLLARIES 1B AND 1C

Section A-A shows that c r-ld (r) can be obtained as the solution to the constrained maximization problem

In some cases though, further knowledge of the error performance of the encoder and decoder, can result in an explicit characterization of the complexity exponent. Take for instance the case of DMT optimal encoding [START_REF] Elia | Explicit space-time codes achieving the diversity-multiplexing gain tradeoff[END_REF], [START_REF] Oggier | Perfect space-time block codes[END_REF] and DMT optimal MMSE-preprocessed lattice decoding [START_REF] Gamal | Lattice coding and decoding achieve the optimal diversity-multiplexing tradeoff of MIMO channels[END_REF], [START_REF] Jaldén | DMT optimality of LR-aided linear decoders for a general class of channels, lattice designs, and system models[END_REF], where the constraint I(µ) ≤ d L (r) in (94a) reverts to the constraint

), which may be recognized to correspond to the no-outage region (cf. [START_REF] Zheng | Diversity and multiplexing: A fundamental tradeoff in multiple-antenna channels[END_REF]). In this case c r-ld (r) can then be explicitly obtained from the optimization problem

which can be solved in a straightforward manner to give that

describing the upper bound on the complexity exponent for MMSE-preprocessed lattice sphere decoding of DMT optimal full-rate codes, which for minimum delay (n T = T ) DMT optimal full-rate codes takes the form

and which further simplifies to Moving on to the universal upper bound, we can see from (71) that, regardless of the fading statistics and the corresponding I(µ), the exponent c r-ld (r) is non-decreasing in d L (r) and is hence maximized when d L (r) is itself maximized, i.e., it is maximized in the presence of DMT optimal encoding and decoding. Combined with the fact that the corresponding maximization problem in (95) does not depend on the fading distribution, other than the natural fact that its tail must vanish exponentially fast, results in the fact that, for any fullrate code and statistical characterization of the channel, the complexity of MMSE-preprocessed lattice SD is universally upper bounded as (cf. [START_REF] Jaldén | Sphere decoding complexity exponent for decoding full rate codes over the quasi-static mimo channel[END_REF])

This proves Corollary 1c.

) 17 , it follows by the bounded orthogonality defect of LLL reduced bases that there is a constant K κ > 0 independent of R r and ρ, for which (cf. [START_REF] Lenstra | Factoring polynomials with rational coefficients[END_REF] and the proof in [START_REF] Taherzadeh | LLL reduction achieves the receive diversity in MIMO decoding[END_REF])

where

denotes the shortest vector in the lattice generated by R r . As a result we have that

Looking to lower bound σ min ( Rr ), we seek a bound on λ(R r ). Towards this let r ′ = rγ for some r ≥ γ > 0, in which case for s being the transmitted symbol vector, and for any ŝ ∈ Z κ such that ŝ = s, it follows that

From (102) it is clear that to find a lower bound on λ(R r ′ ), we need to lower bound r -R r ′ŝ for all ŝ ∈ Z κ and upper bound w . Let us, for now, assume that w 2 ≤ ρ b . To lower bound r -R r ′ŝ , we draw from the equivalence of MMSE preprocessing and the regularized metric (cf. equation (45) in [START_REF] Jaldén | DMT optimality of LR-aided linear decoders for a general class of channels, lattice designs, and system models[END_REF]), and rewrite

where c y H [I -M H r ′ (M H r ′ M r ′ + α 2 r ′ I) -1 M r ′ ]y ≥ 0. We now note that for ŝ = s then y -M r ′ s 2 + α 2 r ′ s 2 .

≤ ρ b ,


Olivier Besson
email: olivier.besson@unine.ch

Martine Picq
email: martine.picq@insa-lyon.fr

Jerome Pousin
email: jerome.pousin@insa-lyon.fr

Jérôme Pousin

A Transport for imaging process

Keywords: AMS Classification 35F40, 35L85, 35R05, 62-99

teaching and research institutions in France or abroad, or from public or private research centers.

Introduction

Modern medical imaging modalities can provide a great amount of information to study the human anatomy and physiological functions in both space and time. In cardiac magnetic resonance imaging (MRI) for example, several slices can be acquired to cover the heart in 3D and at a collection of discrete time samples over the cardiac cycle. From these partial observations, the challenge is to extract the heart's dynamics from these input spatio-temporal data throughout the cardiac cycle [START_REF] Lynch | Segmentation of the Left Ventricle of the Heart in 3D+t MRI Data Using an Optimized Non-Rigid Temporal Model[END_REF], [START_REF] Schaerer | A New Dynamic Elastic Model for Cardiac Image Analysis[END_REF]. Image registration consists in estimating a transformation which insures the warping of one reference image onto another target image (supposed to present some similarity). Continuous transformations are privileged, the sequence of transformations during the estimation process is usually not much considered. Most important is the final resulting transformation and not the way one image will be transformed to the other. Here, we consider a reasonable interpolation process to continuously map the image intensity functions between two images in the context of cardiac motion estimation and modeling. The aim of this paper is to present, in the context of optical flow, an algorithm to compute a time dependent transportation plan without using lagrangian techniques. The paper is organized as follows. The introduction is ended, by recalling the optical flow model (OF) . In section 2, the algorithm is presented, and its convergence is discussed. In section 3 it is shown that the solutions obtained with the proposed algorithm are the solutions minimizing the same energy than the time dependent optimal mass transportation problem. Section 4 is devoted to numerical results. In particular a 2D cardiac medical image is considered.

The optical flow (OF) method

Let ρ be the intensity function, and v be the velocity of the apparent motion of brightness pattern. An image sequence is considered via the gray-value map ρ : Q = (0, 1) × Ω → R where Ω ⊂ R d is a bounded regular domain, the support of images, for d = 1, 2, 3. If image points move according to the velocity field v : Q → R d , then gray values ρ(t, X(t, x)) are constant along motion trajectories X(t, x). One obtains the optical flow equation.

d dt ρ(t, X(t, x)) = ∂ t ρ(t, X(t, x)) + ( v | ∇ X ρ(t, X(t, x)) ) R d = 0. ( 1
)
The previous equations lead to an ill-posed problem for the unknown (ρ, v).

Variational formulations or relaxed minimizing problems for computing jointly (ρ, v) have been first proposed in [START_REF] Aubert | Computing optical flow problem via variational techniques[END_REF] and after by many other authors. Here our concern is somewhat different. Finding (ρ, v) simultaneously is possible by solving a mass transport problem. Similarly to the work developed in [START_REF] Benamou | A computational fluid mechanics solution to the Monge-Kantorovich mass transfer problem[END_REF][START_REF] Benamou | The Monge-Kantorovich mass transfer and its computational fluid Mechanics formulation[END_REF], a characterization of (ρ, v) as solution of a minimizing problem is developed.

Let ρ 0 and ρ 1 be the cardiac images between two times arbitrary fixed to zero and one, the mathematical problem reads: find ρ the gray level function defined from Q with values in [0, 1] verifying

∂ t ρ(t, x) + ( (v(t, x) | ∇ρ(t, x) ) = 0, (t, x) ∈ (0, 1) × Ω ρ(0, x) = ρ 0 (x); ρ(1, x) = ρ 1 (x) x ∈ Ω (2)
The velocity function v, is determined in order to minimize the functional.

inf ρ,v 1 0 Ω ρ(t, x) v(t, x) 2 dtdx. (3)
Thus we get an image sequence through the gray-value map ρ. Let us mention [START_REF] Angenent | Minimizing Flows for the Monge-Kantorovich Problem[END_REF], for example, where the optimal mass transportation approach is used in images processing. In this work, the optimal transport problem in 2D is decomposed in several 1D optimal transport problems which are easier to numerically solved. The algorithm proposed here, is based on a least squares formulation for the transport equation ( [START_REF] Bochev | Least-Squares Finite Element Methods[END_REF] for example), which differs from the methods proposed in [START_REF] Angenent | Minimizing Flows for the Monge-Kantorovich Problem[END_REF], or in [START_REF] Benamou | A computational fluid mechanics solution to the Monge-Kantorovich mass transfer problem[END_REF]. Notice that the proposed transport in this paper differs from the optimal transportation, studied e.g. in the book of C. Villani [START_REF] Villani | Topics in optimal transportation[END_REF].

Algorithm for solving the optical flow

In this section, an algorithm is presented to solve the optical flow given by equations (2), and (3). Let us first specify our hypotheses.

H1 The domain Ω is a bounded C 2,α domain satisfying the exterior sphere condition.

H2 The functions ρ i ∈ C 1,α (Ω) for i = 0, 1, with ρ 0 = ρ 1 on ∂Ω. Moreover there exist two constants such that 0 3) is a velocity field vanishing on ∂Ω. To do so, the following method is used. Assume that

< β ≤ ρ i ≤ β in Ω. Let ρ 0 ∈ C 1,α ([0, 1] × Ω) be given by ρ 0 (t, x) = (1 -t)ρ 0 (x) + tρ 1 (x). We have ∂ t ρ 0 C 0,α ([0,1]×Ω) ≤ C(ρ 0 , ρ 1 ) and ∂ t ρ 0 | ∂Ω = 0. For each t ∈ [0, 1], our need to solve problem (2)-(
ρ n ∈ C 1,α ([0, 1] × Ω) is given with ρ n (0, x) = ρ 0 (x), and ρ n (1, x) = ρ 1 (x) for all x ∈ Ω. • Compute η n+1 such that    -div(ρ n (t, •)∇η n+1 ) = 0 in Ω ρ n (t, •) ∂η n+1 ∂n = 1 on ∂Ω, (4)
and set

C n+1 (t) = 1 |∂Ω| Ω ∂ t ρ n η n+1 dx. • For each t ∈ [0, 1] compute ϕ n+1 solution of -div(ρ n (t, •)∇ϕ n+1 ) = ∂ t ρ n (t, •), in Ω ϕ n+1 = C n+1 (t) on ∂Ω. (5)
• Set v n+1 = ∇ϕ n+1 .

• Compute ρ n+1 , L 2 -least squares solution of

∂ t ρ n+1 (t, x) + v n+1 (t, x) | ∇ρ n+1 (t, x) = 0, (t, x) ∈ (0, 1) × Ω, ρ n+1 (0, x) = ρ 0 (x); ρ n+1 (1, x) = ρ 1 (x) x ∈ Ω. ( 6
)
Remark 2.1. If the requirement ρ 0 ∂Ω = ρ 1 ∂Ω is canceled in hypothesis H2, then the boundary condition of problem ( 5) is replaced by ∂ϕ n+1 ∂n = 0 and we do not need anymore the constants C n+1 .

For each t ∈ [0, 1], since ρ n (t, •), and ∂ t ρ n (t, •) ∈ C 0,α (Ω), theorem 6.14 p. 107 of [START_REF] Gilbarg | Elliptic Partial Differential Equations of Second Order[END_REF] applies, and there exists a unique ϕ n+1 (t, •) ∈ C 2,α (Ω) solution of problem [START_REF] Benamou | A computational fluid mechanics solution to the Monge-Kantorovich mass transfer problem[END_REF]. In problem (5) the time t is a parameter. Since ρ n ∈ C 1,α , and ∂ t ρ n ∈ C 0,α , and C n+1 ∈ C 0,α , the classical C 2,α (Ω) a priori estimates for solutions to elliptic problems allow us to show that ϕ n+1 is a C 0,α function with respect to time. So we have

ϕ n+1 C 0,α ([0,1];C 2,α (Ω)) ≤ M ( C n C 0,α ([0,1]) + ∂ t ρ n C 0,α ([0,1]×Ω) ).
Consider the extension of ϕ n+1 by C n+1 outside of the domain Ω, still denoted byϕ n+1 . Since the right hand side of equation ( 5) vanishes on ∂Ω, this extension is regular, so the function v n+1 vanish outside Ω and belongs to C 0,α ([0, 1];

C 1,α (R 2 )). Define the flow X n+1 + (s, t, x) ∈ C 1,α ([0, 1] × [0, 1] × R 2 ; R 2 ) by d ds X n+1 + (s, t, x) = +v n+1 (s, X n+1 + (s, t, x)) in (0, 1) X n+1 + (t, t, x) = x. (7)
This flows will be constant for x ∈ ∂Ω. Observe that solving the transport equation ( 6) in the least-squares sense is equivalent to solve this equation on each integral curves defined by [START_REF] Besson | Solutions for linear conservation laws with velocity fields in L ∞[END_REF].

Set r(s) = ρ n+1 (s, X n+1 + (s, t, x)), and express equation ( 6) along the integral curves of equation [START_REF] Besson | Solutions for linear conservation laws with velocity fields in L ∞[END_REF]. The equation is reduced to the following ordinary differential equation with initial and final conditions.

d ds r(s) = 0 r(0) = ρ 0 (X n+1 + (0, t, x)); r(1) = ρ 1 (X n+1 + (1, t, x)). (8)
The L 2 least squares solution of ( 8) minimizes d ds r(s), and is given by:

r(s) = (1 -s)ρ 0 (X n+1 + (0, t, x)) + sρ 1 (X n+1 + (1, t, x)).
Therefore the following representation formula for the function ρ n+1 is proved.

Lemma 2.2. The L 2 -least squares solution of problem ( 6) is given by

ρ n+1 (t, x) = (1 -t)ρ 0 (X n+1 + (0, t, x)) + tρ 1 (X n+1 + (1, t, x)). (9)
Remark that the regularity of the function ρ n+1 is a consequence of the regularity of the flow X n+1 + . Let us now consider the convergence of the algorithm ( 4)-( 6).

Theorem 2.3. There exist

(ρ, ϕ) ∈ C 1 ([0, 1] × Ω; R * + ) × C 0 ([0, 1]; C 2 (Ω)), L 2 -least squares solution, respectively solution of ∂ t ρ + ( ∇ϕ | ∇ρ ) = 0, in (0, 1) × Ω ρ(0, x) = ρ 0 (x); ρ(1, x) = ρ 1 (x) in Ω (10) -div(ρ(t, •)∇ϕ) = ∂ t ρ(t, •), in Ω ϕ = C(t) and ∇ϕ = 0 on ∂Ω (11) with C(t) defined as follow.        -div(ρ(t, •)∇η) = 0 in Ω ρ(t, •) ∂ n η = 1 on ∂Ω C = 1 |∂Ω| Ω ∂ t ρ η dx. ( 12
) Proof. Since v 0 C 0,α ([0,1]) + ∂ t ρ 0 C 0,α ([0,1]×Ω) is bounded, ϕ n+1 C 0,α ([0,1];C 2,α (Ω)) and v n+1 C 0,α ([0,1];C 1,α (R 2 )
) are uniformly bounded in n. From lemma 2.2 there exists a unique ρ n+1 , the L 2 -least squares solution of [START_REF] Benamou | The Monge-Kantorovich mass transfer and its computational fluid Mechanics formulation[END_REF]. Let us give an estimate for D 3 X n+1 + . Starting from

D 1 X n+1 + (s, t, x)) = v n+1 (s, X n+1 + (s, t, x)),
we deduce (see [START_REF] Ambrosio | Transport equation and Cauchy problem for BV vector fields[END_REF])

D 3 D 1 X n+1 + (s, t, x) = D 2 v n+1 (s, X n+1 + (s, t, x))D 3 X n+1 + (s, t, x) D 3 X n+1 + (t, t, x) = Id. ( 13
) Since D 3 D 1 X n+1 + (s, t, x) = D 1 D 3 X n+1 + (s, t, x) we get D 3 X n+1 + (s, t, x) = e -s t D 2 (v n+1 (τ,X n+1 + (τ,t,x))) dτ Id. ( 14
)
Thus

D 3 v n+1 + C 0,α ([0,1] 2 ×R 2 ) is uniformly bounded in n. Moreover we have [1] D 2 X n+1 + (s, t, x) = v n+1 (s, t, x) | D 3 X n+1 + (s, t, x) so D 2 v n+1 C 0,α ([0,1] 2 ×R2) is bounded independently of n. From theorem 2.2 we deduce that ρ n+1 C 1,α ([0,1]×Ω) is uniformly bounded. Since the embeddings C 0,α ([0, 1]; C 2,α (Ω)) ֒→ C 0 ([0, 1]; C 2 (Ω)) and C 1,α ([0, 1]×Ω) ֒→ C 1 ([0, 1]×Ω)
are relatively compact there is a subsequence of (ρ n , ϕ n ) solution of ( 4)-( 6), still denoted by and(ρ, ϕ) is the solution of ( 10)-( 12) provided the boundary conditions to be justified. The condition ∇ϕ n | ∂Ω = 0 is valid for the approximations ϕ n (since the functions can be extended by C n outside of Ω). So the convergence in C 0 ([0, 1]; C 2 (Ω)) yields the condition for the gradient of limit function. For the approximations of function ρ, the formula given in lemma 2.2 combined with the regularity result show that the boundary conditions are exactly satisfied. These conditions are thus valid for the limit function due to the convergence in C 1 .

(ρ n , ϕ n ) converging to (ρ, ϕ) in C 1 ([0, 1]×Ω)×C 0 ([0, 1]; C 2 (Ω)),
We will show in the next section that the solution minimizes the same energy as for the time dependent optimal transportation mass.

Interpretation of solutions to problem (10)-(12)

In this section it is shown that the solution to problem (10)-( 12) is a solution to a time dependent mass transportation problem. Zero is a bound from below of the following functional

1 4 1 0 ∂ t u+div(u∇(ψ- C)) 2
L 2 (Ω) dt to be minimized with respect to (u, ψ).

0 = 1 4 1 0 ∂ t ρ + div(ρ∇(ϕ -C)) 2 L 2 (Ω) dt
thus (ρ, ϕ -C) solution of ( 10)-( 12) minimizes Min

{ψ∈L 2 ((0,1);H 1 0 (Ω)), u∈L 2 ((0,1);L 2 (Ω)) ∂tu+( u | ∇ψ )=0 u(0)=ρ 0 ; u(1)=ρ 1 in Ω} 1 4 1 0 ∂ t u + div(u∇ψ) 2 L 2 (Ω) dt.
So the solution (ρ, ϕ) of problem ( 10)-( 12), satisfies

(ρ, ϕ -C) = Argmin {ψ∈L 2 ((0,1);H 1 0 (Ω)), u∈L 2 ((0,1);L 2 (Ω)), ∂tu+( u | ∇ψ )=0, u(0)=ρ 0 ; u(1)=ρ 1 in Ω} 1 4 1 0 ∂ t u + div(u∇ψ) H -1 dt. ( 15
)
Theorem 3.1. Let (ρ, ϕ) be the solution of problem (10)-( 12) given in theorem 2.3, then it satisfies

(ρ, ∇ϕ) = Argmin {v∈L 2 ((0,1);(H 1 (Ω)) 2 ), u∈L 2 ((0,1);L 2 (Ω)), ∂tu+div(uv)=0, ∂tu+( v | ∇u )=0, u(0)=ρ 0 ; u(1)=ρ 1 in Ω} 1 0 Ω u v 2 dxdt. (16)
Proof. Observe that in problem [START_REF] Schaerer | A New Dynamic Elastic Model for Cardiac Image Analysis[END_REF], the transport equation is solved with a L 2 -least square procedure. This is equivalent to find ξ = u -(1 -t)ρ 0 + tρ 1 such that

∂ t ξ + ( v | ∇ξ ) = P R (ρ 1 -ρ 0 + ( v | ∇ ((1 -t)ρ 0 + tρ 1 ) )
where P R is the L 2 projection onto the range of the transport operator for the velocity v.

For u ∈ L ∞ (Ω), the expression v → div(uv) is well defined as a linear continuous on

H 1 0 (Ω). If moreover 0 < β ≤ u ≤ β in Ω, let H = H 1 0 (Ω) be equipped with the following inner product. (θ, ψ) = Ω u ( ∇θ | ∇ψ ) dx,
which induces a semi-norm which is equivalent to the H 1 -norm. The Riesz theorem claims that for the linear continuous form

L v (ψ) =< -div (uv), ψ > H;H ′ , there is a unique θ ∈ H such that L v (ψ) = Ω u ( ∇θ | ∇ψ ) dx, ∀ψ ∈ H.
Therefore v = ∇θ for a θ ∈ H 1 0 (Ω) and problem ( 16) is reduced to

(ρ, ∇ϕ) = Argmin {ψ∈L 2 ((0,1);H 1 0 (Ω)∩H 2 (Ω)), u∈L 2 ((0,1);L 2 (Ω)), ∂tu+div(u∇ψ)=0, ∂tu+( ∇ψ | u )=0, u(0)=ρ 0 ; u(1)=ρ 1 in Ω} 1 0 Ω u ∇ψ 2 dxdt. (17)
or

(ρ, ∇ϕ) = Argmin {-ψ∈L 2 ((0,1);H 1 0 (Ω)∩H 2 (Ω)), u∈L 2 ((0,1);L 2 (Ω)), ∂tu-div(u∇ψ)=0, ∂tu-( u | ∇ψ )=0, u(0)=ρ 0 ; u(1)=ρ 1 in Ω} 1 0 Ω u ∇ψ 2 dxdt. (18)
Since

Ω u ∇ψ 2 dx = -div(u∇ψ) 2 H -1 , problem (17) reads:
(ρ, ∇ϕ) = Argmin

{-ψ∈L 2 ((0,1);H 1 0 (Ω)∩H 2 (Ω)), u∈L 2 ((0,1);L 2 (Ω)), ∂tu-div(u∇ψ)=0, ∂tu-( u | ∇ψ )=0, u(0)=ρ 0 ; u(1)=ρ 1 in Ω} 1 0 div(-u∇ψ) 2 H -1 dt. (19)
From the definition of the linear form L, observe that

1 4 div(-u∇ψ) + ∂ t u 2 H -1 = div(-u∇ψ) 2 H -1 , so problem (19) reads: (ρ, ∇ϕ) = Argmin {ψ∈L 2 ((0,1);H 1 0 (Ω)∩H 2 (Ω)), u∈L 2 ((0,1);L 2 (Ω)), ∂tu+div(u∇ψ)=0, ∂tu+( u | ∇ψ )=0, u(0)=ρ 0 ; u(1)=ρ 1 in Ω} 1 4 1 0 div(u∇ψ) + ∂ t u 2 H -1 dt. ( 20
)
Gathering [START_REF] Picq | Résolution de l'équation du transport sous contraintes[END_REF] with the previous result proves the theorem.

Remark 3.2. The proposed transport is a one which minimizes the divergence of the velocity in a weighted dual norm of H 1 0 Remark 3.3. The Dacorogna-Moser transport [START_REF] Dacorogna | On a partial differential equation involving the Jacobian determinant[END_REF]

ρ(t, x) = (1 -t)ρ 0 (x) + tρ 1 (x) with -∆ϕ = ∂ t ρ in Ω; ∂ n ϕ = 0; on ∂Ω; (21)
and the velocity v = ∇ϕ ρ satisfies (ρ, ∇ϕ) = Argmin {v∈L 2 ((0,1);(H 1 (Ω))

2 ), u∈L 2 ((0,1);L 2 (Ω)),

∂tu+div(v)=0, u(0)=ρ 0 ; u(1)=ρ 1 in Ω} 1 0 Ω v 2 dxdt. (22)
Indeed, as before, consider the space H = H 1 0 equipped with the previous semi-norm, and set

L v (ψ) =< -div (v), ψ > H;H ′ =< ∂ t u, ψ > H;H ′ .
Then the problem is reduced to

(ρ, ∇ϕ) = Argmin {ψ∈L 2 ((0,1);H 1 0 (Ω)∩H 2 (Ω)), u∈L 2 ((0,1);L 2 (Ω)) ∂tu+div(∇ψ)=0 1 4 1 0 ∂ t u -div(∇ψ) 2 H -1 dt. ( 23
)
Using the relation ∂ t u = -div(∇ψ) and Jensen's inequality we get

1 0 div(∇ψ) 2 H -1 dt = 1 0 ∂ t u 2 H -1 dt ≥ 1 0 ∂ t u dt 2 H -1 = ρ 1 -ρ 0 2 H -1 .
Thus the Dacorogna-Moser transport is a minimum of the functional (22), so it is a transport minimizing the divergence of the velocity in H -1 -norm.

Numerical Approximation of the 2D Optimal Extended Optical Flow

The numerical method is based on a finite element time-space L 2 least squares formulation (see [START_REF] Besson | Solutions for linear conservation laws with velocity fields in L ∞[END_REF]) of the transport problem [START_REF] Benamou | The Monge-Kantorovich mass transfer and its computational fluid Mechanics formulation[END_REF]. Define v n+1 as

v n+1 = (1, v n+1 1 , v n+1 2 ) t
and for a sufficiently regular function ϕ defined on Q, set

∇ϕ = ∂ϕ ∂t , ∂ϕ ∂x 1 , ∂ϕ ∂x 2 t , and
div( v n+1 ϕ) = ∂ϕ ∂t + 2 i=1 ∂ ∂x i (v n+1 i ϕ).
Let {ϕ 1 • • • ϕ N } be a basis of a space-time finite element subspace

V h = {ϕ, piecewise regular polynomial functions, with ϕ(0, •) = ϕ(1, •) = 0},
for example, a brick Lagrange finite element of order one ( [START_REF] Besson | Space-time integrated least squares: a time marching approach[END_REF]). Let Π h be the Lagrange interpolation operator. Let also W h be the finite element subspace of H 1 0 (Ω), where the basis functions {ψ 1 • • • ψ M } are the traces at t = 0 of basis functions {ϕ i } N

i=1 . An approximation of problem ( 5) is the following. For a discrete sequence of time t compute

Ω (ρ n h (t, •) ∇(ϕ n+1 h -C n (t)) | ∇ψ h dx = Ω ∂ t ρ n h (t, •)ψ h dx ∀ψ h ∈ W h , (24)
and define v n+1 = ∇ϕ n+1 h . The L 2 least squares formulation of problem ( 6) is defined in the following way. Consider the functional

J(c h ) = 1 2 Q v n+1 | ∇c h + ∇ Π h (1 -t)ρ 0 + tρ 1 2 dx dt.
This functional is convex and coercive in the appropriate anisotropic Sobolev's space

H = {ϕ ∈ L 2 (Q); v n+1 | ∇ϕ ∈ L 2 (Q); ϕ(0, •) = ϕ(1, •) = 0} since the velocity field v n+1 is regular enough. Moreover, v n+1 | ∇ϕ L 2 (Q)
is a norm in H (see [START_REF] Besson | Solutions for linear conservation laws with velocity fields in L ∞[END_REF]). Set

D = β, β , θ h = Π h (1 -t)ρ 0 + tρ 1 ,
and introduce

K h = {ϕ h ∈ V h ; ϕ h + θ h ∈ D}
it is a closed convex subspace. Thus an approximation of ( 6) subject to the constraint ρ n+1 ∈ D is defined with the following minimization problem:

min c∈K h J(c). ( 25
)
The minimizer of problem ( 25) is

c h = ρ n+1 h -θ h which is characterized by the Fermat's rule [2] Q v n+1 | ∇c h v n+1 | ∇ψ h dx dt ≥ Q -v n+1 | ∇θ h v n+1 | ∇ψ h dx dt (26)
for all ψ h ∈ T K h , the contingent cone to K h . Thus an approximation of the solution to problem ( 6) is

ρ n+1 h = c h + Π h (1 -t)ρ 0 + tρ 1 ∈ V h . Remark 4.1. Define the bilinear form a n+1 v (•, •) on V h × V h by a n+1 v (ϕ, ψ) = Q v n+1 | ∇ψ v n+1 | ∇ψ h dx dt.
Then the function c h is the solution of

a n+1 v (c h , ψ) = a n+1 v (-θ h , ψ)
as long as c h is everywhere non negative. Moreover, if V + h consists of the non negative functions in V h , the function c h is the orthogonal projection of θ h onto V + h , according to the inner product induced by the bilinear form a n+1 v .

When the approximation with finite element of the algorithm proposed in section 2 has converged, the computed solution ρ h can be decomposed as

ρ h = θ h + P av V + h θ h ,
that is to say the approximated Dacoragna-Moser transport solution augmented with its orthogonal projection onto non negative functions space according to the inner product a v .

In the next subsections, all the computations are done according to remark 2.1. So it is assumed that the normal velocity satisfy v n+1 n = ∂ϕ n+1 ∂n = 0.

The transport of a bump

As a first example, the displacement of a bump is considered. More precisely, let Ω =]-1, 1[×]-1, 1[, y 0 = 0.5, r 0 = 0.3, β > 0, and α > 0. For (x, y) ∈ Ω, set r = x 2 + (y -y 0 ) 2 , and define

ρ 0 = β + αe -r 2 /(r 2 0 -r 2 ) if r 2 0 > r 2 β elsewhere. ( 27
)
Then for r = x 2 + (y + y 0 ) 2 , define

ρ 1 = β + αe -r 2 /(r 2 0 -r 2 ) if r 2 0 > r 2 β elsewhere, (28)
The domain Ω is subdivided into 40 × 40 elements, and the time interval ]0, T [=]0, 1[ is subdivided into 60 elements, thus the linear system has 102541 unknowns. Is is solved with a preconditioned conjugate gradient.

In figures 1-5, the shape of the bump is presented at the time-steps 0, 12, 24, 36, 48, and 60, the first one is the initial shape, and the last one is the final shape. The next figure 1 describe the bump transport for β = 1. Remark that in this case this transport is very similar to the Dacorogna-Moser one. Remark that when β is small, the linear system becomes very ill-conditioned. Indeed for β ≤ 0.1 usual preconditioners like the classical IC0 one are useless. In this case the parallel Gram-Schmidt least squares preconditioner (DIAG + LS CGS OPT) developed in [START_REF] Straubhaar | Préconditionnement de systèmes linéaires symétriques définis positifs[END_REF][START_REF] Straubhaar | Preconditioners for the conjugate gradient algorithm using Gram-Schmidt and least squares methods[END_REF][START_REF] Straubhaar | Parallel preconditioners for the conjugate gradient algorithm using Gram-Schmidt and least squares methods[END_REF] is used. Then the linear system is solved in parallel using 48 processors. The global tolerance for the iterative scheme developed in section 2 is set to 0.01. The number of iterations, and the CPU time for each value of β is given in table

The left ventricle motion

The iterative strategy described in Section 2 is then used to compute an approximated solution, and to reconstruct the systole to diastole images of a slice of a left ventricle. Ten time steps have been used to compute the solution, and 10000 degrees of freedom for the time-space least squares finite element. The approximated fixed point algorithm converges in about 10 iterations with an accuracy of about 10 -4 . In this case the usual IC0 preconditioner is sufficient; this is essentially due to the fact that there is no large region in the domain with a very low density ρ. In the next figure 6, the initial image and the final image are presented. To summarize, in this work, we present a fixed point algorithm for the computation of the time dependent optimal mass transportation problem, allowing to handle the images tracking problem. The efficiency of the method has been tested with some 2D examples.

Figure 1 :

 1 Figure 1: Bump transport for β = 1

Figure 2 :Figure 3 :Figure 4 :Figure 5 :

 2345 Figure 2: Bump transport for β = 0.5

Figure 6 :

 6 Figure 6: End of diastole of a left ventricular (a), of systole (b)

Figure 7 :

 7 Figure 7: Time step 3 and 6

Table 1 :

 1 4.1. For β < 0.05 the algorithm did not converge. CPU time and number of iterations

	β	nb. iterations CPU time [s]	residue
	1	4	128	8.62 10 -4
	0.5	5	149	9.23 10 -3
	0.2	8	201	4.41 10 -3
	0.1	10	229	9.47 10 -3
	0.05	55	1019	7.77 10 -3