Nektar++
Static Public Member Functions | List of all members
Nektar::LinearSystemSolver Struct Reference

#include <NekLinSys.hpp>

Static Public Member Functions

template<typename BVectorType , typename XVectorType >
static void Solve (const BVectorType &b, XVectorType &x, MatrixStorage m_matrixType, const Array< OneD, const int > &m_ipivot, unsigned int n, const Array< OneD, const double > &A, char m_transposeFlag, unsigned int m_numberOfSubDiagonals, unsigned int m_numberOfSuperDiagonals)
 
template<typename BVectorType , typename XVectorType >
static void SolveTranspose (const BVectorType &b, XVectorType &x, MatrixStorage m_matrixType, const Array< OneD, const int > &m_ipivot, unsigned int n, const Array< OneD, const double > &A, char m_transposeFlag, unsigned int m_numberOfSubDiagonals, unsigned int m_numberOfSuperDiagonals)
 

Detailed Description

Definition at line 69 of file NekLinSys.hpp.

Member Function Documentation

◆ Solve()

template<typename BVectorType , typename XVectorType >
static void Nektar::LinearSystemSolver::Solve ( const BVectorType &  b,
XVectorType &  x,
MatrixStorage  m_matrixType,
const Array< OneD, const int > &  m_ipivot,
unsigned int  n,
const Array< OneD, const double > &  A,
char  m_transposeFlag,
unsigned int  m_numberOfSubDiagonals,
unsigned int  m_numberOfSuperDiagonals 
)
inlinestatic

Definition at line 72 of file NekLinSys.hpp.

References ASSERTL0, Lapack::Dgbtrs(), Lapack::Dgetrs(), Lapack::Dpbtrs(), Lapack::Dpptrs(), Lapack::Dsptrs(), Lapack::Dtptrs(), Nektar::eBANDED, Nektar::eDIAGONAL, Nektar::ErrorUtil::efatal, Nektar::eFULL, Nektar::eLOWER_TRIANGULAR, Nektar::eLOWER_TRIANGULAR_BANDED, Nektar::ePOSITIVE_DEFINITE_SYMMETRIC, Nektar::ePOSITIVE_DEFINITE_SYMMETRIC_BANDED, Nektar::eSYMMETRIC, Nektar::eSYMMETRIC_BANDED, Nektar::eUPPER_TRIANGULAR, Nektar::eUPPER_TRIANGULAR_BANDED, and NEKERROR.

Referenced by Nektar::LinearSystem::Solve().

