template_lapack_geqr2.h

Go to the documentation of this file.
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

Generated on Wed Nov 21 09:32:11 2012 for ergo by  doxygen 1.4.7