ergo
Main Page
Namespaces
Classes
Files
File List
File Members
template_lapack_getf2.h
Go to the documentation of this file.
1
/* Ergo, version 3.2, a program for linear scaling electronic structure
2
* calculations.
3
* Copyright (C) 2012 Elias Rudberg, Emanuel H. Rubensson, and Pawel Salek.
4
*
5
* This program is free software: you can redistribute it and/or modify
6
* it under the terms of the GNU General Public License as published by
7
* the Free Software Foundation, either version 3 of the License, or
8
* (at your option) any later version.
9
*
10
* This program is distributed in the hope that it will be useful,
11
* but WITHOUT ANY WARRANTY; without even the implied warranty of
12
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13
* GNU General Public License for more details.
14
*
15
* You should have received a copy of the GNU General Public License
16
* along with this program. If not, see <http://www.gnu.org/licenses/>.
17
*
18
* Primary academic reference:
19
* KohnâSham Density Functional Theory Electronic Structure Calculations
20
* with Linearly Scaling Computational Time and Memory Usage,
21
* Elias Rudberg, Emanuel H. Rubensson, and Pawel Salek,
22
* J. Chem. Theory Comput. 7, 340 (2011),
23
* <http://dx.doi.org/10.1021/ct100611z>
24
*
25
* For further information about Ergo, see <http://www.ergoscf.org>.
26
*/
27
28
/* This file belongs to the template_lapack part of the Ergo source
29
* code. The source files in the template_lapack directory are modified
30
* versions of files originally distributed as CLAPACK, see the
31
* Copyright/license notice in the file template_lapack/COPYING.
32
*/
33
34
35
#ifndef TEMPLATE_LAPACK_GETF2_HEADER
36
#define TEMPLATE_LAPACK_GETF2_HEADER
37
38
39
template
<
class
Treal>
40
int
template_lapack_getf2
(
const
integer
*m,
const
integer
*n, Treal *a,
const
integer
*
41
lda,
integer
*ipiv,
integer
*info)
42
{
43
/* -- LAPACK routine (version 3.0) --
44
Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,
45
Courant Institute, Argonne National Lab, and Rice University
46
June 30, 1992
47
48
49
Purpose
50
=======
51
52
DGETF2 computes an LU factorization of a general m-by-n matrix A
53
using partial pivoting with row interchanges.
54
55
The factorization has the form
56
A = P * L * U
57
where P is a permutation matrix, L is lower triangular with unit
58
diagonal elements (lower trapezoidal if m > n), and U is upper
59
triangular (upper trapezoidal if m < n).
60
61
This is the right-looking Level 2 BLAS version of the algorithm.
62
63
Arguments
64
=========
65
66
M (input) INTEGER
67
The number of rows of the matrix A. M >= 0.
68
69
N (input) INTEGER
70
The number of columns of the matrix A. N >= 0.
71
72
A (input/output) DOUBLE PRECISION array, dimension (LDA,N)
73
On entry, the m by n matrix to be factored.
74
On exit, the factors L and U from the factorization
75
A = P*L*U; the unit diagonal elements of L are not stored.
76
77
LDA (input) INTEGER
78
The leading dimension of the array A. LDA >= max(1,M).
79
80
IPIV (output) INTEGER array, dimension (min(M,N))
81
The pivot indices; for 1 <= i <= min(M,N), row i of the
82
matrix was interchanged with row IPIV(i).
83
84
INFO (output) INTEGER
85
= 0: successful exit
86
< 0: if INFO = -k, the k-th argument had an illegal value
87
> 0: if INFO = k, U(k,k) is exactly zero. The factorization
88
has been completed, but the factor U is exactly
89
singular, and division by zero will occur if it is used
90
to solve a system of equations.
91
92
=====================================================================
93
94
95
Test the input parameters.
96
97
Parameter adjustments */
98
/* Table of constant values */
99
integer
c__1 = 1;
100
Treal c_b6 = -1.;
101
102
/* System generated locals */
103
integer
a_dim1, a_offset, i__1, i__2, i__3;
104
Treal d__1;
105
/* Local variables */
106
integer
j;
107
integer
jp;
108
#define a_ref(a_1,a_2) a[(a_2)*a_dim1 + a_1]
109
110
111
a_dim1 = *lda;
112
a_offset = 1 + a_dim1 * 1;
113
a -= a_offset;
114
--ipiv;
115
116
/* Function Body */
117
*info = 0;
118
if
(*m < 0) {
119
*info = -1;
120
}
else
if
(*n < 0) {
121
*info = -2;
122
}
else
if
(*lda <
maxMACRO
(1,*m)) {
123
*info = -4;
124
}
125
if
(*info != 0) {
126
i__1 = -(*info);
127
template_blas_erbla
(
"GETF2 "
, &i__1);
128
return
0;
129
}
130
131
/* Quick return if possible */
132
133
if
(*m == 0 || *n == 0) {
134
return
0;
135
}
136
137
i__1 =
minMACRO
(*m,*n);
138
for
(j = 1; j <= i__1; ++j) {
139
140
/* Find pivot and test for singularity. */
141
142
i__2 = *m - j + 1;
143
jp = j - 1 +
template_blas_idamax
(&i__2, &
a_ref
(j, j), &c__1);
144
ipiv[j] = jp;
145
if
(
a_ref
(jp, j) != 0.) {
146
147
/* Apply the interchange to columns 1:N. */
148
149
if
(jp != j) {
150
template_blas_swap
(n, &
a_ref
(j, 1), lda, &
a_ref
(jp, 1), lda);
151
}
152
153
/* Compute elements J+1:M of J-th column. */
154
155
if
(j < *m) {
156
i__2 = *m - j;
157
d__1 = 1. /
a_ref
(j, j);
158
template_blas_scal
(&i__2, &d__1, &
a_ref
(j + 1, j), &c__1);
159
}
160
161
}
else
if
(*info == 0) {
162
163
*info = j;
164
}
165
166
if
(j <
minMACRO
(*m,*n)) {
167
168
/* Update trailing submatrix. */
169
170
i__2 = *m - j;
171
i__3 = *n - j;
172
template_blas_ger
(&i__2, &i__3, &c_b6, &
a_ref
(j + 1, j), &c__1, &
a_ref
(j, j +
173
1), lda, &
a_ref
(j + 1, j + 1), lda);
174
}
175
/* L10: */
176
}
177
return
0;
178
179
/* End of DGETF2 */
180
181
}
/* dgetf2_ */
182
183
#undef a_ref
184
185
186
#endif
source
matrix
template_lapack
lapack
template_lapack_getf2.h
Generated on Wed Nov 21 2012 20:40:23 for ergo by
1.8.1.2