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 72 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 75 of file NekLinSys.hpp.

81  {
82  switch (m_matrixType)
83  {
84  case eFULL:
85  {
86  x = b;
87  int info = 0;
88  Lapack::Dgetrs('N', n, 1, A.get(), n, (int *)m_ipivot.get(),
89  x.GetRawPtr(), n, info);
90  if (info < 0)
91  {
92  std::string message =
93  "ERROR: The " + std::to_string(-info) +
94  "th parameter had an illegal parameter for dgetrs";
95  ASSERTL0(false, message.c_str());
96  }
97  }
98  break;
99  case eDIAGONAL:
100  for (unsigned int i = 0; i < A.size(); ++i)
101  {
102  x[i] = b[i] * A[i];
103  }
104  break;
105  case eUPPER_TRIANGULAR:
106  {
107  x = b;
108  int info = 0;
109  Lapack::Dtptrs('U', m_transposeFlag, 'N', n, 1, A.get(),
110  x.GetRawPtr(), n, info);
111 
112  if (info < 0)
113  {
114  std::string message =
115  "ERROR: The " + std::to_string(-info) +
116  "th parameter had an illegal parameter for dtrtrs";
117  ASSERTL0(false, message.c_str());
118  }
119  else if (info > 0)
120  {
121  std::string message =
122  "ERROR: The " + std::to_string(-info) +
123  "th diagonal element of A is 0 for dtrtrs";
124  ASSERTL0(false, message.c_str());
125  }
126  }
127  break;
128  case eLOWER_TRIANGULAR:
129  {
130  x = b;
131  int info = 0;
132  Lapack::Dtptrs('L', m_transposeFlag, 'N', n, 1, A.get(),
133  x.GetRawPtr(), n, info);
134 
135  if (info < 0)
136  {
137  std::string message =
138  "ERROR: The " + std::to_string(-info) +
139  "th parameter had an illegal parameter for dtrtrs";
140  ASSERTL0(false, message.c_str());
141  }
142  else if (info > 0)
143  {
144  std::string message =
145  "ERROR: The " + std::to_string(-info) +
146  "th diagonal element of A is 0 for dtrtrs";
147  ASSERTL0(false, message.c_str());
148  }
149  }
150  break;
151  case eSYMMETRIC:
152  {
153  x = b;
154  int info = 0;
155  Lapack::Dsptrs('U', n, 1, A.get(), m_ipivot.get(),
156  x.GetRawPtr(), x.GetRows(), info);
157  if (info < 0)
158  {
159  std::string message =
160  "ERROR: The " + std::to_string(-info) +
161  "th parameter had an illegal parameter for dsptrs";
162  ASSERTL0(false, message.c_str());
163  }
164  }
165  break;
166  case ePOSITIVE_DEFINITE_SYMMETRIC:
167  {
168  x = b;
169  int info = 0;
170  Lapack::Dpptrs('U', n, 1, A.get(), x.GetRawPtr(), x.GetRows(),
171  info);
172  if (info < 0)
173  {
174  std::string message =
175  "ERROR: The " + std::to_string(-info) +
176  "th parameter had an illegal parameter for dpptrs";
177  ASSERTL0(false, message.c_str());
178  }
179  }
180  break;
181  case eBANDED:
182  {
183  x = b;
184  int KL = m_numberOfSubDiagonals;
185  int KU = m_numberOfSuperDiagonals;
186  int info = 0;
187 
188  Lapack::Dgbtrs(m_transposeFlag, n, KL, KU, 1, A.get(),
189  2 * KL + KU + 1, m_ipivot.get(), x.GetRawPtr(),
190  n, info);
191 
192  if (info < 0)
193  {
194  std::string message =
195  "ERROR: The " + std::to_string(-info) +
196  "th parameter had an illegal parameter for dgbtrs";
197  ASSERTL0(false, message.c_str());
198  }
199  }
200  break;
201  case ePOSITIVE_DEFINITE_SYMMETRIC_BANDED:
202  {
203  x = b;
204  int KU = m_numberOfSuperDiagonals;
205  int info = 0;
206 
207  Lapack::Dpbtrs('U', n, KU, 1, A.get(), KU + 1, x.GetRawPtr(), n,
208  info);
209 
210  if (info < 0)
211  {
212  std::string message =
213  "ERROR: The " + std::to_string(-info) +
214  "th parameter had an illegal parameter for dpbtrs";
215  ASSERTL0(false, message.c_str());
216  }
217  }
218  break;
219  case eSYMMETRIC_BANDED:
220  NEKERROR(ErrorUtil::efatal, "Unhandled matrix type");
221  break;
222  case eUPPER_TRIANGULAR_BANDED:
223  NEKERROR(ErrorUtil::efatal, "Unhandled matrix type");
224  break;
225  case eLOWER_TRIANGULAR_BANDED:
226  NEKERROR(ErrorUtil::efatal, "Unhandled matrix type");
227  break;
228 
229  default:
230  NEKERROR(ErrorUtil::efatal, "Unhandled matrix type");
231  }
232  }
ASSERTL0
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:215
NEKERROR
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
Definition: ErrorUtil.hpp:209
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:225
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:191
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:288
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:258
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:381
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:242
Nektar::eLOWER_TRIANGULAR_BANDED
@ eLOWER_TRIANGULAR_BANDED
Definition: MatrixStorageType.h:53
Nektar::ePOSITIVE_DEFINITE_SYMMETRIC_BANDED
@ ePOSITIVE_DEFINITE_SYMMETRIC_BANDED
Definition: MatrixStorageType.h:51
Nektar::ePOSITIVE_DEFINITE_SYMMETRIC
@ ePOSITIVE_DEFINITE_SYMMETRIC
Definition: MatrixStorageType.h:48
Nektar::eUPPER_TRIANGULAR_BANDED
@ eUPPER_TRIANGULAR_BANDED
Definition: MatrixStorageType.h:52

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(), and SolveTranspose().

