Hindawi Publishing Corporation EURASIP Journal on Advances in Signal Processing Volume 2010, Article ID 640186, 13 pages doi:10.1155/2010/640186

Research Article

Joint Linear Processing for an Amplify-and-Forward MIMO Relay Channel with Imperfect Channel State Information

Batu K. Chalise1 and Luc Vandendorpe2

1 Center for Advanced Communications, Villanova University, Villanova, PA 19085, USA

2 Communication and Remote Sensing Laboratory, Université catholique de Louvain, Place du Levant, 2, 1348 Louvain la Neuve, Belgium

Correspondence should be addressed to Batu K. Chalise, batu.chalise@villanova.edu

Received 22 March 2010; Accepted 5 August 2010

Academic Editor: Kostas Berberidis

Copyright © 2010 B. K. Chalise and L. Vandendorpe. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

The problem of jointly optimizing the source precoder, relay transceiver, and destination equalizer has been considered in this paper for a multiple-input-multiple-output (MIMO) amplify-and-forward (AF) relay channel, where the channel estimates of all links are assumed to be imperfect. The considered joint optimization problem is nonconvex and does not offer closed-form solutions. However, it has been shown that the optimization of one variable when others are fixed is a convex optimization problem which can be efficiently solved using interior-point algorithms. In this context, an iterative technique with the guaranteed convergence has been proposed for the AF MIMO relay channel that includes the direct link. It has been also shown that, for the double-hop relay case without the receive-side antenna correlations in each hop, the global optimality can be confirmed since the structures of the source precoder, relay transceiver, and destination equalizer have closed forms and the remaining joint power allocation can be solved using Geometric Programming (GP) technique under high signal-to-noise ratio (SNR) approximation. In the latter case, the performance of the iterative technique and the GP method has been compared with simulations to ensure that the iterative approach gives reasonably good solutions with an acceptable complexity. Moreover, simulation results verify the robustness of the proposed design when compared to the nonrobust design that assumes estimated channels as actual channels.

1. Introduction

The application of relays for cooperative communications has received a lot of interest in recent years. It is well known that the channel impairments such as shadowing, multipath fading, distance-dependent path losses, and interference often degrade the link quality between the source and destination in a wireless network. If the link quality degrades severely, relays can be employed between the source and destination nodes for assisting the transmission of data from the source to destination [1]. In the literature, various types of cooperative communications such as amplify-and-forward (AF), decode-and-forward [1], coded-cooperation [2], and compress-and-forward [3] have been presented. In [4], the outage and ergodic capacities have been analyzed for a three-node network where one of the nodes relays the messages of another node towards the third one. Among several cooperation schemes [1-4], the AF scheme is more

attractive due to its simplicity since the relay simply forwards the signal and does not decode it. Recently, space-time coding strategies have been developed for relay networks [5]. In [6], the authors study distributed beamforming for a cooperative network which consists of a transmitter, a receiver, and an arbitrary number of relay nodes. The common things among aforementioned works are that the transmitter, receiver, and the relays are all single-antenna nodes and the channel state information (CSI) (either instantaneous or second-order statistics of the channel) is assumed to be error-free.

The performance of cooperative communications can be further enhanced by employing multiple-input-multiple-output (MIMO) relays [7]. The optimal designs of AF MIMO relays have been investigated in [8, 9] for point-to-point and in [10, 11] for point-to-multipoint communications assuming that the available CSI is perfect. The robust design of MIMO relay for multipoint-to-multipoint

communications has been solved in [12], where the sources and destinations are single antenna nodes. The optimal design of multiple AF MIMO relays in a point-to-point communication scenario has been considered in [13, 14] to minimize the mean-square error (MSE) and satisfy the quality of service (QoS) requirements. These works also assume perfect knowledge of CSI. Recently, the joint robust design of AF MIMO relay and destination equalizer has been investigated in [15] for a double-hop (without direct link) MIMO relay channel. To the best source of our knowledge, the joint optimization of the source precoder, MIMO relay, and the destination equalizer has not been considered in the literature for the case where the CSI is imperfect and the direct link is included. Although the path attenuation for the direct link is much larger than that for the link via relay, due to the fading of the wireless channels, there can be still a significant number of instantaneous channels during which the direct link is better than the relay link. As a result, we consider the direct link in our analysis and exploit the benefit provided by the relay channel in terms of diversity. Moreover, in practice, channel estimation is required to obtain the CSI, where the estimation errors are inevitable due to various factors such as the limited length of the training sequences and the time-varying nature of wireless channels. The performance degradation due to such estimation errors can be mitigated by using robust designs that take into account the possible estimation errors. As a result, robust methods are highly desired for practical applications. The robust techniques can be divided mainly into worst-case and stochastic approaches [16]. The worst-case approach [17, 18] considers that the errors belong to a predefined uncertainty region, where the objective is to optimize the worst system performance for any error in this region. The stochastic approach guarantees a certain system performance averaged over channel realizations [19]. The latter approach has been used in [20] to minimize the power of the transmit beamformer while satisfying the QoS requirements for all users. In the sequel, we use stochastic approach for the robust design.

In this paper, we deal with the joint robust design of source precoder, relay transceiver, and destination equalizer for an AF MIMO relay system where the CSI is considered to be imperfect at all nodes. A stochastic approach is employed in which the objective is to minimize the average sum mean-square error ( If the channel estimation is perfect at the receiver, the minimum mean-square error (MMSE) matrix X can be related to the rate using the relation r = - log det(X). However, if the receiver does not have perfect estimation of the channel, the relation between the rate and MMSE matrix is not straightforward. Consequently, for our current system model where both the estimates of source-relay and relay-destination channels are imperfect, deriving rate expression and solving the optimization problem based on that expression are still an open issue.) under the source and relay power constraints. The considered joint optimization is nonconvex and also does not lead to closed-form solutions. However, it has been shown that the optimization of one parameter when others are fixed is a quadratic convex-optimization problem that can be easily solved within the framework of

convex optimization techniques. We first propose an iterative approach both for the MIMO relay channels with and without the direct link. Although the iterative method guarantees fast convergence for the case with the direct link, the global optimality cannot be proven since the joint optimization problem is nonconvex. As a result, in the second part of this paper, we limit the joint optimization problem for the case without the direct link in which the source-relay and relay-destination MIMO channels have only transmit-side antenna correlations. In the latter case, it is shown that structures of the optimal source precoder and relay transceiver have closed forms, where the remaining joint power allocation problem can be approximately formulated into a Geometric Programming (GP) problem. With the help of computer simulations, we compare the solutions of the iterative technique and the GP approach under high signal-to-noise ratio (SNR) approximation for the case without the direct link. This comparison is helpful to conclude that the iterative approach gives reasonably good solutions with an acceptable complexity.

