Joint source and relay design for two-hop amplify-and-forward relay networks with QoS constraints

In this paper, we consider the joint design of source precoding matrix and the relay precoding matrix in a two-hop multiple-input multiple-output relay network. The goal is to find a pair of matrices in order to minimize the power consumption and at the same time meet pre-selected quality of service constraints that are defined as the mean square error of each data stream. Using majorization theory, we simplify the matrix-valued optimization problem into a scalar-valued one. We then propose a lower bound and an upper bound of the original problem, both in convex forms. Specifically, the latter is solved by a multi-level water-filling algorithm that is much efficient than directly applying the interior point method. Numerical examples corroborate the proposed studies and also demonstrate the tightness of both bounds to the original problem.


Introduction
Relay networks have recently attracted much attention because of their promising capability in achieving reliable communication and wide coverage for the next generation of wireless systems [1,2]. Different types of relaying strategies, e.g., amplify and forward (AF), decode and forward (DF), and compress and forward (CF) were introduced in [3][4][5], respectively. DF and CF decode data before retransmission and are, thus, also known as regenerative strategies; AF relay amplifies the received data only and is known as non-regenerative strategy. The computational simplicity in AF relay makes it highly attractive and a strong candidate for the real-time application. On the other side, multi-input multi-output (MIMO) technique [6] can enhance the data transmission rate by introducing the spatial diversity gain. Therefore, combining relaying and MIMO becomes a natural way to further advance the future wireless communication systems.
Most of the existing works focus on the AF MIMO relay networks, where a certain performance criterion is optimized subject to power constraints at both source and *Correspondence: feifeigao@ieee.org 2 Tsinghua National Laboratory for Information Science and Technology, Beijing, 100084, China Full list of author information is available at the end of the article relay. For example, the mutual information and the total mean square error (MSE) criteria are selected as objective functions in [7,8], respectively. The authors in [9] investigated the diversity multiplexing tradeoff for MIMO relays. Moreover, there are some works on beamforming design for special types of AF MIMO relays; for instance, the authors in [10] considered codebook design for halfduplex AF MIMO relays, while the author in [11] considered beamforming design for rather vast types of MIMO relays. Applying the majorization theory [12], the authors in [13] proposed a unified framework for most performance criteria. The extension of [13] to the multiple relay case was introduced in [14].
All of the mentioned methods above considered enhancing the system performance by maximizing or minimizing a certain objective function constraint on power consumption at some or all nodes. Although this model improves the system performance, it does not guarantee a certain quality of service (QoS) requirement for an individual user. The importance of considering QoS becomes more obvious in practical applications supporting several types of service each with a different reliability requirement. One of the pioneering works considering QoS in MIMO point-to-point systems is [15], where the http://jwcn.eurasipjournals.com/content/2013/1/108 authors optimized the transmission power subject to predefined sets of QoS, e.g., individual MSE, signal-to-noise ratio (SNR), and bit error rate. In [16], the authors investigated the optimization of source beamforming and relay weighting matrix in order to minimize the total power subject to a given set of QoS for multi-input singleoutput broadcast channel. On the aspect of relay network, AF MIMO relay network with QoS constraints has been investigated in [17]. Applying majorization theory, the author in [17] proposed a unified framework for the design of the optimal structure of the source precoding and relay amplifying matrices. The author applied a successive geometric programming method to obtain the optimal power loading among data streams. Unfortunately, the computational complexity of the proposed solution in [17] compromises its suitability for practical implementation. With similar assumptions, the authors in [18] considered a simplified version of the problem in which only the relay power is minimized. Then, the minimization is executed over a convex lower bound of the objective function. In a rather more general setup, the authors in [19] studied the joint relay and source power minimization and applied the majorization theory to reduce the problem to a scalar one. Then, using QoS convex relaxation, an upper bound and a lower bound on the optimal results are presented.
In this paper, building upon the results in [19], we take a specific look into a dual-hop a AF MIMO relay network. We first jointly design the source and the relay precoding matrices such that the overall transmission power is minimized subject to a given set of QoS constraints. Applying the majorization theory, we reduce the original matrixvalued problem to a scalar-valued one and then propose two new convex optimization problems whose objective values serve as the lower bound and the upper bound of the original problem. While both new problems can be handled by the existing convex optimization tools, e.g., CVX [20], we specifically design a multi-level water-filling (MLWF) algorithm to solve the upper bound problem that can further reduce the computational complexity. Compared with the successive geometric programming approach developed in [17], the MLWF algorithm does not require any optimization tool and thus is easier to implement for practical relay systems. Numerical results corroborate the proposed studies and clearly demonstrate the tightness of the proposed lower bound and upper bound, especially over low MSE region.
The rest of this paper is arranged as follows: Section 2 presents the system model and formulates the optimization problem in matrix form. In Section 3, the optimization is simplified to a scalar-valued problem from the majorization theory. Two suboptimal problems whose objectives serve as the upper bound and the lower bound of the original problems are derived in the same section.
In Section 4, the upper bound problem is solved from a multi-level water-filling algorithm coupled with decomposition methods. The simulation results are presented in Section 5, and conclusions are drawn in Section 6.