77  {
78  switch(m_matrixType)
79  {
80  case eFULL:
81  {
82  x = b;
83  int info = 0;
84  Lapack::Dgetrs('N',n,1,A.get(),n,(int *)m_ipivot.get(),x.GetRawPtr(),n,info);
85  if( info < 0 )
86  {
87  std::string message = "ERROR: The " + std::to_string(-info) + "th parameter had an illegal parameter for dgetrs";
88  ASSERTL0(false, message.c_str());
89  }
90 
91  }
92  break;
93  case eDIAGONAL:
94  for(unsigned int i = 0; i < A.num_elements(); ++i)
95  {
96  x[i] = b[i]*A[i];
97  }
98  break;
99  case eUPPER_TRIANGULAR:
100  {
101  x = b;
102  int info = 0;
103  Lapack::Dtptrs('U', m_transposeFlag, 'N', n, 1, A.get(), x.GetRawPtr(), n, info);
104 
105  if( info < 0 )
106  {
107  std::string message = "ERROR: The " + std::to_string(-info) + "th parameter had an illegal parameter for dtrtrs";
108  ASSERTL0(false, message.c_str());
109  }
110  else if( info > 0 )
111  {
112  std::string message = "ERROR: The " + std::to_string(-info) + "th diagonal element of A is 0 for dtrtrs";
113  ASSERTL0(false, message.c_str());
114  }
115  }
116  break;
117  case eLOWER_TRIANGULAR:
118  {
119  x = b;
120  int info = 0;
121  Lapack::Dtptrs('L', m_transposeFlag, 'N', n, 1, A.get(), x.GetRawPtr(), n, info);
122 
123  if( info < 0 )
124  {
125  std::string message = "ERROR: The " + std::to_string(-info) + "th parameter had an illegal parameter for dtrtrs";
126  ASSERTL0(false, message.c_str());
127  }
128  else if( info > 0 )
129  {
130  std::string message = "ERROR: The " + std::to_string(-info) + "th diagonal element of A is 0 for dtrtrs";
131  ASSERTL0(false, message.c_str());
132  }
133  }
134  break;
135  case eSYMMETRIC:
136  {
137  x = b;
138  int info = 0;
139  Lapack::Dsptrs('U', n, 1, A.get(), m_ipivot.get(), x.GetRawPtr(), x.GetRows(), info);
140  if( info < 0 )
141  {
142  std::string message = "ERROR: The " + std::to_string(-info) + "th parameter had an illegal parameter for dsptrs";
143  ASSERTL0(false, message.c_str());
144  }
145  }
146  break;
147  case ePOSITIVE_DEFINITE_SYMMETRIC:
148  {
149  x = b;
150  int info = 0;
151  Lapack::Dpptrs('U', n, 1, A.get(), x.GetRawPtr(), x.GetRows(), info);
152  if( info < 0 )
153  {
154  std::string message = "ERROR: The " + std::to_string(-info) + "th parameter had an illegal parameter for dpptrs";
155  ASSERTL0(false, message.c_str());
156  }
157  }
158  break;
159  case eBANDED:
160  {
161  x = b;
162  int KL = m_numberOfSubDiagonals;
163  int KU = m_numberOfSuperDiagonals;
164  int info = 0;
165 
166  Lapack::Dgbtrs(m_transposeFlag, n, KL, KU, 1, A.get(), 2*KL+KU+1, m_ipivot.get(), x.GetRawPtr(), n, info);
167 
168  if( info < 0 )
169  {
170  std::string message = "ERROR: The " + std::to_string(-info) + "th parameter had an illegal parameter for dgbtrs";
171  ASSERTL0(false, message.c_str());
172  }
173  }
174  break;
175  case ePOSITIVE_DEFINITE_SYMMETRIC_BANDED:
176  {
177  x = b;
178  int KU = m_numberOfSuperDiagonals;
179  int info = 0;
180 
181  Lapack::Dpbtrs('U', n, KU, 1, A.get(), KU+1, x.GetRawPtr(), n, info);
182 
183  if( info < 0 )
184  {
185  std::string message = "ERROR: The " + std::to_string(-info) + "th parameter had an illegal parameter for dpbtrs";
186  ASSERTL0(false, message.c_str());
187  }
188  }
189  break;
190  case eSYMMETRIC_BANDED:
191  NEKERROR(ErrorUtil::efatal, "Unhandled matrix type");
192  break;
193  case eUPPER_TRIANGULAR_BANDED:
194  NEKERROR(ErrorUtil::efatal, "Unhandled matrix type");
195  break;
196  case eLOWER_TRIANGULAR_BANDED:
197  NEKERROR(ErrorUtil::efatal, "Unhandled matrix type");
198  break;
199 
200  default:
201  NEKERROR(ErrorUtil::efatal, "Unhandled matrix type");
202  }
203 
204  }
ASSERTL0
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:216
NEKERROR
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mod...
Definition: ErrorUtil.hpp:209
Lapack::Dsptrs
static void Dsptrs(const char &uplo, const int &n, const int &nrhs, const double *ap, const int *ipiv, double *b, const int &ldb, int &info)
Solve a real symmetric matrix problem using Bunch-Kaufman pivoting.
Definition: Lapack.hpp:115
Lapack::Dpbtrs
static void Dpbtrs(const char &uplo, const int &n, const int &kd, const int &nrhs, const double *ab, const int &ldab, double *b, const int &ldb, int &info)
Solve a real, Positive definite banded symmetric matrix problem using Cholesky factorization.
Definition: Lapack.hpp:157
Lapack::Dgetrs
static void Dgetrs(const char &trans, const int &n, const int &nrhs, const double *a, const int &lda, int *ipiv, double *b, const int &ldb, int &info)
General matrix LU backsolve.
Definition: Lapack.hpp:190
Lapack::Dtptrs
static void Dtptrs(const char &uplo, const char &trans, const char &diag, const int &n, const int &nrhs, const double *a, double *b, const int &ldb, int &info)
Definition: Lapack.hpp:245
Lapack::Dpptrs
static void Dpptrs(const char &uplo, const int &n, const int &nrhs, const double *ap, double *b, const int &ldb, int &info)
Solve a real Positive defiinte symmetric matrix problem using Cholesky factorization.
Definition: Lapack.hpp:139
Lapack::Dgbtrs
static void Dgbtrs(const char &trans, const int &n, const int &kl, const int &ku, const int &nrhs, const double *a, const int &lda, const int *ipiv, double *b, const int &ldb, int &info)
Solve general banded matrix using LU factorisation.
Definition: Lapack.hpp:174