The remainder of this paper is organized as follows. The system model for MIMO relay channel is presented in Section 2. In Section 3, the iterative approach is described for jointly optimizing the source precoder, relay transceiver, and destination equalizer for the MIMO relay channel with the direct link. The closed-form solutions and the approximate GP problem formulation are provided in Section 4 for the MIMO relay channel without the direct link where singleside antenna correlations have been considered for source-relay and relay-destination channels. In Section 5, simulation results are presented to show the performance of the proposed robust and nonrobust methods, and in Section 6, conclusions are drawn.

Notations . Upper (lower) bold face letters will be used for matrices (vectors); (■)*,(■ )T, (• )H,E{-}, In, and || • || denote conjugate, transpose, Hermitian transpose, mathematical expectation, n X n identity matrix, and Frobenius norm, respectively. tr(-), vec(-), CMxM, 0 denote the matrix trace operator, vectorization operator, space of M X M matrices with complex entries, and the Kronecker product, respectively.

2. System Model

We consider a cooperative communication system that consists of a source, a relay, and a destination which are all multiantenna nodes. The block diagram is shown in Figure 1. Notice that the direct link between the source and destination is taken into account, so that the diversity order of the cooperative system can be maintained. The source has M antennas, the relay has NR receiving antennas and NT transmitting antennas, and the destination has ND antennas. The relay protocol consists of two timeslots. In the first timeslot, the source sends a symbol vector to the destination and relay. The relay linearly processes the source symbol vector and sends it to the destination in the second timeslot. The source remains idle during the second timeslot. At the

Figure 1: Cooperative MIMO relay channel.

end of two time slots, the destination linearly combines the symbol vectors received from the source and relay [1]. It is assumed that the estimates of the source-relay, relay-destination, and source-destination channels are available instead of their exact knowledge. The MIMO channels are considered to be spatially correlated block-fading frequency-flat Rayleigh channels. The signal received by the relay is given by

yr = H1Fs + nr,

where s e CNsX1 is the complex source signal of length NS, F e cmxns is the precoder that maps the NS X 1 symbol vector into M source antennas, H1 e cnrxm is the MIMO channel between the source and relay, and nr e CNrx1 is the additive Gaussian noise vector at the relay. We also assume that the elements of s are statistically independent with the zero-mean and unit variance, that is, E{ssH} = INs. The precoder F at the source operates under the power constraint PS = tr(FFH) < Pmax, where PSmsx is the maximum power of the source. We consider nr ~ nc(0, ffr2INR), that is, the entries of nr are zero-mean circularly symmetric complex Gaussian (ZMCSCG) with the variance ffr2. In order to ensure that the symbol s can be recovered at the destination, it is assumed that ND, NR, and NT are greater than or equal to NS. The signal received by the destination in first timeslot can be expressed as follows:

yd,i = HoFs + na,i,

where H0 e CNdxM is the MIMO channel between the source and destination and nd,i e CNdX1 is the additive Gaussian noise vector at the destination and follows nd,1 ~ Nc(0, o^Ind). The MIMO relay processes the signal yr using the linear operator Z e CNtXNr and forwards the following signal to the destination in the second timeslot:

yo = ZH1Fs + Znr,

where the relay transceiver Z operates under the power constraint PR = E{yHyo} < PR33* with total relay power of PRmax. The signal received by the destination in the second timeslot is

yd,2 = H2ZH1 Fs + H2Znr + nd,2,

where H2 e cNdXNt is the MIMO channel between the relay and destination and nd,2 is the additive noise described

as nd,2 ~ Nc(0, o^Ind). The double-sided spatially correlated source-relay, relay-destination and source-destination are modelled according to Kronecker model as follows:

h1 = e1/2HwY1/2

1/2ttww1/2

-1 H1 T1