Notations
Vectors and matrices are boldface small and capital letters, respectively; the transpose, complex conjugate, Hermitian, inverse, and pseudo-inverse of A are denoted by A T , A * , A H , A −1 , and A † , respectively; a denotes the Euclidean norm of the vector a; diag{a} is the diagonal matrix with diagonal elements given by a, while diag{A} is a vector with entries taken from the diagonal elements of A ; I is the identity matrix; and E{·} is the statistical expectation. Moreover, basic notations and definitions of majorization theory can be found in Appendix 1.

System model
As shown in Figure 1, we consider a three-node multiantenna relay network that consists of a source node, a relay node, and a destination node, equipped with M, N, and P antennas, respectively. We assume that the direct link between the source and the destination is weak and, therefore, can be neglected. Denote the baseband MIMO channels between the source and relay as the N×M matrix H s , while that between the relay and destination as the P × N matrix H r . We further assume that the channel state information is perfectly known at all nodes. Suppose that L data streams, denoted by x, are precoded by an M × L precoding matrix B at source node. We require L ≤ min{M, N, P} so that the data streams can be detected with linear method at the destination. With an N × N precoding matrix F at relay, the received signal at destination is as follows: where v r and v d are additive white complex Gaussian noise at relay and destination, respectively, i.e., v r ∈ CN (0, σ 2 r I N ) and v d ∈ CN (0, σ 2 d I N ). Without loss of generality, we set σ 2 d = σ 2 r = 1. Since the correlation and power can be designed through B, the data streams from the source can be assumed independent from each other, i.e., E{xx H } = I.
We consider the minimum mean square error detection at destination with a P × L decoding matrix T. The estimated data can be expressed as follows: with the error matrix The MSE is defined as tr(C). For the given B and F, one can easily find the optimal T from the standard approach [21] and the corresponding error matrix is As in [15], the QoS measurement here is taken as the MSE of each individual data stream. Let us denote the QoS vector as ρ =[ ρ 1 , . . . , ρ L ] T , i.e., the MSE of the ith data stream is required to be smaller than ρ i . Note that ρ i < 1 is necessary to avoid trivial solutions. From (5), the QoS constraints are stated as follows: where '≤' denotes the element-wise operation if used between vectors. The average power consumed by the source is computed from while that spent by the relay is

Optimization problem
The goal is to find the optimal B and F to minimize the overall power spent by the source and relay and at the same time meet the QoS requirements. Define The optimization problem is then expressed as follows: Unfortunately, the problem is non-convex and cannot be solved in an efficient way.

Equivalent problem
Let us define a new optimization: where w means weak majorization, and the details can be found in Appendix 1.

Theorem 1. Problems (P1) and (P2) are equivalent.
Proof. The idea is to prove that for each feasible point of (P1), there is a corresponding feasible point in (P2) that yields the same objective value and vice versa.
where Z (I + M H RM) −1 . From Z ii ≤ ρ i , one can conclude that diag{Z} w ρ. We further know from Lemma 2 in Appendix 1 that λ{Z} w diag{Z}. Therefore, λ{Z} w http://jwcn.eurasipjournals.com/content/2013/1/108 ρ and we achieve the conclusion that for any feasible point (B, F) in (P1), there is always a corresponding feasible point (B, F) in (P2) that gives the same objective value.
Moreover, the objective function of (P1) with (B, F) is the same as the objective function of (P2).