◆ SolveTranspose()

template<typename BVectorType , typename XVectorType >
static void Nektar::LinearSystemSolver::SolveTranspose ( const BVectorType &  b,
XVectorType &  x,
MatrixStorage  m_matrixType,
const Array< OneD, const int > &  m_ipivot,
unsigned int  n,
const Array< OneD, const double > &  A,
char  m_transposeFlag,
unsigned int  m_numberOfSubDiagonals,
unsigned int  m_numberOfSuperDiagonals 
)
inlinestatic

Definition at line 207 of file NekLinSys.hpp.

References ASSERTL0, Lapack::Dgbtrs(), Lapack::Dgetrs(), Lapack::Dtptrs(), Nektar::eBANDED, Nektar::eDIAGONAL, Nektar::ErrorUtil::efatal, Nektar::eFULL, Nektar::eLOWER_TRIANGULAR, Nektar::eLOWER_TRIANGULAR_BANDED, Nektar::ePOSITIVE_DEFINITE_SYMMETRIC, Nektar::ePOSITIVE_DEFINITE_SYMMETRIC_BANDED, Nektar::eSYMMETRIC, Nektar::eSYMMETRIC_BANDED, Nektar::eUPPER_TRIANGULAR, Nektar::eUPPER_TRIANGULAR_BANDED, NEKERROR, and Xxt::Solve().

Referenced by Nektar::LinearSystem::SolveTranspose().

