00001 /* Ergo, version 3.2, a program for linear scaling electronic structure 00002 * calculations. 00003 * Copyright (C) 2012 Elias Rudberg, Emanuel H. Rubensson, and Pawel Salek. 00004 * 00005 * This program is free software: you can redistribute it and/or modify 00006 * it under the terms of the GNU General Public License as published by 00007 * the Free Software Foundation, either version 3 of the License, or 00008 * (at your option) any later version. 00009 * 00010 * This program is distributed in the hope that it will be useful, 00011 * but WITHOUT ANY WARRANTY; without even the implied warranty of 00012 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00013 * GNU General Public License for more details. 00014 * 00015 * You should have received a copy of the GNU General Public License 00016 * along with this program. If not, see <http://www.gnu.org/licenses/>. 00017 * 00018 * Primary academic reference: 00019 * KohnâSham Density Functional Theory Electronic Structure Calculations 00020 * with Linearly Scaling Computational Time and Memory Usage, 00021 * Elias Rudberg, Emanuel H. Rubensson, and Pawel Salek, 00022 * J. Chem. Theory Comput. 7, 340 (2011), 00023 * <http://dx.doi.org/10.1021/ct100611z> 00024 * 00025 * For further information about Ergo, see <http://www.ergoscf.org>. 00026 */ 00027 00028 /* This file belongs to the template_lapack part of the Ergo source 00029 * code. The source files in the template_lapack directory are modified 00030 * versions of files originally distributed as CLAPACK, see the 00031 * Copyright/license notice in the file template_lapack/COPYING. 00032 */ 00033 00034 00035 #ifndef TEMPLATE_LAPACK_GEQR2_HEADER 00036 #define TEMPLATE_LAPACK_GEQR2_HEADER 00037 00038 00039 template<class Treal> 00040 int template_lapack_geqr2(const integer *m, const integer *n, Treal *a, const integer * 00041 lda, Treal *tau, Treal *work, integer *info) 00042 { 00043 /* -- LAPACK routine (version 3.0) -- 00044 Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., 00045 Courant Institute, Argonne National Lab, and Rice University 00046 February 29, 1992 00047 00048 00049 Purpose 00050 ======= 00051 00052 DGEQR2 computes a QR factorization of a real m by n matrix A: 00053 A = Q * R. 00054 00055 Arguments 00056 ========= 00057 00058 M (input) INTEGER 00059 The number of rows of the matrix A. M >= 0. 00060 00061 N (input) INTEGER 00062 The number of columns of the matrix A. N >= 0. 00063 00064 A (input/output) DOUBLE PRECISION array, dimension (LDA,N) 00065 On entry, the m by n matrix A. 00066 On exit, the elements on and above the diagonal of the array 00067 contain the min(m,n) by n upper trapezoidal matrix R (R is 00068 upper triangular if m >= n); the elements below the diagonal, 00069 with the array TAU, represent the orthogonal matrix Q as a 00070 product of elementary reflectors (see Further Details). 00071 00072 LDA (input) INTEGER 00073 The leading dimension of the array A. LDA >= max(1,M). 00074 00075 TAU (output) DOUBLE PRECISION array, dimension (min(M,N)) 00076 The scalar factors of the elementary reflectors (see Further 00077 Details). 00078 00079 WORK (workspace) DOUBLE PRECISION array, dimension (N) 00080 00081 INFO (output) INTEGER 00082 = 0: successful exit 00083 < 0: if INFO = -i, the i-th argument had an illegal value 00084 00085 Further Details 00086 =============== 00087 00088 The matrix Q is represented as a product of elementary reflectors 00089 00090 Q = H(1) H(2) . . . H(k), where k = min(m,n). 00091 00092 Each H(i) has the form 00093 00094 H(i) = I - tau * v * v' 00095 00096 where tau is a real scalar, and v is a real vector with 00097 v(1:i-1) = 0 and v(i) = 1; v(i+1:m) is stored on exit in A(i+1:m,i), 00098 and tau in TAU(i). 00099 00100 ===================================================================== 00101 00102 00103 Test the input arguments 00104 00105 Parameter adjustments */ 00106 /* Table of constant values */ 00107 integer c__1 = 1; 00108 00109 /* System generated locals */ 00110 integer a_dim1, a_offset, i__1, i__2, i__3; 00111 /* Local variables */ 00112 integer i__, k; 00113 Treal aii; 00114 #define a_ref(a_1,a_2) a[(a_2)*a_dim1 + a_1] 00115 00116 00117 a_dim1 = *lda; 00118 a_offset = 1 + a_dim1 * 1; 00119 a -= a_offset; 00120 --tau; 00121 --work; 00122 00123 /* Function Body */ 00124 *info = 0; 00125 if (*m < 0) { 00126 *info = -1; 00127 } else if (*n < 0) { 00128 *info = -2; 00129 } else if (*lda < maxMACRO(1,*m)) { 00130 *info = -4; 00131 } 00132 if (*info != 0) { 00133 i__1 = -(*info); 00134 template_blas_erbla("GEQR2 ", &i__1); 00135 return 0; 00136 } 00137 00138 k = minMACRO(*m,*n); 00139 00140 i__1 = k; 00141 for (i__ = 1; i__ <= i__1; ++i__) { 00142 00143 /* Generate elementary reflector H(i) to annihilate A(i+1:m,i) 00144 00145 Computing MIN */ 00146 i__2 = i__ + 1; 00147 i__3 = *m - i__ + 1; 00148 template_lapack_larfg(&i__3, &a_ref(i__, i__), &a_ref(minMACRO(i__2,*m), i__), &c__1, & 00149 tau[i__]); 00150 if (i__ < *n) { 00151 00152 /* Apply H(i) to A(i:m,i+1:n) from the left */ 00153 00154 aii = a_ref(i__, i__); 00155 a_ref(i__, i__) = 1.; 00156 i__2 = *m - i__ + 1; 00157 i__3 = *n - i__; 00158 template_lapack_larf("Left", &i__2, &i__3, &a_ref(i__, i__), &c__1, &tau[i__], & 00159 a_ref(i__, i__ + 1), lda, &work[1]); 00160 a_ref(i__, i__) = aii; 00161 } 00162 /* L10: */ 00163 } 00164 return 0; 00165 00166 /* End of DGEQR2 */ 00167 00168 } /* dgeqr2_ */ 00169 00170 #undef a_ref 00171 00172 00173 #endif