Processing matrix variables
Based on Theorem 1, we can solve (P2) instead of (P1).
Define two new positive semi-definite matrices as well as their singular value decomposition (SVD) as follows: where U X and U Y are the N × L and P × L orthonormal matrices, respectively. Throughout this paper, we always sort the singular values and eigenvalues in an increasing order. Note that it is still possible that X or Y contains some zero entries.
The matrices H sB and H r F can be represented as follows: where V X can be any L × L unitary matrix, and V Y can be any N × L orthonormal matrices. These two matrices will be designed later to fulfill the optimality requirement. Let be the SVD of H s and H r respectively. We further express the singular matrices as U s [ U s,2 , U s,1 ] and U r [ U r,2 , U r,1 ], in which U s,1 and U r,1 contain L first columns.
One can rewrite (15) and (16) as follows: There are three cases to be discussed: • If N = M = L, thenB can be re-expressed from (18) asB • If M > N = L, (18) can be simplified as The solution forB is not unique, and the one which minimizes the objective function should be chosen. Basically, we need to solve Note that the second term in the objective function is not included in (21) since it will be taken cared by F. From Lemma 9 in [14], the optimal structure ofB is as follows: ×N) . In this case, (22) is still the optimal structure forB. Please refer to [14] for the detailed derivation.
Following the same procedure, the optimal structure of F can be computed as follows: We then proceed to the optimization over new variables U X , V X , U Y , V Y , X , and Y . Substituting (22) and (23) into (P2) yields the new objective function where the pseudo-inverses of † s and † r are We apply the inversion lemma [22], for the constraint and obtain where G is defined as the corresponding term. Substituting (22) and (23) into G yields Remark 1. Representing the problem in terms of these new variables decreases the dependency between the objective function and constraints and thus facilitates the http://jwcn.eurasipjournals.com/content/2013/1/108 optimization procedure. For example, U Y is only involved in objective function, but not in the constraints. Hence, we can adjust U Y to change the objective without affecting the constraint. Reversely, V X and V Y only affect the constraint, but not the objective function.

Remark 2.
On the other hand, the diagonality constraint in (P2) can be satisfied by only adjusting V X so that the other constraint and the objective functions will not be affected.

Simplification
In the following subsections, we derive the optimum structures for V X ,V Y , U X , and U Y , respectively.

Optimal V X and V Y
The functionality of the optimal V X is to make G diagonal. Moreover, if there is a V Y such that all components in diag{I − G} are minimized simultaneously, then this V Y must be optimal because it provides the largest possible freedom to optimize over the rest variables, i.e., U X , U Y , X , and Y . Applying (9.H.2 in [12]) on G, we obtain is the multiplication of the diagonal matrices in (27).
According to the definition of submajorization and supermajorization in Appendix 1, (28) can also be written as follows: and is equivalent to where 1 is an all one vector. Moreover, (31) can be reexpressed as follows: Since one can always find V X that makes G diagonal, there is λ{I − G} = diag{I − G}. Moreover, from Lemma 2, we can rewrite (32) as follows: which indicates that diag{I −Ĝ} is a simultaneous lower bound for all the elements in diag{I − G}. Therefore, G = G must hold at the optimal point, and this can be achieved if V Y = U X and V X = I.

Optimal U X and U Y
With the optimal V Y , the constraint is independent from both U Y and U X . Therefore, one can find the optimal structure of U X and U Y purely from We need the following matrix inequality (9.H.1.h in [12]) to proceed: Given two L × L positive semi-definite matrices A 1 and A 2 with eigenvalues λ l (A 1 ) and λ l (A 2 ) arranged in the increasing order, there is Then, the first part in (34) can be lower bounded by the following: Obviously, the minimum value can be achieved when U X = U s , since the diagonal elements of s and X are all arranged in increasing order. Similar discussion holds for the second term in (34), namely and the minimum is achieved when U Y = U r,1 . Substituting the optimal U X , U Y , V X , and V Y into the objective function and the constraints of (P2), we obtain the following: Define a i and b i as the ith diagonal entries of ( † s ) 2 and ( † r ) 2 , respectively. Further, define x i and y i as the ith diagonal entries of X and Y , respectively. Then, problem (P2) is converted to a scalar form: x i ≥ 0, y i ≥ 0, ∀i.
Unfortunately, the constraint (40) is non-convex, and the problem cannot be solved efficiently. We then propose http://jwcn.eurasipjournals.com/content/2013/1/108 the following two convex bounds for each summand on the left-hand side of (40): Replacing the corresponding term in (40) by the righthand side (RHS) in (41) or (42) while keeping the same objective function, we can obtain the upper or the lower bound of the original problem, respectively. Both the lower bound and the upper bound problems can be solved by the existing convex optimization tools based on the interior point method, e.g., CVX, [20,23]. http://jwcn.eurasipjournals.com/content/2013/1/108 Figure 3 Proposed master algorithm.

