1. Introduction
Cooperative communications has been proposed as a positive anti-interference scheme [1]-[5]. As a key technique in 3GPP LTE-Advanced, Coordinated Multi-Point (CoMP) is the most effective method to enhance the performance of the cell-edge user (CEU). Therefore, the CoMP technique research is of important significance.
1.1. Previous Research
Ever since proposed, CoMP has been paid much attention to. It covers downlink precoding design, cell selection, limited feedback codebook design, and robust design and so on. The following gives simple summary of the exiting researches.
Precoding design of CoMP has been studied in [6]-[7], based on the assumption that local channel knowledge is known by the base station (BS). They all choose virtual signal-to-interference-plus-noise ratio (SINR) as the objective function, and achieve the precoding vector by beamforming (BF) and a generalized zero forcing (ZF) algorithm respectively.
The above precoding is under the perfect channel information (CSI). However, the provision of perfect CSI is often a formidable task in wireless systems. A non-ideal channel can typically be obtained due to different errors in practice. In order to ensure the network performance under non-ideal conditions, robust optimization theory was proposed in [8] and [9]. The SINR and minimum mean square error (MMSE) of a single cell downlink transmission was studied in [8]. Paper [9] applied the general estimation model assumption that error lies in hyper-sphere body and designs a linear precoding matrix through maximizing weighted sum-rate and the minimum rate under the worst conditions individually.
There are usually two strategies in cooperation set selection. Previous researches mainly focused on cooperation set selection based on receiving state (e.g., outage probability, the average signal-to-noise ratio or the receiving power), whose disadvantage includes: (1) it may lead to frequent switch of coordination set; (2) it may cause interference, e.g., selecting further BS to coordinate will influence neighboring cell; (3) the reliability may be affected, e.g., the delay of further BS is quite large; (4) outage may occur due to load unbalancing, i.e., the busy BS participates in the coordination set. In addition, there are three ways for CoMP collaboration set selection [11]-[15].
Static collaboration selects several fixed base stations to cooperate according to certain criteria. Although this approach is simple, if the cooperated set for users at different locations is the same, it may not be able to eliminate the interference effectively. Dynamic collaboration is to select the base station dynamically based on the feedback. Although this method can maximize the elimination of the inter-cell interference, the complexity increases rapidly at the same time. Therefore, there is a compromise between the semi-dynamic cooperation combining static collaboration and dynamic collaboration. In such a way, one large static cooperating set is predetermined, and the users select base stations involved in the set according to the criteria dynamically.
How to set predetermined set is a problem of semi-dynamic collaboration. The system tends to be more complicated if the set is too big. If the set is too small, the system is tantamount to static collaboration.
For the above reasons, this paper provided and investigated region-partitioning problem of coordinated multi-point transmission based on the semi-dynamic cooperation idea and MBSFN regional planning. The problems are in different regions of the cell whether the user should open the CoMP working mode or not and how to select cooperation set. Of course, the core work of region-partitioning problem is involved in the user’s location information, channel information, SINR and so on.
1.2. Outlines of the paper
Based on the ideas mentioned above, robust design of coordinated set planning considering the feedback delay and link error will be provided and studied in this paper. The main contributions of this paper can be briefly summarized as follows:
The rest of the paper is organized as follows. In Section II, the various downlink models of a coordinated multi-cell multi-user system with different non-ideal cases are introduced, and the robust net ergodic capacity optimization problem is proposed. Upper and lower bounds of the robust capacity are investigated and a practical scheme is proposed to determine the optimal number of cooperative BSs in Section III. In section IV, simulation results are shown and analyzed. The paper is concluded in Section V.
1.3. Notations
Bold fonts in both lower and upper cases are used to denote vectors and matrices, respectively. If not explicitly stated, the dimensions will be clear from the context. I is the identity matrix and 0 is the zero-matrix. The trace of a matrix is denoted by denote the spectral norm and the Frobenius norm respectively. The conjugate (Hermitian) transpose is written as . diag(a1, L, ak ) is a diagonal matrix with elements a1, L, ak on the main diagonal. is the expectation operator.
2. System Model
This section first gives the various downlink models of a coordinated multi-cell multi-user system with different non-ideal cases and then simplifies the models.
2.1. Model with delayed feedback link
A coordinated multi-cell multi-user downlink system is considered, as shown in Fig. 1. This system includes M cooperative cells, each cell including K users and one base station. To simplify the analysis, we still assume that each user only has one antenna and each base station is equipped with Nt antennas.
Fig. 1.The model of multi-cell cooperation considering feedback delay and error
Assume that the channel remains unchanged within a single time slot and it is subject to flat Rayleigh fading, so the received signal of the user k in the n th channel slots is
where, xm [n ] represents the transmission data from base station m in the n th time slot; wm [n ] ∈ CNt ×1 means corresponding precoding vector; zk[n ] indicates that the additive white Gaussian noise at user equipment (UE) k, with zk[n ] ~ CN (0, σ2I); shows the base station m to the k th user's channel; ak,m is the path loss and hk,m [n ] stands for the small-scale fading. In order to facilitate the subsequent analysis, these channels are mutually independent upon the assumption that the base station antenna distance is large enough.
Group the UE received signals into a vector and denote transmission signals after precoding to arrive at the send antenna port by then the received signal of the user k can be rewritten as
in which, Assuming that the channel Hk[n ] is well modeled as a spatially white Gaussian channel, with entries hk,m[n ] ~ CN (0, I), and the channels are i.i.d. over different users. The average transmit power of BS is P, so the power constraint simplification is
In order to facilitate the numerical analysis, we employ Gaussian Markov stationary ergodic block fading channel. Denote the length of one slot by Ts, and it remains unchanged in a time slot interval. Besides, feedback delay is referred to τ. τ = DTs indicates that the channel information obtained by the transmitting side is delayed D symbol period. The channel vector can be rewritten as [18]
where, the correlation coefficient ρ is ρ = J0(2πfdτ) based on the classic Clerk equi-direction scattering model. It is found that ρ is decided by the product of Doppler shift and feedback delay, which is referred to as Doppler delay product. J0(⋅) is the zero-order Bessel function. The variance of channel error vector, ek,i[n ], is 1 - ρ2 with the distribution ek,i[n ] ~ CN (0, (1 - ρ2) I ), which is mutually independent of hk,i[n - D]. Noteworthy, τ = 0 corresponds to no feedback delay, i.e. D = 0, ρ = 1, and at the moment CSI is perfectly known.
2.2. Channel estimation model
The length of the channel resource block can be divided into the length of the training pilot and the length of the transmission data. In the training period, BSs send the orthogonal pilot symbol, and users use MMSE or other error estimation methods to estimate channel coefficient. Actual channel can be decomposed into the estimated channel vector and error vector resulted from feedback delay e1k,m[n ], so there is
Suppose that at the beginning of each slot, the user can apply the pilot signal estimation to acquire channel information . After obtaining the estimated channel state information, each user quantifies channel quality information (CQI) first, and then transmits the limited bits through a feedback channel to base station controller.
2.3. Channel quantitative model
For the limited feedback frequency division duplex (FDD) system, we generally choose the closest quantized codeword from the unit vector codebook set by the inner product [16] measure
where h is feedback channel, cl is codeword.
Similarly, due to the fact that optimal quantization vector is generally unknown, [17] proved random vector quantization (RVQ) theoretical analysis can provide a performance close to optimal quantization vector. Therefore, in this section only the simplest isotropic distribution of the random vector quantization model is considered. can be decomposed into CQI and channel direction information (CDI). The important criterion to measure the channel quantization error is mean square angle distortion. And the quantization error is defined as is the quantized version of .
2.4. Backhaul link model
The channel vector needs to be transmitted through the backhaul link from BS to base station control (BSC). In addition to link delay, there are link errors, such as bit 0 being misjudged as 1. This model only takes errors in the link between BS to BSC into consideration. After signals from the source node of the relay system passing through a series of relay nodes [19], then arriving at the destination node, the receiving signal at the BSC can be expressed as
where hBS [n ] is the receiving information of the UE channel at the BS, zCS[n ] ~ CN (0, σ2CSI). The andom error matrix Ek is related to the random failure of backhaul links, as in [19]
where ek,m is identically distributed Bernoulli random variables, and its distribution function is, P(ek,m = 0) = ε.
3. Problem Formulation of the Coordinated Set Planning
On the principle of cooperative cell clustering in [7] and for single cluster collaboration model, a typical double-cell cellular collaboration system composed by seven hexagonal cell is considered in this paper. The multi-user downlink model is shown in Fig.2. Referring to zoning standard of the traditional relay systems and distributed antenna systems in [21], a radius r is first provided to draw the boundaries of CoMP and non-CoMP area. Then, in CoMP area we select the cooperative BSs by the principle of proximity, and the boundary is divided by connection line between base stations.
Fig. 2.The model of CoMP based on the multi-user scenario
Due to cellular system with good symmetry and assuming all users and base stations are uniformly distributed, the approximately equal probability density function for (ρ, θ) is
Considering the above error models, coordinated set planning issues take the robust ergodic capacity for the worst case as optimization objectives, and the mathematical model can be expressed as
where shows the subscriber capacity located at the region of non-CoMP. is the robust capacity upon max-minimum principle, that is, according to the error model of non-ideal channel, optimize the precoding vector of the worst case to maximum coordinated capacity. After is gotten, coordinated set planning model is discussed as follow.
Taking hardware complexity and bit overhead of CoMP into account, a penalty factor ξ(ξ≤ 1) is proposed to express the equivalent capacity loss. The robust net ergodic capacity can be expressed as
According to Ref. [10], , where T is the coherence time interval, T0 is additional bit overhead and transmission delay interval, Mc is the number of the cooperated base station.
Due to cellular system with good symmetry and assuming all users and base stations are uniformly distributed, the ergodic capacity for the user status in all position just requires to study a single one-twelfth of the triangle area covered. Removing the constant factor from (11), finally the optimization goal can be simplified to
4. Analysis of the Coordinated Set Planning based on the Robust Capacity
In this section, we investigate the robust capacity and analyze upper and lower bounds of the robust capacity. Based on upper bound of the robust capacity, we define the cooperative gain, and then propose a practical scheme to determine the optimal number of cooperative BSs.
4.1. Analysis of robust capacity
Before transmitting the data in n time slot, we need to describe the channel information in detail. According to channel estimate model, quantized model and feedback delay model, gotten by BS is
Based on the formulation of and random error model of backhaul link, whose noise is negligible, the channel information obtained by controller is
Apply the above results(13) and (14), the relationship between the final channel information for precoding and actual channel information is
where e1k,m[n ] is error vector resulted from feedback delay, e2k,m[n ] contains quantized error and estimate error. Denote the small scale fading channel matrix as and the large scale fading matrix as . Combine (15) into matrix form
where E2k[n ] = [e2k,1[n ], L ,e2k,M [n ]]H , E1k[n ] = [e1k,1[n ],L ,e1k,M [n ]]H.
To facilitate the numerical analysis, an analog feedback with prediction is as shown in [18]. Typically, for analog feedback with d step MMSE predictor and the Gauss-Markov model, the error variance is , where ρ is correlation coefficient and ε0 is the Kalman filtering mean-square error.
According to (16), the received signal at user k can be rewritten as
Furthermore, (17) can be written as
Let Zk [n ] = (ρE1k[n ] + E2k[n ])x[n ] + zk[n ]. The received signal to interference plus noise ratio (SINR) at user k is thus equal to:
When the distributed functions of error and channel are given, can be written as
We assume that each user feeds back its index to the BS through a zero-delay and error-free feedback channel with B bits. Based on robust design on delay and quantization case, the precoding matrix [20] is
where, Then the achievable robust ergodic capacity is
According to SINR of the received signal (19), ergodic capacity under robust design can be obtained
is the minimum capacity of the worst channel conditions, and is the maximum capacity with precoding in the worst case. Because capacity is positive linear with SINR, then the equivalent objective function can also be referred to as
Even in a single cell downlink transmission, the optimization problem for maxmin SINR robust design has no solution mathematically [17]. Here, the upper and lower bounds of the above objective function (24) will be given .
1) Analysis of the lower bound
At first, taking a relaxation of minimization of inside (24), the resulting lower bound is
Because the three error variables are independent, so the lower bound is compact. According to the lemma given in literature [9], the equivalent expression is
where, ε is the probability of the link error. Similarly, it is possible to obtain
And suppose that the error vector is bounded noise model [22], that is ║E1k[n ]║2F ≤ ε1, ║E2k[n ]║2F ≤ ε2. Therefore, the problem (27) can be rewritten as
With the high-order spread spectrum modulation, the probability of the link error ε is minimal with respect to the channel coefficient. Then and the problem (27) further becomes
Obviously, this problem can be converted to semi-definite programming (SDP) problem for solving by introducing slack variables. Conversion and solving process is as follows.
Introduce slack variables τ, that is
Add another slack variable δ, constraints can be converted to
The above constraints satisfy the convex constraints, so equation (31) is second-order cone programming (SOCP) problem, the standard form for (31) is
This section gives a binary search algorithm to solve this problem. The concrete steps are as follows:
2) Analysis of the upper bound
Theorem 1 : For any function f(x, y), is found.
Proof : If ,
Through the minimax inequality conversion [9] , which is proved by theorem 1 , the upper bound for the problem (24) can be the
Therefore, the upper bound of the problem is
It is easy to know the maximization problem within upper bound can be solved by generalized Rayleigh quotient [23],
Without considering power constraints and denote then the optimal precoding vector is
The maximum value of the target SINR is the maximum generalized eigenvalue for TkH Tk relative to QkH Qk + σ2I/MP . If the power constraint is considered, the power constraint for each base station can be relaxed to the total power constraint, then Wkopt[n ] is need to be multiplied by a power control factor.
After then, the upper bound problem can be simplified as
That means minimizing the maximum generalized eigenvalue of the relative matrix. This matrix is
Due to the optimization eigenvalue problem is more abstract, in order to get a precise mathematical expression, a single-user single data stream is considered. Let R = qqH + σ2I/MP , and the problem (38) is converted into
According to the relative matrix knowledge [23],
Again apply the eigenvalue nature of the singular matrix,
In order to construct a unitary matrix , the matrix qqH can be written as
where, the matrix is unitary matrix. Then,
Inverse matrix R-1 is
Theorem 2 : From(45), there is
Proof:
Since tH t is a constant, the objective function (42) can be further written as
Obviously, if can reach the smallest eigenvalue of the matrix R-1, you can get the optimal solution of the above problem. and the smallest eigenvalues of R-1 can be seen from the formula (45), and it is
Therefore, only with the same direction can satisfy the conditions of the optimal solution, that is and α0 ≤ ρε1 + ε2. The problem then is transformed into
Clearly, when ║q║ = ρε1 + ε2 is maximum which means α0 is taken to the upper limit, the objective function is the minimum. Therefore, the optimal solution is
The equation(46) is established.
4.2. The optimal number of cooperative BSs
The coordinated set planning based on the robust net ergodic capacity is to be considered in a theoretical analysis due to that the solution of optimization problem (22) is very difficult. Based on the robust capacity above, we select capacity gain as a performance comparison indicator.
Robust coordinated set planning analysis is not meaningful for the two - tier cell scene. According to Wyner model [19], the large scale fading elements of A is
where, M is as the dimensions of A, and the scaling factor α (α ∈ [0,1]) is the distance between the position of UE and the cell center, normalized to the maximum distance within a cell.
Specific implementation process of coordinated set planning in this section is as follows:
5. Simulation Results and Discussion
In this section, we start with simulation results that compares upper and lower bounds of robust capacity and ergodic capacity with SNR of robust precoding and other traditional precoding, respectively. Then, the CoMP capacity and non-CoMP capacity with path loss factor are simulated. Finally, the optimal number of cooperative BSs with large scale fading factor is simulated and analyzed.
5.1. Simulation scenario and parameter
The system simulation scenario and parameters are shown as following.
Table 1.Simulation scenario and parameter
According to the COST231 Hata model , the path loss model is
where, hbs, hms are height of BS’s and MS’s antenna, fc is carrier frequency, in units of MHz; d is the straight distance between the BS and MS, C is a constant; parameters for urban macrocell are hbs = 32m, hms = 1.5m, fc = 1900MHz, C = 3dB, the correction model of the path loss is
5.2 Simulation results and discussion
The upper and lower bounds of robust capacity are changing with the SNR in the two cases, which is given in Fig. 3. Upper1, Lower1 corresponding to the case of ρ = 1, ε = ε1 = ε2 = 0.05, and Upper2, Lower2 corresponds to ρ = 0.9, ε = ε1 = ε2 = 0.1. As can be seen from Fig. 3, the upper bound becomes more close to lower bound with the increment of SNR. Hence, the analysis of robust capacity given in this paper is reasonable and meaningful.
Fig. 3.The upper and lower bound of roubst capacity with SNR for two cases (Nt = 2, M = 2, alpha = 0.8)
Fig. 4 shows ergodic capacity with SNR of several traditional precoding and robust precoding. And, the number of feedback bits is 10, the number of cooperative BSs is 2, and the number of antennas of each base station is 2, and feedback delay Doppler product is 0.1. From this figure, the robust design based on maxmin principle is better than ZF and MMSE precoding scheme. In addition, we can also observe that when SNR is high, MMSE precoding converges to ZF precoding, and this drawback can be overcomed by robust precoding.
Fig. 4.The capacity of several tradtional precoding and robust precoding with SNR (Nt = 2, M = 2, alpha = 0.8)
Fig. 5 shows the change curve of the optimal number of cooperative BSs with large scale fading factor in three cases. Long dashed line corresponds to the results of the coordinated set planning only considering the quantization error corresponding to the case of B = 10. Short dashed line corresponds to the coordinated set planning result considering the quantization error and feedback delay corresponding to the case of B = 10, ρ = 0.9. Solid line corresponds to the coordinated set planning result considering the feedback delay, the quantization error, the estimated error and link error corresponding to the case of B = 10 , ρ = 0.9, ε = ε1 = ε2 = 0.1. From this three curves, the greater α, the larger the number of cooperative base stations. This is because the larger α, farther away is the user from the cell center, and the worse of the channel correlation from the adjacent cell, so the inter-cell interference is more obvious, the greater the capacity gain we can obtain CoMP relative to non-CoMP mode. Therefore, the number of cooperative base stations increases. When α is the same , the optimal number of cooperative BSs of the case 3 is greater than the case 2, which is greater than the case 1. Since more practical channel condition is considered in case 3, the cooperative gain is relatively decreased, so that larger numbers of BSs are needed to cooperate to maximize the ergodic capacity under case 3.
Fig. 5.The optimal number of coordinated BSs with path loss factor for three cases (Nt = 2, SNR = 0)
6. Conclusion
CoMP has the potential to realize significant gains in throughput and reliability. In practical systems, perfect BSs cooperation or global processing is very difficult, if not impossible, to achieve. This paper has studied robust design of the multi-cell collaborative coordinated set planning considering the feedback delay and link error. To solve the tradeoff between the advantage and disadvantage of CoMP, a penalty factor was introduced to express the equivalent capacity loss. The net ergodic capacity optimization problem, whose optimization variables were the number of the coordinated BSs and the dividing-area radius, was derived and simplified. By employing the maximum minimum criteria, upper and lower bounds of the robust capacity were investigated. And the gap between two bounds gets smaller as transmission power increases. Based on upper bound of the robust capacity, we defined the cooperative gain, and then have proposed a practical scheme to determine the optimal number of cooperative BSs under the case of each fixed path loss factor.
The robust design based on maxmin principle is better than BF, ZF and MMSE precoding scheme. Besides, in the case that the large scale fading is higher, the greater the capacity gain we can obtain CoMP relative to Non CoMP mode and the number of the cooperated BSs shall be greater; in the case that the large scale fading is the same and the channel is less reliable, the cooperative gain is relatively decreased and larger numbers of BSs are needed to cooperate to maximize the ergodic capacity.
참고문헌
- 3GPP TR 36.814 V1.1.1, "Further Advancements for E-UTRA-Physical Layer Aspects," June 2009.
- 3GPP TS 36.331 V8.2.0, "Evolved Universal Terrestrial Radio Access (E-UTRA) Radio Resource Control," 2008.
- 3GPP TSG RAN WG1 Meeting #57, R1-091919. "Updates on cell clustering for CoMP transmission/reception," 2009.
- R. Irmer, H. Droste, P. Marsch, et al. "Coordinated multipoint: Concepts, performance, and field trial results," Communications Magazine, IEEE, 49(2), pp.102-111, 2011.
- G.Y. Kim, J.A.C. Lee, S.J. Hong, "Analysis of Macro-Diversity in LTE-Advanced," KSII Transactions on Internet and Information Systems, vol. 5(9), pp.1596-1612, 2011.
- E. Bjornson, R. Zakhour, "Cooperative Multicell Precoding : Rate Region Characterization and Distributed Strategies With Instantaneous and Statistical CSI," Transactions on Signal Process., vol. 58(8), pp. 4298-4310, 2010. https://doi.org/10.1109/TSP.2010.2049996
- J. Zhang, R.-h. Chen, J. G. Andrews, A. Ghosh, & R. W., Jr. Heath, "Networked MIMO with Clustered Linear Precoding," IEEE Transactions on Wireless Communications, 8(4), pp.1910-1921,2009. https://doi.org/10.1109/TWC.2009.080180
- N. Vucic and H. Boche, "Robust QoS-constrained optimization of downlink multiuser MISO systems," IEEE Transactions on Signal Process., vol.57 (2), pp. 714-725, 2009. https://doi.org/10.1109/TSP.2008.2008553
- A. Tajer, N. Prasad, and X. Wang, "Robust linear precoder design for multi-cell downlink transmission," IEEE Transactions on Signal process., vol. 59( 1), pp. 235-251, 2011. https://doi.org/10.1109/TSP.2010.2082537
- J.X. Dai, M. Chen, M. Zhao, et al. "Study of Coordinated Set of Coordinated Multi-Point Transmission with Limited Feedback," IEICE Transactions on Communications, 97(1), pp. 171-181, 2014.
- D. Lecompte, F. David, "Evolved multimedia broadcast/multicast service in LTE-advanced: overview and Rel-11 enhancements," Communications Magazine, vol. 50(11), pp. 68-74, 2012.
- D. Gesbert, S. Hanly, H. Huang, S. Shamai, O. Simeone, & W. Yu, "Multi-cell MIMO cooperative networks: a new look at interference," IEEE Journal on Selected Areas in Communications, vol. 28(9), 1380-1408, 2010. https://doi.org/10.1109/JSAC.2010.101202
- J.Q. Ni, Z.S. Fei, C.W. Xing, "Terminal-based Dynamic Clustering Algorithm in Multi-Cell Cellular System," KSII Transactions on Internet and Information Systems, vol. 6(9), pp.2086-2097, 2012.
- S.-h. Xiao, "Discussion on Strategies for Adaptive Dynamical Clustering in Cooperative Multi-point Downlink Transmission Systems," Wireless Personal Communications, vol. 67(3), 525-539, 2012. https://doi.org/10.1007/s11277-011-0394-4
- R.A. Abdelaal, K. Elsayed, & M.H. Ismail, "Cooperative scheduling, precoding, and optimized power allocation for LTE-advanced CoMP systems," In Proceedings 2012 IFIP Wireless Days, pp.1-6, Dublin, 2012.
- C.K. Au-yeung, S. Member, and D.J. Love, "On the Performance of Random Vector Quantization Limited Feedback Beamforming in a MISO System," IEEE Transactions on Wireless Commun., vol. 6(2), pp. 458-462, 2007. https://doi.org/10.1109/TWC.2007.05351
- N. Jindal, "MIMO Broadcast Channels With Finite-Rate Feedback," IEEE Transactions on Information Theory,vol.52(11), pp. 5045-5060, 2006. https://doi.org/10.1109/TIT.2006.883550
- J. Zhang, J.G. Andrews, "Mode Switching for MIMO Broadcast Channel Based on Delay and Channel Quantization," EURASIP Journal on Advances in Signal Processing, vol. 2009, pp. 1-15, 2009.
- J. Hoydis, M. Kobayashi, "On the optimal number of cooperative base stations in network mimo," http://arxiv.org/abs/1003.0332v2.
- C. Zhang, W.Xu, M. Chen, "Robust MMSE Beamforming for Multiuser MISO Systems with Limited Feedback," IEEE Transactions on Signal Process. Letter, vol. 16 (7), pp. 588-591, 2009. https://doi.org/10.1109/LSP.2009.2020455
- X.-Z. Wang, P.-C. Zhu, & M. Chen, "Antenna location design for generalized distributed antenna systems," IEEE Communications Letters, pp. 13(5), 315-317, 2009. https://doi.org/10.1109/LCOMM.2009.090123
- X. Zhang, D. P. Palomar, "Statistically robust design of linear MIMO transceivers," IEEE Transactions on Signal Process. , vol. 56(8), pp. 3678-3689, 2008. https://doi.org/10.1109/TSP.2008.919384
- G. H. Golub and C. F. V. Loan, "Matrix Computations", 3rd Edition, Baltimore, MD: The John Hopkins University Press, 1996.