TBTK
LUSolver.h
Go to the documentation of this file.
1 /* Copyright 2017 Kristofer Bj√∂rnson
2  *
3  * Licensed under the Apache License, Version 2.0 (the "License");
4  * you may not use this file except in compliance with the License.
5  * You may obtain a copy of the License at
6  *
7  * http://www.apache.org/licenses/LICENSE-2.0
8  *
9  * Unless required by applicable law or agreed to in writing, software
10  * distributed under the License is distributed on an "AS IS" BASIS,
11  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
12  * See the License for the specific language governing permissions and
13  * limitations under the License.
14  */
15 
23 #ifndef COM_DAFER45_TBTK_LU_SOLVER
24 #define COM_DAFER45_TBTK_LU_SOLVER
25 
26 #include "TBTK/Communicator.h"
27 #include "TBTK/Matrix.h"
28 #include "TBTK/SparseMatrix.h"
29 
30 #include <complex>
31 
32 #include "slu_zdefs.h"
33 
34 namespace TBTK{
35 
37 class LUSolver : public Communicator{
38 public:
42  enum class DataType {None, Double, ComplexDouble};
43 
45  LUSolver();
46 
48  ~LUSolver();
49 
51  void setMatrix(const SparseMatrix<double> &sparseMatrix);
52 
54  void setMatrix(const SparseMatrix<std::complex<double>> &sparseMatrix);
55 
58 
60  void solve(Matrix<double> &b);
61 
63  void solve(Matrix<std::complex<double>> &b);
64 private:
66  SuperMatrix *L;
67 
69  SuperMatrix *U;
70 
72  int *rowPermutations;
73 
75  int *columnPermutations;
76 
78  SuperLUStat_t *statistics;
79 
81  DataType matrixDataType;
82 
84  void allocatePermutationMatrices(
85  unsigned int numRows,
86  unsigned int numColumns
87  );
88 
90  void initStatistics();
91 
93  void allocateLUMatrices();
94 
96  void initOptionsAndPermutationMatrices(
97  superlu_options_t &options,
98  SuperMatrix &matrix
99  );
100 
102  void performLUFactorization(SuperMatrix &matrix);
103 
105  void checkSolveAssert(unsigned int numRows);
106 
108  void checkXgstrfErrors(
109  int info,
110  std::string functionName,
111  int numColumns
112  );
113 
115  void checkXgstrsErrors(int info, std::string functionName);
116 };
117 
119  return matrixDataType;
120 }
121 
122 }; //End of namespace TBTK
123 
124 #endif
void setMatrix(const SparseMatrix< double > &sparseMatrix)
Sparse matrix.
DataType getMatrixDataType() const
Definition: LUSolver.h:118
Custom matrix.
Definition: Matrix.h:33
Definition: ModelFactory.h:35
Definition: SparseMatrix.h:35
Base class that communicate their current status during execution.
Solves Hx = y for x, where H is a SparseMatrix.
Definition: LUSolver.h:37
Definition: Communicator.h:28
void solve(Matrix< double > &b)
DataType
Definition: LUSolver.h:42