00001 /* Ergo, version 3.7, a program for linear scaling electronic structure 00002 * calculations. 00003 * Copyright (C) 2018 Elias Rudberg, Emanuel H. Rubensson, Pawel Salek, 00004 * and Anastasia Kruchinina. 00005 * 00006 * This program is free software: you can redistribute it and/or modify 00007 * it under the terms of the GNU General Public License as published by 00008 * the Free Software Foundation, either version 3 of the License, or 00009 * (at your option) any later version. 00010 * 00011 * This program is distributed in the hope that it will be useful, 00012 * but WITHOUT ANY WARRANTY; without even the implied warranty of 00013 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00014 * GNU General Public License for more details. 00015 * 00016 * You should have received a copy of the GNU General Public License 00017 * along with this program. If not, see <http://www.gnu.org/licenses/>. 00018 * 00019 * Primary academic reference: 00020 * Ergo: An open-source program for linear-scaling electronic structure 00021 * calculations, 00022 * Elias Rudberg, Emanuel H. Rubensson, Pawel Salek, and Anastasia 00023 * Kruchinina, 00024 * SoftwareX 7, 107 (2018), 00025 * <http://dx.doi.org/10.1016/j.softx.2018.03.005> 00026 * 00027 * For further information about Ergo, see <http://www.ergoscf.org>. 00028 */ 00029 00030 /* This file belongs to the template_lapack part of the Ergo source 00031 * code. The source files in the template_lapack directory are modified 00032 * versions of files originally distributed as CLAPACK, see the 00033 * Copyright/license notice in the file template_lapack/COPYING. 00034 */ 00035 00036 00037 #ifndef TEMPLATE_LAPACK_GETF2_HEADER 00038 #define TEMPLATE_LAPACK_GETF2_HEADER 00039 00040 00041 template<class Treal> 00042 int template_lapack_getf2(const integer *m, const integer *n, Treal *a, const integer * 00043 lda, integer *ipiv, integer *info) 00044 { 00045 /* -- LAPACK routine (version 3.0) -- 00046 Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., 00047 Courant Institute, Argonne National Lab, and Rice University 00048 June 30, 1992 00049 00050 00051 Purpose 00052 ======= 00053 00054 DGETF2 computes an LU factorization of a general m-by-n matrix A 00055 using partial pivoting with row interchanges. 00056 00057 The factorization has the form 00058 A = P * L * U 00059 where P is a permutation matrix, L is lower triangular with unit 00060 diagonal elements (lower trapezoidal if m > n), and U is upper 00061 triangular (upper trapezoidal if m < n). 00062 00063 This is the right-looking Level 2 BLAS version of the algorithm. 00064 00065 Arguments 00066 ========= 00067 00068 M (input) INTEGER 00069 The number of rows of the matrix A. M >= 0. 00070 00071 N (input) INTEGER 00072 The number of columns of the matrix A. N >= 0. 00073 00074 A (input/output) DOUBLE PRECISION array, dimension (LDA,N) 00075 On entry, the m by n matrix to be factored. 00076 On exit, the factors L and U from the factorization 00077 A = P*L*U; the unit diagonal elements of L are not stored. 00078 00079 LDA (input) INTEGER 00080 The leading dimension of the array A. LDA >= max(1,M). 00081 00082 IPIV (output) INTEGER array, dimension (min(M,N)) 00083 The pivot indices; for 1 <= i <= min(M,N), row i of the 00084 matrix was interchanged with row IPIV(i). 00085 00086 INFO (output) INTEGER 00087 = 0: successful exit 00088 < 0: if INFO = -k, the k-th argument had an illegal value 00089 > 0: if INFO = k, U(k,k) is exactly zero. The factorization 00090 has been completed, but the factor U is exactly 00091 singular, and division by zero will occur if it is used 00092 to solve a system of equations. 00093 00094 ===================================================================== 00095 00096 00097 Test the input parameters. 00098 00099 Parameter adjustments */ 00100 /* Table of constant values */ 00101 integer c__1 = 1; 00102 Treal c_b6 = -1.; 00103 00104 /* System generated locals */ 00105 integer a_dim1, a_offset, i__1, i__2, i__3; 00106 Treal d__1; 00107 /* Local variables */ 00108 integer j; 00109 integer jp; 00110 #define a_ref(a_1,a_2) a[(a_2)*a_dim1 + a_1] 00111 00112 00113 a_dim1 = *lda; 00114 a_offset = 1 + a_dim1 * 1; 00115 a -= a_offset; 00116 --ipiv; 00117 00118 /* Function Body */ 00119 *info = 0; 00120 if (*m < 0) { 00121 *info = -1; 00122 } else if (*n < 0) { 00123 *info = -2; 00124 } else if (*lda < maxMACRO(1,*m)) { 00125 *info = -4; 00126 } 00127 if (*info != 0) { 00128 i__1 = -(*info); 00129 template_blas_erbla("GETF2 ", &i__1); 00130 return 0; 00131 } 00132 00133 /* Quick return if possible */ 00134 00135 if (*m == 0 || *n == 0) { 00136 return 0; 00137 } 00138 00139 i__1 = minMACRO(*m,*n); 00140 for (j = 1; j <= i__1; ++j) { 00141 00142 /* Find pivot and test for singularity. */ 00143 00144 i__2 = *m - j + 1; 00145 jp = j - 1 + template_blas_idamax(&i__2, &a_ref(j, j), &c__1); 00146 ipiv[j] = jp; 00147 if (a_ref(jp, j) != 0.) { 00148 00149 /* Apply the interchange to columns 1:N. */ 00150 00151 if (jp != j) { 00152 template_blas_swap(n, &a_ref(j, 1), lda, &a_ref(jp, 1), lda); 00153 } 00154 00155 /* Compute elements J+1:M of J-th column. */ 00156 00157 if (j < *m) { 00158 i__2 = *m - j; 00159 d__1 = 1. / a_ref(j, j); 00160 template_blas_scal(&i__2, &d__1, &a_ref(j + 1, j), &c__1); 00161 } 00162 00163 } else if (*info == 0) { 00164 00165 *info = j; 00166 } 00167 00168 if (j < minMACRO(*m,*n)) { 00169 00170 /* Update trailing submatrix. */ 00171 00172 i__2 = *m - j; 00173 i__3 = *n - j; 00174 template_blas_ger(&i__2, &i__3, &c_b6, &a_ref(j + 1, j), &c__1, &a_ref(j, j + 00175 1), lda, &a_ref(j + 1, j + 1), lda); 00176 } 00177 /* L10: */ 00178 } 00179 return 0; 00180 00181 /* End of DGETF2 */ 00182 00183 } /* dgetf2_ */ 00184 00185 #undef a_ref 00186 00187 00188 #endif