ergo
template_lapack_getrs.h
Go to the documentation of this file.
1 /* Ergo, version 3.3, a program for linear scaling electronic structure
2  * calculations.
3  * Copyright (C) 2013 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_GETRS_HEADER
36 #define TEMPLATE_LAPACK_GETRS_HEADER
37 
38 
39 template<class Treal>
40 int template_lapack_getrs(const char *trans, const integer *n, const integer *nrhs,
41  const Treal *a, const integer *lda, const integer *ipiv, Treal *b, const integer *
42  ldb, integer *info)
43 {
44 /* -- LAPACK routine (version 3.0) --
45  Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,
46  Courant Institute, Argonne National Lab, and Rice University
47  March 31, 1993
48 
49 
50  Purpose
51  =======
52 
53  DGETRS solves a system of linear equations
54  A * X = B or A' * X = B
55  with a general N-by-N matrix A using the LU factorization computed
56  by DGETRF.
57 
58  Arguments
59  =========
60 
61  TRANS (input) CHARACTER*1
62  Specifies the form of the system of equations:
63  = 'N': A * X = B (No transpose)
64  = 'T': A'* X = B (Transpose)
65  = 'C': A'* X = B (Conjugate transpose = Transpose)
66 
67  N (input) INTEGER
68  The order of the matrix A. N >= 0.
69 
70  NRHS (input) INTEGER
71  The number of right hand sides, i.e., the number of columns
72  of the matrix B. NRHS >= 0.
73 
74  A (input) DOUBLE PRECISION array, dimension (LDA,N)
75  The factors L and U from the factorization A = P*L*U
76  as computed by DGETRF.
77 
78  LDA (input) INTEGER
79  The leading dimension of the array A. LDA >= max(1,N).
80 
81  IPIV (input) INTEGER array, dimension (N)
82  The pivot indices from DGETRF; for 1<=i<=N, row i of the
83  matrix was interchanged with row IPIV(i).
84 
85  B (input/output) DOUBLE PRECISION array, dimension (LDB,NRHS)
86  On entry, the right hand side matrix B.
87  On exit, the solution matrix X.
88 
89  LDB (input) INTEGER
90  The leading dimension of the array B. LDB >= max(1,N).
91 
92  INFO (output) INTEGER
93  = 0: successful exit
94  < 0: if INFO = -i, the i-th argument had an illegal value
95 
96  =====================================================================
97 
98 
99  Test the input parameters.
100 
101  Parameter adjustments */
102  /* Table of constant values */
103  integer c__1 = 1;
104  Treal c_b12 = 1.;
105  integer c_n1 = -1;
106 
107  /* System generated locals */
108  integer a_dim1, a_offset, b_dim1, b_offset, i__1;
109  /* Local variables */
110  logical notran;
111 
112 
113  a_dim1 = *lda;
114  a_offset = 1 + a_dim1 * 1;
115  a -= a_offset;
116  --ipiv;
117  b_dim1 = *ldb;
118  b_offset = 1 + b_dim1 * 1;
119  b -= b_offset;
120 
121  /* Function Body */
122  *info = 0;
123  notran = template_blas_lsame(trans, "N");
124  if (! notran && ! template_blas_lsame(trans, "T") && ! template_blas_lsame(
125  trans, "C")) {
126  *info = -1;
127  } else if (*n < 0) {
128  *info = -2;
129  } else if (*nrhs < 0) {
130  *info = -3;
131  } else if (*lda < maxMACRO(1,*n)) {
132  *info = -5;
133  } else if (*ldb < maxMACRO(1,*n)) {
134  *info = -8;
135  }
136  if (*info != 0) {
137  i__1 = -(*info);
138  template_blas_erbla("GETRS ", &i__1);
139  return 0;
140  }
141 
142 /* Quick return if possible */
143 
144  if (*n == 0 || *nrhs == 0) {
145  return 0;
146  }
147 
148  if (notran) {
149 
150 /* Solve A * X = B.
151 
152  Apply row interchanges to the right hand sides. */
153 
154  template_lapack_laswp(nrhs, &b[b_offset], ldb, &c__1, n, &ipiv[1], &c__1);
155 
156 /* Solve L*X = B, overwriting B with X. */
157 
158  template_blas_trsm("Left", "Lower", "No transpose", "Unit", n, nrhs, &c_b12, &a[
159  a_offset], lda, &b[b_offset], ldb);
160 
161 /* Solve U*X = B, overwriting B with X. */
162 
163  template_blas_trsm("Left", "Upper", "No transpose", "Non-unit", n, nrhs, &c_b12, &
164  a[a_offset], lda, &b[b_offset], ldb);
165  } else {
166 
167 /* Solve A' * X = B.
168 
169  Solve U'*X = B, overwriting B with X. */
170 
171  template_blas_trsm("Left", "Upper", "Transpose", "Non-unit", n, nrhs, &c_b12, &a[
172  a_offset], lda, &b[b_offset], ldb);
173 
174 /* Solve L'*X = B, overwriting B with X. */
175 
176  template_blas_trsm("Left", "Lower", "Transpose", "Unit", n, nrhs, &c_b12, &a[
177  a_offset], lda, &b[b_offset], ldb);
178 
179 /* Apply row interchanges to the solution vectors. */
180 
181  template_lapack_laswp(nrhs, &b[b_offset], ldb, &c__1, n, &ipiv[1], &c_n1);
182  }
183 
184  return 0;
185 
186 /* End of DGETRS */
187 
188 } /* dgetrs_ */
189 
190 #endif
int integer
Definition: template_blas_common.h:38
#define maxMACRO(a, b)
Definition: template_blas_common.h:43
int template_lapack_laswp(const integer *n, Treal *a, const integer *lda, const integer *k1, const integer *k2, const integer *ipiv, const integer *incx)
Definition: template_lapack_laswp.h:40
int template_blas_erbla(const char *srname, integer *info)
Definition: template_blas_common.cc:144
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:40
int template_blas_trsm(const char *side, const char *uplo, const char *transa, const char *diag, const integer *m, const integer *n, const Treal *alpha, const Treal *a, const integer *lda, Treal *b, const integer *ldb)
Definition: template_blas_trsm.h:41
bool logical
Definition: template_blas_common.h:39
logical template_blas_lsame(const char *ca, const char *cb)
Definition: template_blas_common.cc:44