Algorithm for the upper bound problem
Power minimization is ensured by considering the minimization of the upper bound of the objective value. Therefore, we take a detailed look into the upper bound problem and target at a more efficient solution. The upper bound problem is restated as follows: It is hard to obtain the closed form solution even with the Karush-Kuhn-Tucker (KKT) conditions. Observing the symmetry, we can apply the primal decomposition method [23] to break down the problem into two simpler subproblems.

The decomposition method
Problem (43) can be decomposed into two parallel subproblems together with a master problem connecting as a bridge [23]. Defining a new auxiliary vector t = [ t 1 , t 2 , . . . , t N ] T , the subproblems are as follows: x i ≥ 0, ∀i http://jwcn.eurasipjournals.com/content/2013/1/108 In many cases, the subproblems do not have closedform solution. Therefore, searching the optimum t for the master problem can be iteratively done from the subgradient method. Note that the convergence of the decomposed problem to the global optimal point can be guaranteed since (43) is convex.

Solving the subproblems
The Lagrangian corresponding to x (t) is as follows: and the corresponding KKT conditions are the following: For simplicity, we defineμ k = N i=k μ i . Multiplying both sides of (46) by x k and combining with (49) give the following: After a straightforward calculation, we obtain the following: Based on (52), the KKT conditions are reduced to the following: Next, we propose an efficient algorithm that could find the solution that satisfies (53), (54), and (55) simultaneously. The standard water-filling algorithm embedded in step 3 can be found in Appendix 2. For an explicit demonstration, we present the algorithm in a diagram in Figure 2.

Algorithm 1 Multi-level water-filling algorithm
We next need to prove the optimality of the above MLWF. Since the problem is convex, the output of the algorithm is the optimal solution if and only if all of the KKT conditions are satisfied. The conditions (55) are simultaneously satisfied from the step 4 in the algorithm. Condition (53) is itself satisfied by the nature of the water-filling algorithm, that is the water levels are always non-negative. The above MLWF algorithm satisfies (54) as well, as seen from the following lemma. http://jwcn.eurasipjournals.com/content/2013/1/108 Lemma 1. Successive water levels achieved by the proposed MLWF algorithm are ordered decreasingly.
Proof. The algorithm has two loops: one is the inner loop 'steps 3 → 4 → 5 → 3' and the other is the outer loop 'steps 3 → 4 → 5 → 6 → 3. ' Each time when the inner loop finishes, one water levelμ(L, H) will be achieved. Then, we proceed to compute other water levels. In the inner loop, we first adopt the water level given by standard water-filling algorithm. Our hypothesis is that the aforementioned water level satisfies all of the conditions (55). Then, we check our hypothesis by searching whether there is a k = l 0 for which the corresponding condition in (55) is violated.
Meanwhile, the subproblem y (t) can be solved with the same MLWF, which will not be restated here.

Algorithm 2 The master algorithm
Input: The number of all positive eigenvalues (N ), all inverted eigenvalues 1: Choose a positive increasing vector for t, e.g., Set the initial objective value T (0) = 0; 2: Obtain the solution and water levels of x (t) and y (t) from MLWF algorithm; 3: Calculate the objective function in the mth iteration as This algorithm is presented by a diagram in Figure 3.