212  {
213  switch(m_matrixType)
214  {
215  case eFULL:
216  {
217  x = b;
218  int info = 0;
219  Lapack::Dgetrs('T',n,1,A.get(),n,(int *)m_ipivot.get(),x.GetRawPtr(), n,info);
220 
221  if( info < 0 )
222  {
223  std::string message = "ERROR: The " + std::to_string(-info) + "th parameter had an illegal parameter for dgetrs";
224  ASSERTL0(false, message.c_str());
225  }
226  }
227 
228  break;
229  case eDIAGONAL:
230  Solve(b, x, m_matrixType, m_ipivot, n, A, m_transposeFlag, m_numberOfSubDiagonals, m_numberOfSuperDiagonals);
231  break;
232  case eUPPER_TRIANGULAR:
233  {
234  char trans = m_transposeFlag;
235  if( trans == 'N' )
236  {
237  trans = 'T';
238  }
239  else
240  {
241  trans = 'N';
242  }
243 
244  x = b;
245  int info = 0;
246  Lapack::Dtptrs('U', trans, 'N', n, 1, A.get(), x.GetRawPtr(), n, info);
247 
248  if( info < 0 )
249  {
250  std::string message = "ERROR: The " + std::to_string(-info) + "th parameter had an illegal parameter for dtrtrs";
251  ASSERTL0(false, message.c_str());
252  }
253  else if( info > 0 )
254  {
255  std::string message = "ERROR: The " + std::to_string(-info) + "th diagonal element of A is 0 for dtrtrs";
256  ASSERTL0(false, message.c_str());
257  }
258  }
259 
260  break;
261  case eLOWER_TRIANGULAR:
262  {
263  char trans = m_transposeFlag;
264  if( trans == 'N' )
265  {
266  trans = 'T';
267  }
268  else
269  {
270  trans = 'N';
271  }
272  x = b;
273  int info = 0;
274  Lapack::Dtptrs('L', trans, 'N', n, 1, A.get(), x.GetRawPtr(), n, info);
275 
276  if( info < 0 )
277  {
278  std::string message = "ERROR: The " + std::to_string(-info) + "th parameter had an illegal parameter for dtrtrs";
279  ASSERTL0(false, message.c_str());
280  }
281  else if( info > 0 )
282  {
283  std::string message = "ERROR: The " + std::to_string(-info) + "th diagonal element of A is 0 for dtrtrs";
284  ASSERTL0(false, message.c_str());
285  }
286  }
287  break;
288  case eSYMMETRIC:
289  case ePOSITIVE_DEFINITE_SYMMETRIC:
290  case ePOSITIVE_DEFINITE_SYMMETRIC_BANDED:
291  Solve(b, x, m_matrixType, m_ipivot, n, A, m_transposeFlag, m_numberOfSubDiagonals, m_numberOfSuperDiagonals);
292  break;
293  case eBANDED:
294  {
295  x = b;
296  int KL = m_numberOfSubDiagonals;
297  int KU = m_numberOfSuperDiagonals;
298  int info = 0;
299 
300  Lapack::Dgbtrs(m_transposeFlag, n, KL, KU, 1, A.get(), 2*KL+KU+1, m_ipivot.get(), x.GetRawPtr(), n, info);
301 
302  if( info < 0 )
303  {
304  std::string message = "ERROR: The " + std::to_string(-info) + "th parameter had an illegal parameter for dgbtrs";
305  ASSERTL0(false, message.c_str());
306  }
307  }
308  break;
309  case eSYMMETRIC_BANDED:
310  NEKERROR(ErrorUtil::efatal, "Unhandled matrix type");
311  break;
312  case eUPPER_TRIANGULAR_BANDED:
313  NEKERROR(ErrorUtil::efatal, "Unhandled matrix type");
314  break;
315  case eLOWER_TRIANGULAR_BANDED:
316  NEKERROR(ErrorUtil::efatal, "Unhandled matrix type");
317  break;
318 
319  default:
320  NEKERROR(ErrorUtil::efatal, "Unhandled matrix type");
321  }
322  }
ASSERTL0
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:216
NEKERROR
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mod...
Definition: ErrorUtil.hpp:209
Lapack::Dgetrs
static void Dgetrs(const char &trans, const int &n, const int &nrhs, const double *a, const int &lda, int *ipiv, double *b, const int &ldb, int &info)
General matrix LU backsolve.
Definition: Lapack.hpp:190
Lapack::Dtptrs
static void Dtptrs(const char &uplo, const char &trans, const char &diag, const int &n, const int &nrhs, const double *a, double *b, const int &ldb, int &info)
Definition: Lapack.hpp:245
Nektar::LinearSystemSolver::Solve
static void Solve(const BVectorType &b, XVectorType &x, MatrixStorage m_matrixType, const Array< OneD, const int > &m_ipivot, unsigned int n, const Array< OneD, const double > &A, char m_transposeFlag, unsigned int m_numberOfSubDiagonals, unsigned int m_numberOfSuperDiagonals)
Definition: NekLinSys.hpp:72
Lapack::Dgbtrs
static void Dgbtrs(const char &trans, const int &n, const int &kl, const int &ku, const int &nrhs, const double *a, const int &lda, const int *ipiv, double *b, const int &ldb, int &info)
Solve general banded matrix using LU factorisation.
Definition: Lapack.hpp:174
Generated on Wed Jan 20 2021 03:33:53 for Nektar++ by   doxygen 1.8.13