◆ 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 235 of file NekLinSys.hpp.

243  {
244  switch (m_matrixType)
245  {
246  case eFULL:
247  {
248  x = b;
249  int info = 0;
250  Lapack::Dgetrs('T', n, 1, A.get(), n, (int *)m_ipivot.get(),
251  x.GetRawPtr(), n, info);
252 
253  if (info < 0)
254  {
255  std::string message =
256  "ERROR: The " + std::to_string(-info) +
257  "th parameter had an illegal parameter for dgetrs";
258  ASSERTL0(false, message.c_str());
259  }
260  }
261 
262  break;
263  case eDIAGONAL:
264  Solve(b, x, m_matrixType, m_ipivot, n, A, m_transposeFlag,
265  m_numberOfSubDiagonals, m_numberOfSuperDiagonals);
266  break;
267  case eUPPER_TRIANGULAR:
268  {
269  char trans = m_transposeFlag;
270  if (trans == 'N')
271  {
272  trans = 'T';
273  }
274  else
275  {
276  trans = 'N';
277  }
278 
279  x = b;
280  int info = 0;
281  Lapack::Dtptrs('U', trans, 'N', n, 1, A.get(), x.GetRawPtr(), n,
282  info);
283 
284  if (info < 0)
285  {
286  std::string message =
287  "ERROR: The " + std::to_string(-info) +
288  "th parameter had an illegal parameter for dtrtrs";
289  ASSERTL0(false, message.c_str());
290  }
291  else if (info > 0)
292  {
293  std::string message =
294  "ERROR: The " + std::to_string(-info) +
295  "th diagonal element of A is 0 for dtrtrs";
296  ASSERTL0(false, message.c_str());
297  }
298  }
299 
300  break;
301  case eLOWER_TRIANGULAR:
302  {
303  char trans = m_transposeFlag;
304  if (trans == 'N')
305  {
306  trans = 'T';
307  }
308  else
309  {
310  trans = 'N';
311  }
312  x = b;
313  int info = 0;
314  Lapack::Dtptrs('L', trans, 'N', n, 1, A.get(), x.GetRawPtr(), n,
315  info);
316 
317  if (info < 0)
318  {
319  std::string message =
320  "ERROR: The " + std::to_string(-info) +
321  "th parameter had an illegal parameter for dtrtrs";
322  ASSERTL0(false, message.c_str());
323  }
324  else if (info > 0)
325  {
326  std::string message =
327  "ERROR: The " + std::to_string(-info) +
328  "th diagonal element of A is 0 for dtrtrs";
329  ASSERTL0(false, message.c_str());
330  }
331  }
332  break;
333  case eSYMMETRIC:
334  case ePOSITIVE_DEFINITE_SYMMETRIC:
335  case ePOSITIVE_DEFINITE_SYMMETRIC_BANDED:
336  Solve(b, x, m_matrixType, m_ipivot, n, A, m_transposeFlag,
337  m_numberOfSubDiagonals, m_numberOfSuperDiagonals);
338  break;
339  case eBANDED:
340  {
341  x = b;
342  int KL = m_numberOfSubDiagonals;
343  int KU = m_numberOfSuperDiagonals;
344  int info = 0;
345 
346  Lapack::Dgbtrs(m_transposeFlag, n, KL, KU, 1, A.get(),
347  2 * KL + KU + 1, m_ipivot.get(), x.GetRawPtr(),
348  n, info);
349 
350  if (info < 0)
351  {
352  std::string message =
353  "ERROR: The " + std::to_string(-info) +
354  "th parameter had an illegal parameter for dgbtrs";
355  ASSERTL0(false, message.c_str());
356  }
357  }
358  break;
359  case eSYMMETRIC_BANDED:
360  NEKERROR(ErrorUtil::efatal, "Unhandled matrix type");
361  break;
362  case eUPPER_TRIANGULAR_BANDED:
363  NEKERROR(ErrorUtil::efatal, "Unhandled matrix type");
364  break;
365  case eLOWER_TRIANGULAR_BANDED:
366  NEKERROR(ErrorUtil::efatal, "Unhandled matrix type");
367  break;
368 
369  default:
370  NEKERROR(ErrorUtil::efatal, "Unhandled matrix type");
371  }
372  }
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:75

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 Solve().

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

Generated on Sat Mar 25 2023 12:35:55 for Nektar++ by doxygen 1.9.1