H2 = i)/2w[^2

Ho = s0/2hw<2,

where Z1 e CNRxNR, I2 e CNdxNd, and !0 e CNdxNd are the receive-side spatial correlation matrices, and e CMxM, Y2 e CNtxNt, and Y0 e CMxM are the transmit-side correlation matrices for the channels H1, H2, and H0, respectively. The elements of Hw, Hw, and Hw are ZMCSCG random variables with the unit variance. Note that the transmit and receive spatial correlation matrices are positive semidefinite matrices and are a function of the antenna spacing, average direction of arrival/departure of the wavefronts at/from the transmitter/receiver, and the corresponding angular spread (see [21] and the references therein). The spatial correlation matrices represent the second-order statistics of the channels which vary slowly and can be precisely estimated. However, estimation of the fast fading parts Hw, Hw, and Hw of the spatially correlated MIMO channels can lead to a significant amount of estimation error. For the linear minimum mean-square error (MMSE) estimation, we can write the following error model [22]:

Hw = hw + Ew,

Hw _ TTW . rw

2 = H2 + E2 ,

Hw _ Tjw I cw

0 = H0 + E0>

where Hw, Hw, and Hw are the estimated CSI, and Ew, Ew, and Ew are the corresponding channel estimation errors whose elements are ZMCCSG random variables with the variances ffe2,1, ffe2,2, and ffe2, 0, respectively. Substituting (6) into (5), the error modelling for the actual channels H1, H2, and H0 can be simply given by

H1 = Zl/2Hw^1/2 + £1/2 Ew^1/2 = H1 + E1, H2 = £2/2Hw^1/2 +1!!2 Ew^1/2 = ÍI2 + E2,

Ho = s¿/2Hw^1,/2 + 4/2 ew^0/2 = Ho + Eo,

which shows that the errors Ei, E2, and E0 are also double-sided correlated like the MIMO channels. The destination recovers the source signal s by linearly combining the signals yd,1 (2) and yd,2 (4) of two time slots as follows:

S = W1yd,1 + W2yd,2,

where W1, W2 e CNs xNd denote the linear operators for the signals received from the direct and relay links, respectively. The MSE between s and s can be defined as follows:

M(F,Z, W1, W2) = e{(s - s)(s - s)H|

where the MSE matrix depends on F, Z, W1 and W2 and the mathematical expectation is only taken with respect to noise and signal realizations. Considering that nr, nd,1, nd,2 and s are statistically independent and applying (2) and (4) into (8), we can write (9) as follows:

M(F, Z, W1, W2) = W1

HqffH hH + a2ni iNd

(W1H0 + W2H2ZH1)F - (W1H0F)

- (W2H2ZH1F)H + (w1H0FFH HH Zh HHwH

W1H0 FFH HH Zh HH WH

+ W2 H2ZH1F(H2ZH1F)H

+< H2Z(H2Z)H + < Ind

WH + Ins .

For the given channel estimates H1, H2, and HQ, the MSE matrix of (10) is random due to the random errors E1, E2, and Eq. Since the exact errors are not known and only the covariance matrices of these errors are known, we need to derive the average MSE matrix. This can be done by averaging the MSE matrix (10) over the independent errors E1, E2, and Eq. Hence, we can write

Eeq {I} = H QFFH HH + EEq { EQFFH EHj

= HH QFFH HH + tr( FFH Y0) Z Q,

where ZQ = ff2,QZQ, and we have applied (7) and used the facts that EQ is zero-mean and EX{XAXH} = tr(A)ox;I for X ~

NC(0, oxJI). Similarly, since EQ, E1, and E2 are independent, we can easily show that

ee 0 ,e1 ,e2

EEQ,E1,EJ (Hq + Eq) FFH (H 1 + E1)HZH (H + E2

= H qFFh HHH ZH HH HH.

Furthermore, we can write

Ee,E2 {III} =EE^ HZE^) H1FFHHHj ZHHHj, (13) where the inner expectation is

Ee^ H1FFH HH } = EE J (H1 + E1) H FFH (H1 + E:

H 1FFH H H + tr FFH Y1 )Z 1

and Z1 = ff,2,1Z1. Substituting the result of (14) into (13), we have

Ee,E2 {III} = EE2{ H2ZAZh HHj

= H 2ZAZh HH + tr (ZAZH Y2) Z2,

where Z2 = ^¿Z^ Applying similar steps, we can also get

Ee2 {IV} = H 2ZZH ff H + tr( ZZH Y2) Z2. (16)

Using the results of (11) to (16), the average MSE matrix can be written as follows:

M(F, Z, W1, W2) = Ee!,E2,E0 {M(F, Z, W1, W2)}

1QFF n q

Q Z Q + °rn iNd

+ W1HQFFH HH ZH ff H wH

+( W1 H qFFh ftH Zh HHwH

iv t7- tr / (v \ r^j ~

H2ZAZHHH +tr ZAZHY2 Z2 + < Ind

/ \ W1ii QF + W2H 2ZH 1F

V Dm )

W1H0F + W2H2ZH1F + Ins,

where A = A + a^ INr. The instantaneous relay power can be obtained as follows:

Pr = tr(E{yoyHj) = tr(zn FFHHHZH) + al tr(ZZH),

where expectation is taken w.r.t. noise and signal realizations. After including the estimation error E1 in (18), the relay power averaged over E1 can be expressed with the help of (14) as follows:

3. Joint Optimization: Iterative Approach

The objective of joint optimization is to minimize the sum of the average MSE (17) under power constraints of the source and relay. This optimization problem can be expressed as follows:

min /mse = tr(M(F, Z, W1, W2)) s.t. tr( FF

F,Z,W1 ,W2

< PS03*,

tr ZAZH ) < P

The constraints of the optimization problem (20) do not depend on W1 and W2. As a result, the optimal W1 and W2 can be easily obtained in terms of F and Z. Unfortunately, after substituting such optimal W1 and W2 into the objective function of (20), the resulting objective function in terms of F and Z appears to be a nontractable nonconvex problem. This fact will be later shown in this section. The joint optimization problem (20) is a nonconvex problem and does not offer closed-form solutions. However, it can be easily observed that the considered problem is a convex problem over one optimization variable when others are fixed. Hence, we propose to solve this optimization problem using iterative technique, where each optimization variable is updated at a time considering others as fixed. The iterative algorithm may be implemented as follows. The destination estimates the source-destination and relay-destination channels and the relay estimates the source-relay channel, separately with the help of training sequence. The relay sends the estimated source-relay channel to the destination where the iterative algorithm is executed. The destination feedbacks optimally designed F and Z to the source and relay, respectively. The channel is considered to remain constant within a block but vary from one block to another, where the block consists of training signal and useful data ( Notice that the adaptation of the source precoder and relay transceiver matrices in fast fading scenario can be impractical if the design is based on instantaneous channels [23]. Therefore, robust designs based on channel covariance information [20] can be appropriate for such a scenario. ) .

Remark 1. It is worthwhile to mention here that the minimization of the sum of the source and relay powers under the MSE constraint can also be solved by using the iterative framework that we are proposing in the sequel. Moreover,

the quality of fairness approach such as minimizing the sum of the source and relay powers while fulfilling the SNR/MSE requirements of each symbol stream can also be handled by the proposed iterative method. For conciseness, the latter two methods are not considered in this paper.

After solving the first-order partial derivative of the objective function of (20) w.r.t to W1, we get

W1 = [FHHH - W2Bm] A-J.

Substituting (21) into the objective of (20), the latter can be expressed in terms of W2, Z, and F as follows:

tr(M(F,Z,W2)) = tr(W2 (Cm - BmAm1Bm)WH

DH WH) + N

w7d„

-tr( W2BmAm1H qF + FH H H A-1

X BmWH - HqF

Now, solving the derivative of (22) w.r.t to W2, we get the optimal W2 as follows:

W2 = (DH - FHHHA^Bm) (C„

BmAm1Bj \ (23)

The optimal W1 can be obtained by substituting W2 from (23) into (21). Using the results of (21) and (23) and then resubstituting Am, Bm, and Cm, (22) can be written in terms ofF and Z as follows:

tr(M(F,Z)) = tr(G) - trl H2ZH 1 FGGiniZn F

H2ZH1FG H2ZH1F + r

G = INs - FH HH

HqFFhffH +tr FFH YQ ZQ+and N

FH HH ZN1fiQF + INS

r = H 2Z (Z1 tr (FFH Y1) + a2r InJ ZH ff H

■td ZA ZH Y^ Z 2 + a2d INd .

It is interesting to observe that G is the MMSE matrix of the direct link, where the sum MMSE is simply given by /mmse,DL = tr(G). The second equality for G in (25) is obtained by using the fact that XH(XXH + I)-1X =

PR = tr( ZAZH I.

I - (XHX +1) 1. Applying the same fact and after some manipulations, we get

tr( M(F, Z))

INs + Gh/2Fh H h Zh H h r-1H 2ZH1 FG1/2

tr( [G-1 + Y]-1),

where second equality is obtained after simple steps using the fact that G is a positive definite square matrix. Notice that Y is only related to the MMSE of the double-hop channel, where the sum MMSE is 7mmse,DH = tr (INs + Y)-1. With these observations, we can formulate the following lemma:

Lemma 1. The sum MMSE of the MIMO relay system with the direct link is upper bounded by the sum MMSE of the direct link and source-relay-destination link.

Proof. This Lemma can be easily proven by using the properties of the positive (semi) definite matrices. Since G-1 is positive definite and Y and rt are positive semidefinite, we can show that

G-1 + Y h G-1 — (G-1 + Y) 1 ^ G — tr( (G-1 + Y)

< tr(G) 4 f

mmse,DL>

G-1+ Y = rt + Ins+ Y h Ins+ Y — (rt+ (Ins+ Y))-1 « (Ins+ Y)-1 tr((G-1+ Y)-1) < tr((lNs+ Y)-1) 4 f

The results of (27) and (28) prove the Lemma.

mmse,DH' (28)

It can be seen that the minimization of (26) under source and relay power constraints is a nontractable problem. We have noticed that even in the case of the nonrobust design, such an objective is difficult to handle. This difficulty has motivated us to use the iterative optimization based on the MSE function (17) for which W1 and W2 have been already determined in terms of Z and F (see (21) and (23)). In the following, we show the optimizations over Z and F when other variables are fixed.

(1) Optimization over Z. With some straightforward manipulations of (17) and using the fact that tr(XXH) = ||X|| 2,

the average sum MSE (objective function of (20)) can be alternatively expressed as follows:

fmse = II W1H0+ W2H2ZH1 F - I

-tr(FFHY1)tr(BZZ1Z^ + o2n tr(ZHBZ

■ tr (ZAZHY2) tr (W2Z2WH

+ on2d t^w1 wh + W2WH)

+ t^W1Z0WH) tr(FFHY0),

where B = HHWHW2H2. Applying the following results [24]:

vec(XWY) = (yt (g) X vec(W), tr (XHYXW^ = vec (X)H (Wr (g) y) vec(X) and denoting zL 4 vec(Z) e cNtNrx1, we can write

fmse =

vec( W1H 0Fj + ^ (H1F + ZH DzL + S3 + s4,

)W2H2 )zl-vec(Ins) (31)

D = sA Z^ B) + o2r (InrQS> B + s2{A

s1 4 tr ( FFHY1), s0 4 tr ( FFHY0),

S2 4 tr ( W2Z2WH), S3 4 o2 tr (W2WH ),

54 ^ < tr(W1WHJ +tr(W1ZoWHjso.

The optimization problem w.r.t. Z can be thus written as follows:

p1 : min f mse s t-

[AT 0 Int]

- nmax < PR •

Noting that zHDzl = ||D1/2zl||2, the optimization problem (33) can be written as follows:

P1 : min if + tf s.t.

Zl,(1,(2

vec W1H0F + HiF

XZl - vec(Ins) [At 0 In^1/2

< tu D1/2zL < t2,

< , pmax

Using the notation t 4 [ 11, t2 ]T, the fact that t2 +t2 = tTt and introducing an auxiliary variable t > 0, the problem (34) can be formulated as the following standard convex optimization problem:

Pi : min t s. t.

ZL,i,t

vec WiHoF + HiF 0W2H2 Zl

- vec( INJ || < art,

||D1/2zL|| < brt,

[~T "i 1 A1 0 INtJ

^ /prt,

where aT = [1,0], bT = [0,1], and the quadratic inequality constraint tTt < t is converted to a linear matrix inequality constraint (LMI) using Schur-Complement theorem [16].

Remark 2. Notice that when other variables are fixed, Z can be optimized by solving the Karush-Kuhn-Tucker (KKT) conditions, where the Lagrangian multiplier that arises due to the relay power constraint can be obtained by using the bisection algorithm like in [15]. However, in order to make the proposed iterative approach applicable for other related problems briefly discussed in the beginning of this section and also for the optimization over F in the sequel, we propose to formulate our optimization problem in the convex form P1, which has been proven to be computationally efficient and flexible to accommodate even a large number of convex constraints [16].

(2) Optimization over F. First, we define the following scalars that do not depend on F:

55 4 tr BZ£iZH , 55 4 tr BZZH ,

57 4 tr ZZH ¥2 ,

4 tr Z£iZH¥2 .

After some simple steps and again using the fact that tr(XXH) = ||X||2 in (17), the average sum MSE (20) can also be expressed as follows:

f mse = || ( WiH 0 + W2H2ZH i) F - I№ 11 + S2tr( FFH ¥0) + a2ni 56+ 53 X 52 tr (FHEF) + (55 + 5258) tr (FFH¥i)

+ 5257^n2r+ <tr(WiWH) ,

where s2 = tr(W1ZQWH) andE = fHHZH Y2ZH 1. Noting that (I ® X)vec(Y) [24] and applying (30), we can

vec(XY) write

f mse = || I® (WiHo + W2H2 ZH i fL - vec( Ins )

+ fH [5^Ins® E) + S2 (Ins ® ¥<,) + (55 + 5258)

X (Ins 0 ¥i) ] fL + 53 + < tr (WiWH) ,

where fL = vec(F) e cmn^a. The relay power in terms of fL can be expressed as follows:

Pr = 59 + fH[5i^ Ins® ¥

i J + (INs

Q]fL, (39)

where we have used s9 4 a2 tr(ZZH),s1Q 4 tr(ZZ 1ZH), and Q 4 H— ZHZH 1. The optimization problem w.r.t. F can be now formulated as follows:

P2 : minf mse

IHlN2 < Psmax,

INs ¥i+ Q)] fi

- nmax < PR .

Finally, like in the case of the optimization over Z, we can transform (40) into the following convex problem:

P2 : min t s. t.

fl,t,t

WiII0 + W2H2ZHi) fL - vec(Ins)

< aTt, 11Si/2 fi 11 < bTt,

11 fL N ^Psmax,

iNs (g)(5io¥i+ Q)] fi

' nmax ~ PR - 59,

whereS = [s2(Ins ® E)+S2(In^ Yq)+(s5+S2Ss)(In^ Y1)] and the LMI is equivalent to the quadratic inequality constraint like in the case of P1. Note that the optimization problems P1 and P2 can also be solved first by reformulating them as quadratic matrix programming [25] problems and then applying the semidefinite programming (SDP) relaxation technique. However, since our problems are already convex and second-order cone programming (SOCP) formulation is possible, applying SDP relaxation only increases the computational complexity as we know that the SDP has higher computational burden than the SOCP [26]. The iterative optimization technique can be now summarized as follows:

+ 5257% + ffn, 56

Algorithm 1. The iterative algorithm for joint optimization of W1, W2, Z and F.

(1) Initialize the algorithm with Z0 and F0 such that the source and relay power constraints are satisfied.

(a) repeat

(b) Update W2 using (23).

(c) Update W1 using (21) and (23).

(d) Update Z; by solving the convex problem P1.

(e) Update F; by solving the convex problem P2.

(f) i = i + 1;

(2) until both | f mse,; - 7mse,i+11 is smaller than a threshold e, where the index i denotes the ith iteration.

Remark 3. In the case of the relay channel without the direct link, the optimization problem over destination equalizer, Z and F can be iteratively solved by omitting all terms containing H0 and W1 in (23), P1, and P2.

Remark 4. If the channel estimates are perfect, (21), (23), P1, and P2 can be changed to shorter forms by using the fact that o^ = o22 = o20 = 0. The resulting equations and optimization problems correspond to the perfect CSI case.

3.1. Computational Complexity. The computational complexity of Algorithm 1 mainly depends on the work loads of the convex optimization problems P1 and P2 which consist of second-order cone (SOC) as well as SDP constraints. An enormous amount of effort will be required to compute the exact computational complexity of P1 and P2, and thus their exact complexity analysis is beyond the scope of this paper. However, using the results of [26], we determine the worst-case computational complexity of P1 and P2. In P1, there are 3 SOC constraints where the first SOC constraint has a dimension (also the size of the cone) with 2N| + 1 real variables and the remaining two SOC constraints have the same size of 2NTNR + 1. The SOC constraints in P1 consist of 2NtNr + 2 real optimization variables. The only one SDP constraint of P1 has a size of 3 with 3 real optimization variables. Therefore, according to [26], the computational load of P1 per iteration is O((2NTNR + 2)2(2NS2 + 4NTNR + 3) + 81). The number of iterations required can be upper bounded by O(2 V3). Hence, the overall worst-case complexity of P1 is 0(2V3((2NtNr + 2)2(2Ns2 + 4NtNr + 3) + 81)). Similarly, we can compute the worst-case computational load for P2. The first SOC constraint in P2 has a size of 2N| + 1 while the other SOC constraints are of the same size of 2MNS + 1. The SOC constraints consist of 2MNS + 2 real optimization variables. Like in the case of P1, the single SDP constraint is of size 3 with 3 real optimization variables. Thus (see also [26]), we find that P2 has a work load of O((2MNs + 2)2(2NS2 + 6MNs + 4) + 81) per iteration, where the required number of iterations is upper bounded by O(V3 + 2). As a result, the worst-case complexity of P2 is O((V3 + 2)((2MNS + 2)2(2NS + 6MNS + 4) + 81)). Notice that in practice the interior-point algorithms used for

solving P1 and P2 behave much better than predicted by the aforementioned worst-case analysis [27].

3.2. Convergence Analysis. It can be shown that the proposed iterative method converges. It has been already discussed that the optimization problem is convex w.r.t. each variable when the others are fixed. For the given F and Z, the solutions given by (21) and (23) correspond to the MMSE receiver. As a result, we have tr(M(F;, Z;,W1+1,W2+1)) < tr(M(F;,Z;, W1, W2)). Similarly, for the given W1, W2, and F, the optimization problem (20) is convex w.r.t. Z and, thus, the problem P1 provides optimal solution for Z which means that tr(M(F;,Z;+1, W1, W2)) < tr(M(F;,Z;, W1, W2)). Finally, for the given W1, W2, and Z, problem (20) is convex w.r.t. F and, hence, the problem P2 gives optimal solution for F thereby confirming tr(M(F;+1,Z;, W1, W2)) < tr(M(F;,Z;, W1, W2)). Therefore, it can be found that with each update of W1, W2, F, and Z, the objective function decreases and the iterative method converges. It will be later shown via numerical results that the iterative algorithm gives satisfactory performance with acceptable convergence speed. However, the global optimality of the solutions of the iterative method for the relay channel with the direct link cannot be guaranteed as the joint problem is nonconvex. In the next section, the joint optimization problem will be restricted to a double-hop MIMO relay, where receive antenna correlations are assumed to be negligible for each hop. In this case, the optimization problem turns to the joint power allocation, for which the global optimality can be guaranteed under high-SNR approximation.

4. Joint Power Allocation with GP: Without Direct Link

Due to severe shadowing, hotspots, and so forth, the destination may be out of the coverage area of the source. In such a scenario, it is reasonable to consider that the direct link between the source and destination does not exist. In the latter case of the relay channel without the direct link, the sum MSE can be expressed in terms of F and Z (see also (28)) as follows:

mmse,DH

tr INs + FHHHZhHHr-1H2ZH 1F

The objective is to minimize the average sum MSE (42) under the constraints of source and relay powers. Thus, the optimization problem is

min tr(MMSE(Z, F)) s.t. tr(FFH) < PS™,

tr( ZAZH) < PR™,

where the MMSE matrix is MMSE(Z, F) = [INs +

FHIIHZHHHr-1il2Zii 1F]-1. Let A(MMSE(Z, F)) be the vector of eigenvalues and d(MMSE(Z, F)) be the vector of the diagonal elements of MMSE(Z, F) in decreasing order. In

this case, d is said to be majorized by A, and the Schur-concave function can be defined as /(A(MMSE(Z,F))) < /(d(MMSE(Z, F))), where /(x) stands for the function of x. It is clear that the minimum of this function is obtained if the diagonal elements of MMSE(Z, F) become its eigenvalues [8]. This can occur when MMSE(Z, F) is a diagonal matrix with its elements in decreasing order. Let the singular value decompositions of H1 and H2 be

H1 = U1A1/2VH, H2 = U2 A2/2V;

1/2 A H

where the diagonal elements of A1 and A2 are considered to be in the decreasing order. If these elements are not in decreasing order, some permutation matrices can be applied to make the diagonal elements of A1 and A2 in the decreasing order [8]. This means that, in (44), the permutation matrices are implicitly included. The closed-form expressions for the optimal Z and F are difficult to obtain in general. However, for the double-hop channels without receive-side antenna correlations, that is, when Z1 = INr and Z2 = INd, the optimal Z and F that diagonalize the MMSE matrix can be given by

V2 aZ/2uH ,

F = V,A,

where AZ and AF are the diagonal matrices with the elements in decreasing order. Substituting (44) and (45) into the MMSE matrix and after some straightforward steps, we get

/ AT/2AT/2AT/2AT/2r-1

J mmse = tr I AF A1 AZ A2 rD

X (AT/2AT/2AT/2AT/2) 2 + INs

Td = (ah trfBAF/2AT72) + aDA}/2AZ/2AT/2AT/2

D = \Oe, 1tqBaf AF ) + °nx) a2 AZ AZ a2

+ tr( C AZ/2A}/2AF/2AFT/2AT/2AT/2) ae2>2lND+ ah'

tr(CAZ/2AT/2) [a21tr(BaF/2aT/^ + al]Ind + a2niInd.

In (47), we use B = V-Y1V1 and C = V-Y2V2. Similarly, the source and relay powers become

Ps = tr(AF/2AT/2 ,

Pr = tr( AZ/2A1/2AF/2AT/2AT/2AZ/2)

-tr(AZ/2AT/^ [a2,1t^(BaF/2aT/2) + al].

Let q = min(M,NS) and v = min(NT,NR), and {yfxj]j=1 and let 1 be the nonzero diagonal elements of A

1 /2 F

and A^2, respectively, in the decreasing order. For brevity, we also define

d 4 a^tr B AF/2AT/2 + a2) 4 I a^ B ^

'lF + anr

e 4 t^CAZ/2A1/2AF/2AT/2AT/2AT/2)ae2,2 4 a22.

/ 4 ae2,2 tr(CAZ/2AZ/2)[a21 tr(BAF/2AT/2) + an,]

4 ae2,2zL ^^¡,i^Z

IB ¡,iAFa21 + an2r ¡=1

where B¡,i and C¡,i are the ¡¡th diagonal elements of B and C, respectively, and p = min(NT,NR,M,NS). Since B and C are positive semidefinite, we have B¡j > 0 and C¡j > 0 for all i. Using (47) and (49), the sum MMSE can be finally expressed as follows:

11+ AFAZA1A2) / dAZA$ + e + / + a2

1 + SNRr

where SNRr = XrFXrZXr1Xr2/(dXrZXr2 + e + / + a^ is the SNR of each data stream provided that NS is also smaller than or equal to M and R = min(NT,NR,M,NS, and ND). It is interesting to note that when the channel estimates are perfect, that is, when a^ = a^,2 = 0, the sum MMSE of (50) reduces to the objective function of [28]. Applying (48), the joint power allocation problem can be formulated as follows:

{xF}q=1,{xz}:=1

/mmse s.t. ¿>F < P;

max F < PS

Y^x^xy + d£xZ < Pmax,

m=1 k=1

which is a nonlinear and nonconvex problem. Since this problem is nonconvex, the global optimal solution is difficult to obtain. Considering the fact that the global optimal solutions of the problems similar to the nonrobust version of (51) can be obtained only with very high computational cost (see [29, 30]), the authors of [31] use an iterative waterfilling technique to solve the nonrobust form of (51). However, we have noticed that it is hard to solve (51) using the waterfilling method of [31]. The major difficulty arises from the fact that the first-order partial derivatives of the corresponding Lagrangian function w.r.t. xF for the fixed {xrZ}r=1 and w.r.t. xrZ for the given {xrF}R=1 do not lead to equations that are decoupled in xF and xZ, respectively. It is easier to see that this difficulty appears due to the reason that d, e, and / in (51) consist of not only xF and xZ but also xF and

A|, for all k e {1,...,R}. Furthermore, although iterative waterfilling method is computationally efficient, it does not guarantee the global optimal solution. In the following, we use an alternative approach based on GP technique. Note that the optimization problem (51) is not a GP problem but a signomial programming (SP) problem [16] which can be iteratively solved as a GP problem after approximating the required posynomial terms by monomial terms. It is known that the SPs do not guarantee global optimality and the computational cost is high. Thus, using the high-SNR approximation, (51) can be solved as a GP whose global optimality can be confirmed. In this regard, we have

f ^12 + e + f + <

J mmse ~ Z-, iririnr •

r = 1 1F1Z1112

Using the upper bound (52), the power optimization problem can be expressed as follows:

AF}L,{4h=i =

]wr s.t.

d\rZXr2 + e + f + o2d < trXrFXrZX[Xr2, V r q

XaF < Psmax,

XAzAZA? + d^>Z < Pjfx,

m = 1 k=1

which is a GP problem that can be solved efficiently to guarantee the global optimality.

Remark 5. Notice that the joint power allocation problem is nonconvex even for the case without channel estimation errors [32]. In such a nonrobust design, it has been shown in [28] that the lower and upper bounds can be established for the MSE for each data stream. Unfortunately, this is not the case for the proposed robust design due to the fact that the terms e, f, and d are again functions of AZ, for all k = 1,..., v and AF, for all j = 1,..., q. It is also worthwhile to note that several optimization problems which can be solved using the GP method and still provide solutions close to the optimal solution of the sum MSE minimization problem are; (a) maximization of the minimum of the SNRs of the data streams and (b) maximization of the geometric mean of the SNRs of the data streams, both under the source and relay power constraints.

5. Numerical Results and Discussions

In this section, the performance of the proposed methods will be investigated. The proposed robust designs are also compared with the nonrobust case, where the channel estimation errors are not taken into account. Notice that the nonrobust design corresponds to the so called naive design, where the optimization problems of interest are

solved assuming that there are no errors in channel estimates (see also Remark 4) although in reality the estimates are erroneous. For all numerical simulations, we take NS = 3, M = 4, Nr = 3, NT = 4, and ND = 3. The spatial covariance matrices for source-relay, relay-destination, and source-destination channels are modelled according to the widely used exponential correlation model. In our examples, we take

Zl = Z2 = Zo =

= Y2 = Yo

1 P P' P 1 P

P2 P 1

1 a a2 a3 a 1 a a2

a3 a2 a 1

In all cases, the optimization problems p1, P2, and (53) are solved using the CVX software [33]. The SNRs for source-relay, source-destination, and relay-destination channels are defined as SNRsr = P^Vo2, SNRsd = P^Vo^, and SNRrd = PRlax/o2d, respectively. Throughout all examples, we take Pmax = Pr3x = 0 dBw and vary the values of o2r and o2d to change SNRsr, SNRsd, and SNRrd. The estimated channels are generated according to (7), so that the elements of actual channels H1, H2, and H0 have the variance of 1. For all results, we compute the average MSE by taking 200 realizations of the estimated channels.

The convergence behaviour of the proposed iterative method as a function of iteration index is illustrated in Figure 2 for the relay channel with the direct link. The parameters in this figure are set as a = 0.2, 3 = 0.1, oe2,1 = o22 = 0.01, and oe2,0 = 0.03. We take three sets of SNRs as SNRsr = SNRsd = SNRrd = 10 dB, SNRsr = SNRsd = SNRrd = 20 dB and SNRsr = 10 dB, and SNRsd = SNRrd = 0 dB. It can be seen from this figure that the iterative method converges in about 15 iterations. The convergence is faster for the lower values of the SNR. The effect of different initializations on the convergence behaviours of the iterative method is also displayed in Figure 2. The convergence speed for the cases where F and Z are initialized to randomly generated matrices of ZMCSCG random variables is similar to the cases where F and Z are initialized to the matrices proportional to identity (i.e., F oc I and Z oc I). Moreover, it can be noticed from Figure 2 that different initializations lead to the similar solutions. In Figure 3, the performance of the iterative method as a function of the iteration index is illustrated for the relay channel without the direct link. We take a = 0.3, 3 = 0, and o^ = o22 = 0.01 in this figure. As a reference, the performance of the GP problem which gives optimal solution under high-SNR assumption is also displayed in Figure 3. It can be noticed from Figure 3 that the difference between the solutions of the iterative method and GP method is negligible after 10-15 iterations.

The performance of the proposed iterative method for the MIMO relay channel with the direct link is shown in Figure 4 for different values of o2. The performance of the

1.4 1.2 1

10 15 20

Iteration index

SNRSr,rd>sd SNRSr,rd>sd SNRSr,rd,sd SNRSr,rd,sd SNRsr SNRsr

; 20dB , F, Z m I : 10 dB , F, Z M I : 10 dB, F, Z are ZMCSCG ; 20dB , F, Z are ZMCSCG : 10 dB,SNRrd,sd = 0dB , F, Z are ZMCSCG : 10 dB,SNRrd,sd = 0dB , F, Z m I

Figure 2: Convergence behaviour of the iterative approach for different initializations with the direct link (a = 0.2, ji = 0.1,

Q 5 1Q 15 2Q

SNRrd = SNRsd (dB)

-e— Prop.robust al = Q.QQ5 -0— Non-robust a2e= Q.QQ5

Robust with F m I, a} = Q.QQ5 -B- Prop.robust a] = Q.QQ8 - <i - Non-robust al = Q.QQ8 ■ - Robust with F m I, a] = Q.QQ8

Figure 4: Sum MSE as a function of SNR (SNRrd relay channel with the direct link (a = 0.6, ji =

SNRsd) for the

o\2 = 0.01, o2e0 = 0.03, andPf" = Pjf" = 0 dB).

= a2e, and Pfax = Pg

0 dB).

0.4, а2еЛ

H 0.5 S

3 0.4 0.35 0.3 0.25 0.2

0 j> 0 0 >00000000<>»VV00i ff Щ g g У W W % 1

0 5 10 15 20

Iteration index

—e— SNRsr,rd = 20 dB , iterative method —в— SNRsr,rd, = 20 dB, GP-power allocation -»- SNRsi' ' SNRsr

: 15 dB, SNRrd= 20 dB, GP-power allocation : 15 dB,SNRrd= 20 dB , ilterative method

Figure 3: Convergence performance of the iterative approach and the GP power allocation method for the case without the direct link (a = 0.3, ji = 0, a\x = ae2>2 = 0.01, andPf* = PJT1 = 0 dB).

nonrobust method which considers the estimated channels as actual channels and the performance of the robust method without a source precoder [15] (the case with F o cI where c is the positive scaling factor chosen for satisfying the source power constraint) are also displayed in this figure. We keep a = 0.4, в = 0.6, and SNRsr = 30 dB and change SNRsd

and SNRrd from 0 to 30 dB. The threshold e for stopping the iterative process is set to 1e - 4. It can be noticed from this figure that in all cases, the MSE decreases when the SNR increases and when the variance of the channel estimation error al decreases. Furthermore, the proposed robust design outperforms both the nonrobust method and the robust method with F o cI [15]. In Figure 5, the performance between the proposed iterative method and the GP power allocation method for the relay channel without the direct link is compared. This figure displays the sum MSE as a function of the correlation coefficient a for different values of SNRsd and SNRrd. We keep в = 0, ae2 = 0.01 and e = 1e - 4 for Figure 5. It can be observed from this figure that the proposed robust methods significantly outperform the nonrobust method. Moreover, the performance gap between the iterative approach and the GP power allocation method is negligible for a < 0.4. For a > 0.4, the iterative approach outperforms the GP method. In all cases, the sum MSE increases with the increasing values ofcorrelation coefficient.

6. Conclusions

The problem of jointly optimizing the source precoder, relay transceiver, and destination equalizer for a cooperative MIMO relay system has been treated in this paper. An iterative approach with the guaranteed convergence has been proposed to solve the nonconvex problem. It has been shown that for the case without the direct link where the two-hop MIMO channels have only transmit-side antenna correlations, the joint optimization turns to the joint source and relay power allocation problem which has been solved by

0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1

—в— Robust-GP, SNRsr = 30 dB, SNRrd = 20 dB

Non-robust-iterative, SNRsr = 30 dB, SNRrd = 20 dB —v— Robust-iterative, SNRsr,rd = 30 dB -tr- Robust-GP, SNRsr,rd = 30 dB —«— Non-robust-iterative, SNRsr,rd = 30 dB

Figure 5: Sum MSE as a function of a for the relay channel without the direct link (ß = 0, a\x = = 0.01, and Pf" = P^ = 0dB).

using GP technique under high-SNR approximation. Simulation results confirm the efficiency and good performance of the iterative approach for the MIMO relay channel with and without the direct link. Furthermore, the proposed robust methods significantly outperform the nonrobust method when the channel estimates are imperfect.

Acknowledgments

This work is supported by the European project NEW-COM++, the Wallone region project COSMOS, and the concerted action SCOOP. B. K. Chalise is currently with the Center for Advanced Communications, Villanova University, Villanova, PA 19085, USA, (Phone: +1-610-519-7371, Fax: +1-610-519-6118, Email: batu.chalise@villanova.edu). L. Vandendorpe is with the Communication and Remote Sensing Laboratory, Universite catholique de Louvain, Place du Levant, 2, B-1348 Louvain la Neuve, Belgium, (Phone: +32-10-472312, Fax: +32-10-472089, Email: luc.vandendorpe@uclouvain.be). The first author completed his part of contribution when he was with the Communication and Remote Sensing Laboratory, Universite Catholique de Louvain.

References

[1] J. N. Laneman and G. W. Wornell, "Exploiting distributed spatial diversity in wireless networks," in Proceedings of the Allerton Conference on Communication, Control, and Computing, Monticello, Ill, USA, October 2000.

[2] M. Janani, A. Hedayat, T. E. Hunter, and A. Nosratinia, "Coded cooperation in wireless communications: space-time

transmission and iterative decoding," IEEE Transactions on Signal Processing, vol. 52, no. 2, pp. 362-371, 2004.

[3] G. Kramer, M. Gastpar, and P. Gupta, "Cooperative strategies and capacity theorems for relay networks," IEEE Transactions on Information Theory, vol. 51, no. 9, pp. 3037-3063, 2005.

[4] R. U. Nabar, H. Bolcskei, and F. W. Kneubuhler, "Fading relay channels: performance limits and space-time signal design," IEEE Journal on Selected Areas in Communications, vol. 22, no. 6, pp. 1099-1109,2004.

[5] Y. Jing and H. Jafarkhani, "Using orthogonal and quasiorthogonal designs in wireless relay networks," IEEE Transactions on Information Theory, vol. 53, no. 11, pp. 4106-4118,

[6] V. Havary-Nassab, S. Shahbazpanahi, A. Grami, and Z.-Q. Luo, "Distributed beamforming for relay networks based on second-order statistics of the channel state information," IEEE Transactions on Signal Processing, vol. 56, no. 9, pp. 4306-4316,

[7] B. Wang, J. Zhang, and A. H0st-Madsen, "On the capacity of MIMO relay channels," IEEE Transactions on Information Theory, vol. 51, no. 1, pp. 29-43, 2005.

[8] O. Munoz-Medina, J. Vidal, and A. Agustín, "Linear transceiver design in nonregenerative relays with channel state information," IEEE Transactions on Signal Processing, vol. 55, pp. 2593-2604, 2007.

[9] X. Tang and Y. Hua, "Optimal design of non-regenerative MIMO wireless relays," IEEE Transactions on Wireless Communications, vol. 6, no. 4, pp. 1398-1406, 2007.

[10] C.-B. Chae, T. Tang, R. W. Heath Jr., and S. Cho, "MIMO relaying with linear processing for multiuser transmission in fixed relay networks," IEEE Transactions on Signal Processing, vol. 56, no. 2, pp. 727-738, 2008.

[11] R. Zhang, C. C. Chai, and Y.-C. Liang, "Joint beamforming and power control for multiantenna relay broadcast channel with QoS constraints," IEEE Transactions on Signal Processing, vol. 57, no. 2, pp. 726-737, 2009.

[12] B. K. Chalise and L. Vandendorpe, "MIMO relay design for multipoint-to-multipoint communications with imperfect channel state information," IEEE Transactions on Signal Processing, vol. 57, no. 7, pp. 2785-2796, 2009.

[13] A. S. Behbahani, R. Merched, and A. M. Eltawil, "Optimizations of a MIMO relay network," IEEE Transactions on Signal Processing, vol. 56, no. 10, pp. 5062-5073, 2008.

[14] W. Guan, H. Luo, and W. Chen, "Linear relaying scheme for MIMO relay system with QoS requirements," IEEE Signal Processing Letters, vol. 15, pp. 697-700, 2008.

[15] X. Chengweng, M. Shaodan, and W. Yik-Chung, "Robust joint design of linear precoder and destination equalizer for dual-hop amplifyand-forward MIMO relay systems," IEEE Transactions on Signal Processing, vol. 58, no. 4, pp. 2273-2283, 2010.

[16] S. Boyd and L. Vandenberghe, Convex Optimization, Cambridge University Press, Cambridge, UK, 2004.

[17] S. A. Vorobyov, A. B. Gershman, and Z.-Q. Luo, "Robust adaptive beamforming using worst-case performance optimization: a solution to the signal mismatch problem," IEEE Transactions on Signal Processing, vol. 51, no. 2, pp. 313-324, 2003.

[18] A. Pascual-Iserte, D. Perez Palomar, A. I. Perez-Neira, and M. Aí . Lagunas, "A robust maximin approach for MIMO communications with imperfect channel state information based on convex optimization," IEEE Transactions on Signal Processing, vol. 54, no. 1, pp. 346-360, 2006.

Robust-iterative,SNRsr = 30 dB, SNRrd = 20 dB

[19] X. Zhang, D. P. Palomar, and B. Ottersten, "Statistically robust design of linear MIMO transceivers," IEEE Transactions on Signal Processing, vol. 56, pp. 3678-3689, 2008.

[20] B. K. Chalise, S. Shahbazpanahi, A. Czylwik, and A. B. Gershman, "Robust downlink beamforming based on outage probability specifications," IEEE Transactions on Wireless Communications, vol. 6, no. 10, pp. 3498-3503, 2007.

[21] J. Luo, J. R. Zeidler, and S. McLaughlin, "Performance analysis of compact antenna arrays with MRC in correlated Nakagami fading channels," IEEE Transactions on Vehicular Technology, vol. 50, no. 1, pp. 267-277, 2001.

[22] T. Yoo and A. Goldsmith, "Capacity of fading MIMO channels with channel estimation error," in Proceedings of the IEEE International Conference on Communications (ICC '04), pp. 808-813, Paris, France, June 2004.

[23] S. Zhou and G. B. Giannakis, "Optimal transmitter eigen-beamforming and space-time block coding based on channel correlations," IEEE Transactions on Information Theory, vol. 49, no. 7, pp. 1673-1690, 2002.

[24] R. A. Horn and C. R. Johnson, Topics in Matrix Analysis, Cambridge University Press, Cambridge, UK, 1991.

[25] A. Beck, "Quadratic matrix programming," SIAM Journal on Optimization, vol. 17, no. 4, pp. 1224-1238, 2007.

[26] M. S. Lobo, L. Vandenberghe, S. Boyd, and H. Lebret, "Applications of second-order cone programming," Linear Algebra and Its Applications, vol. 284, no. 1-3, pp. 193-228, 1998.

[27] Y. C. Eldar, A. Ben-Tal, and A. Nemirovski, "Robust mean-squared error estimation in the presence of model uncertainties," IEEE Transactions on Signal Processing, vol. 53, no. 1, pp. 168-181,2005.

[28] C. Li, X. Wang, L. Yang, and W.-P. Zhu, "A joint source and relay power allocation scheme for a class of MIMO relay systems," IEEE Transactions on Signal Processing, vol. 57, no. 12, pp. 4852-4860, 2009.

[29] W. Zhang, U. Mitra, and M. Chiang, "Optimization of amplify-and-forward multicarrier two-hop transmission," submitted to IEEE Transactions on Communications.

[30] R. Cendrillon, W. Yu, M. Moonen, J. Verlinden, and T. Bostoen, "Optimal multiuser spectrum balancing for digital subscriber lines," IEEE Transactions on Communications, vol. 54, no. 5, pp. 922-933, 2006.

[31] Y. Rong, X. Tang, and Y. Hua, "A unified framework for optimizing linear nonregenerative multicarrier MIMO relay communication systems," IEEE Transactions on Signal Processing, vol. 57, no. 12, pp. 4837-4851, 2009.

[32] I. Hammerstrom and A. Wittneben, "Power allocation for amplify and forward MIMO-OFDM relay links," IEEE Transactions on Wireless Communications, vol. 6, no. 5, pp. 17741786, 2007.

[33] M. Grant and S. Boyd, "CVX: Matlab software for disciplined convex programming," June 2009, http://stanford.edu/~boyd/ cvx.