ergo
template_lapack_gesv.h
Go to the documentation of this file.
1 /* Ergo, version 3.7, a program for linear scaling electronic structure
2  * calculations.
3  * Copyright (C) 2018 Elias Rudberg, Emanuel H. Rubensson, Pawel Salek,
4  * and Anastasia Kruchinina.
5  *
6  * This program is free software: you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License as published by
8  * the Free Software Foundation, either version 3 of the License, or
9  * (at your option) any later version.
10  *
11  * This program is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  * GNU General Public License for more details.
15  *
16  * You should have received a copy of the GNU General Public License
17  * along with this program. If not, see <http://www.gnu.org/licenses/>.
18  *
19  * Primary academic reference:
20  * Ergo: An open-source program for linear-scaling electronic structure
21  * calculations,
22  * Elias Rudberg, Emanuel H. Rubensson, Pawel Salek, and Anastasia
23  * Kruchinina,
24  * SoftwareX 7, 107 (2018),
25  * <http://dx.doi.org/10.1016/j.softx.2018.03.005>
26  *
27  * For further information about Ergo, see <http://www.ergoscf.org>.
28  */
29 
30  /* This file belongs to the template_lapack part of the Ergo source
31  * code. The source files in the template_lapack directory are modified
32  * versions of files originally distributed as CLAPACK, see the
33  * Copyright/license notice in the file template_lapack/COPYING.
34  */
35 
36 
37 #ifndef TEMPLATE_LAPACK_GESV_HEADER
38 #define TEMPLATE_LAPACK_GESV_HEADER
39 
40 
41 template<class Treal>
42 int template_lapack_gesv(const integer *n, const integer *nrhs, Treal *a, const integer
43  *lda, integer *ipiv, Treal *b, const integer *ldb, integer *info)
44 {
45 /* -- LAPACK driver routine (version 3.0) --
46  Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,
47  Courant Institute, Argonne National Lab, and Rice University
48  March 31, 1993
49 
50 
51  Purpose
52  =======
53 
54  DGESV computes the solution to a real system of linear equations
55  A * X = B,
56  where A is an N-by-N matrix and X and B are N-by-NRHS matrices.
57 
58  The LU decomposition with partial pivoting and row interchanges is
59  used to factor A as
60  A = P * L * U,
61  where P is a permutation matrix, L is unit lower triangular, and U is
62  upper triangular. The factored form of A is then used to solve the
63  system of equations A * X = B.
64 
65  Arguments
66  =========
67 
68  N (input) INTEGER
69  The number of linear equations, i.e., the order of the
70  matrix A. N >= 0.
71 
72  NRHS (input) INTEGER
73  The number of right hand sides, i.e., the number of columns
74  of the matrix B. NRHS >= 0.
75 
76  A (input/output) DOUBLE PRECISION array, dimension (LDA,N)
77  On entry, the N-by-N coefficient matrix A.
78  On exit, the factors L and U from the factorization
79  A = P*L*U; the unit diagonal elements of L are not stored.
80 
81  LDA (input) INTEGER
82  The leading dimension of the array A. LDA >= max(1,N).
83 
84  IPIV (output) INTEGER array, dimension (N)
85  The pivot indices that define the permutation matrix P;
86  row i of the matrix was interchanged with row IPIV(i).
87 
88  B (input/output) DOUBLE PRECISION array, dimension (LDB,NRHS)
89  On entry, the N-by-NRHS matrix of right hand side matrix B.
90  On exit, if INFO = 0, the N-by-NRHS solution matrix X.
91 
92  LDB (input) INTEGER
93  The leading dimension of the array B. LDB >= max(1,N).
94 
95  INFO (output) INTEGER
96  = 0: successful exit
97  < 0: if INFO = -i, the i-th argument had an illegal value
98  > 0: if INFO = i, U(i,i) is exactly zero. The factorization
99  has been completed, but the factor U is exactly
100  singular, so the solution could not be computed.
101 
102  =====================================================================
103 
104 
105  Test the input parameters.
106 
107  Parameter adjustments */
108  /* System generated locals */
109  integer a_dim1, a_offset, b_dim1, b_offset, i__1;
110  /* Local variables */
111 
112  a_dim1 = *lda;
113  a_offset = 1 + a_dim1 * 1;
114  a -= a_offset;
115  --ipiv;
116  b_dim1 = *ldb;
117  b_offset = 1 + b_dim1 * 1;
118  b -= b_offset;
119 
120  /* Function Body */
121  *info = 0;
122  if (*n < 0) {
123  *info = -1;
124  } else if (*nrhs < 0) {
125  *info = -2;
126  } else if (*lda < maxMACRO(1,*n)) {
127  *info = -4;
128  } else if (*ldb < maxMACRO(1,*n)) {
129  *info = -7;
130  }
131  if (*info != 0) {
132  i__1 = -(*info);
133  template_blas_erbla("GESV ", &i__1);
134  return 0;
135  }
136 
137 /* Compute the LU factorization of A. */
138 
139  template_lapack_getrf(n, n, &a[a_offset], lda, &ipiv[1], info);
140  if (*info == 0) {
141 
142 /* Solve the system A*X = B, overwriting B with X. */
143 
144  template_lapack_getrs("No transpose", n, nrhs, &a[a_offset], lda, &ipiv[1], &b[
145  b_offset], ldb, info);
146  }
147  return 0;
148 
149 /* End of DGESV */
150 
151 } /* dgesv_ */
152 
153 #endif
int integer
Definition: template_blas_common.h:40
#define maxMACRO(a, b)
Definition: template_blas_common.h:45
int template_blas_erbla(const char *srname, integer *info)
Definition: template_blas_common.cc:146
int template_lapack_getrs(const char *trans, const integer *n, const integer *nrhs, const Treal *a, const integer *lda, const integer *ipiv, Treal *b, const integer *ldb, integer *info)
Definition: template_lapack_getrs.h:42
int template_lapack_getrf(const integer *m, const integer *n, Treal *a, const integer *lda, integer *ipiv, integer *info)
Definition: template_lapack_getrf.h:42
int template_lapack_gesv(const integer *n, const integer *nrhs, Treal *a, const integer *lda, integer *ipiv, Treal *b, const integer *ldb, integer *info)
Definition: template_lapack_gesv.h:42