ergo
|
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_GETF2_HEADER 00036 #define TEMPLATE_LAPACK_GETF2_HEADER 00037 00038 00039 template<class Treal> 00040 int template_lapack_getf2(const integer *m, const integer *n, Treal *a, const integer * 00041 lda, integer *ipiv, 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 June 30, 1992 00047 00048 00049 Purpose 00050 ======= 00051 00052 DGETF2 computes an LU factorization of a general m-by-n matrix A 00053 using partial pivoting with row interchanges. 00054 00055 The factorization has the form 00056 A = P * L * U 00057 where P is a permutation matrix, L is lower triangular with unit 00058 diagonal elements (lower trapezoidal if m > n), and U is upper 00059 triangular (upper trapezoidal if m < n). 00060 00061 This is the right-looking Level 2 BLAS version of the algorithm. 00062 00063 Arguments 00064 ========= 00065 00066 M (input) INTEGER 00067 The number of rows of the matrix A. M >= 0. 00068 00069 N (input) INTEGER 00070 The number of columns of the matrix A. N >= 0. 00071 00072 A (input/output) DOUBLE PRECISION array, dimension (LDA,N) 00073 On entry, the m by n matrix to be factored. 00074 On exit, the factors L and U from the factorization 00075 A = P*L*U; the unit diagonal elements of L are not stored. 00076 00077 LDA (input) INTEGER 00078 The leading dimension of the array A. LDA >= max(1,M). 00079 00080 IPIV (output) INTEGER array, dimension (min(M,N)) 00081 The pivot indices; for 1 <= i <= min(M,N), row i of the 00082 matrix was interchanged with row IPIV(i). 00083 00084 INFO (output) INTEGER 00085 = 0: successful exit 00086 < 0: if INFO = -k, the k-th argument had an illegal value 00087 > 0: if INFO = k, U(k,k) is exactly zero. The factorization 00088 has been completed, but the factor U is exactly 00089 singular, and division by zero will occur if it is used 00090 to solve a system of equations. 00091 00092 ===================================================================== 00093 00094 00095 Test the input parameters. 00096 00097 Parameter adjustments */ 00098 /* Table of constant values */ 00099 integer c__1 = 1; 00100 Treal c_b6 = -1.; 00101 00102 /* System generated locals */ 00103 integer a_dim1, a_offset, i__1, i__2, i__3; 00104 Treal d__1; 00105 /* Local variables */ 00106 integer j; 00107 integer jp; 00108 #define a_ref(a_1,a_2) a[(a_2)*a_dim1 + a_1] 00109 00110 00111 a_dim1 = *lda; 00112 a_offset = 1 + a_dim1 * 1; 00113 a -= a_offset; 00114 --ipiv; 00115 00116 /* Function Body */ 00117 *info = 0; 00118 if (*m < 0) { 00119 *info = -1; 00120 } else if (*n < 0) { 00121 *info = -2; 00122 } else if (*lda < maxMACRO(1,*m)) { 00123 *info = -4; 00124 } 00125 if (*info != 0) { 00126 i__1 = -(*info); 00127 template_blas_erbla("GETF2 ", &i__1); 00128 return 0; 00129 } 00130 00131 /* Quick return if possible */ 00132 00133 if (*m == 0 || *n == 0) { 00134 return 0; 00135 } 00136 00137 i__1 = minMACRO(*m,*n); 00138 for (j = 1; j <= i__1; ++j) { 00139 00140 /* Find pivot and test for singularity. */ 00141 00142 i__2 = *m - j + 1; 00143 jp = j - 1 + template_blas_idamax(&i__2, &a_ref(j, j), &c__1); 00144 ipiv[j] = jp; 00145 if (a_ref(jp, j) != 0.) { 00146 00147 /* Apply the interchange to columns 1:N. */ 00148 00149 if (jp != j) { 00150 template_blas_swap(n, &a_ref(j, 1), lda, &a_ref(jp, 1), lda); 00151 } 00152 00153 /* Compute elements J+1:M of J-th column. */ 00154 00155 if (j < *m) { 00156 i__2 = *m - j; 00157 d__1 = 1. / a_ref(j, j); 00158 template_blas_scal(&i__2, &d__1, &a_ref(j + 1, j), &c__1); 00159 } 00160 00161 } else if (*info == 0) { 00162 00163 *info = j; 00164 } 00165 00166 if (j < minMACRO(*m,*n)) { 00167 00168 /* Update trailing submatrix. */ 00169 00170 i__2 = *m - j; 00171 i__3 = *n - j; 00172 template_blas_ger(&i__2, &i__3, &c_b6, &a_ref(j + 1, j), &c__1, &a_ref(j, j + 00173 1), lda, &a_ref(j + 1, j + 1), lda); 00174 } 00175 /* L10: */ 00176 } 00177 return 0; 00178 00179 /* End of DGETF2 */ 00180 00181 } /* dgetf2_ */ 00182 00183 #undef a_ref 00184 00185 00186 #endif