Simulation results
In this section, we numerically examine the performance of the proposed method. For all examples, we assume that the channel matrices (H r and H s ) have independent and identically distributed Gaussian entries with zero mean and variance 1. We consider 6 × 6 × 6 MIMO relay system, and 1000 Monte-Carlo runs are taken for average. The upper bound optimization will be solved both by the proposed decomposition algorithms and by the CVX convex optimization toolbox [20], while the lower-bounded optimization will only be solved by CVX. Moreover, we also compare the proposed solutions with a suboptimal method and a trivial method listed in the following: 1. Diagonalization method. We reduce the problem into a scalar problem using SVD of channel matrices. If we consider the structure of B = U s B and F = U r F V H s , we can reduce the problem into a simple scalar problem. In this method we basically choose the precoding matrices to match the channel matrices. The optimization problem of (10) is simplified to the following: subject to  where a i , b i , x i , and y i are the i th diagonal entries of ( † s ) 2 , ( † r ) 2 , B , and F , respectively. Note that the problem in (60) is a special case of (P3) in (39). 2. Naive method. We consider a solution that satisfies all of the constraints in (43) with equality. In this method, x and y are given by the following: In the first example, we take equal QoS requirements for all data streams. The total consumed powers versus MSE are depicted in Figure 4 for the upper bound with the proposed method, the upper bound with CVX, the lower bound, and the suboptimal diagonalization method with CVX, as well as the naive method. It can be observed that the performance of the decomposition method is exactly the same as the one obtained from CVX. We can also observe that the lower and the upper bounds are quite tight at high SNR when they meet each other. The numerical results from the diagonalization method depicts weaker performance than the proposed method with roughly 2 dB more power consumption. The naive method, on the other hand, has the poorest performance and requires more power than the proposed method (roughly 3 dB).    The performances with unequal QoS constraints are presented in Figure 5, where the QoS vector is chosen as ρ = ρ[ 1, 1.5, 2, 2.5, 3, 3] T . One may notice that all curves are about 2 dB higher than those in Figure 4. This is because the smallest QoS constraint drags down the overall performance and thus demands for more power. Once again, we find that the upper bound and the lower bound are quite close to each other especially at high SNR values. Meanwhile, the proposed method coincides with the upper bound from the CVX. Moreover, the advantage of the proposed method over the suboptimal diagonalization method and the naive method is fairly explicit.
We also measure the average running time for the proposed method and the upper bound with CVX, measured in MATLAB. The results are shown in Table 1.

Remark 3.
Note that resorting to simulation cannot serve as a rigorous measurement of the complexity analysis. Nevertheless, the large difference in the running time fairly indicates the much higher efficiency of the proposed decomposition method over the regular convex optimization tool.
The effect of the number of antennas on each node is also investigated here. We set the QoS requirements for all data streams as 0.2 and vary the number of antennas from 2 to 9. The results are displayed in Figure 6.

Conclusion
In this paper, we considered the joint design of source precoding matrix and relay precoding matrix in a two-hop AF MIMO relay network. We minimized the power consumption subject to a set of predefined QoS constraints of each data stream. Using matrix calculus and majorization theory, we simplified the original matrix-valued problem to a relatively simpler scalar one and proposed two bounding problems that are convex and can be solved efficiently. We specifically designed a primary decomposition method to solve the upper bound problem that has less complexity than directly applying the interior point method. Numerical examples are provided to corroborate the proposed studies.

Basics of majorization theory
Here, we briefly introduce the basics of majorization theory while the more comprehensive discussion can be found in [12]. represents the elements of x in an increasing order. Similarly, assume x [1] ≥ . . . ≥ x [n] represents the elements of x in a decreasing order. For any x, y ∈ IR n , x is majorized by y if and is denoted by x ≺ y or y x. Definition 3. [12] For any x, y ∈ IR n , x is weakly supermajorized by y if We denote this with x ≺ w y (or equivalently y w x). Also, x is weakly submajorized by y if We denote this with x ≺ w y (or equivalently y w x).
Lemma 2. (9.B.1 in [12]) Let M be an n × n Hermitian matrix with a vector diag{M} denoting its diagonal elements and let vector λ{M} contain its eigenvalues. Then, diag{M} ≺ diag{M} ≺ λ{M} wherediag{M} i = mean(diag{M}). Reversely, given vectors a and b with a ≺ b, then there exists an n × n Hermitian matrix M whose diagonal elements are a and eigenvalues are b. Lemma 3. (5.A.9 in [12]) For any a and b satisfying a ≺ w b, there exists a vector x such that a ≥ x.