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