10 #include"Teuchos_LAPACK.hpp"
13 #ifndef __LINEAR_ALGEBRA__
14 #define __LINEAR_ALGEBRA__
24 void trisolve(ROL::Ptr<Teuchos::LAPACK<int,Real> > lapack,
25 const std::vector<Real>& a,
26 const std::vector<Real>& b,
27 const std::vector<Real>& c,
28 const std::vector<Real>& r,
29 std::vector<Real>& x ) {
31 const char TRANS =
'N';
32 const int N = b.size();
36 std::vector<Real> dl(a);
37 std::vector<Real> d(b);
38 std::vector<Real> du(c);
39 std::vector<Real> du2(N-2,0.0);
40 std::vector<int> ipiv(N);
43 lapack->GTTRF(N,&dl[0],&d[0],&du[0],&du2[0],&ipiv[0],&info);
48 lapack->GTTRS(TRANS,N,NRHS,&dl[0],&d[0],&du[0],&du2[0],&ipiv[0],&x[0],N,&info);
59 void lusolve(ROL::Ptr<Teuchos::LAPACK<int,Real> > lapack,
60 const std::vector<Real> &A,
61 const std::vector<Real> &B,
62 std::vector<Real> &X) {
64 const char TRANS =
'N';
65 const int N2 = A.size();
66 const int N = sqrt(N2);
67 const int M = int(B.size()/N);
70 std::vector<int> ipiv(N);
71 std::vector<Real> PLU(A);
75 lapack->GETRF(N,N,&PLU[0],LDA,&ipiv[0],&info);
78 lapack->GETRS(TRANS,N,M,&PLU[0],LDA,&ipiv[0],&X[0],LDA,&info);
90 void cholsolve(ROL::Ptr<Teuchos::LAPACK<int,Real> > lapack,
91 const std::vector<Real> &A,
92 const std::vector<Real> &B,
93 std::vector<Real> &X) {
95 const char UPLO =
'U';
96 const int N2 = A.size();
97 const int N = sqrt(N2);
98 const int M = int(B.size()/N);
101 std::vector<Real> C(A);
105 lapack->POTRF(UPLO,N,&C[0],LDA,&info);
108 lapack->POTRS(UPLO,N,&C[0],LDA,&X[0],LDA,&info);
void trisolve(ROL::Ptr< Teuchos::LAPACK< int, Real > > lapack, const std::vector< Real > &a, const std::vector< Real > &b, const std::vector< Real > &c, const std::vector< Real > &r, std::vector< Real > &x)
Solve a tridiagonal system.
void cholsolve(ROL::Ptr< Teuchos::LAPACK< int, Real > > lapack, const std::vector< Real > &A, const std::vector< Real > &B, std::vector< Real > &X)
void lusolve(ROL::Ptr< Teuchos::LAPACK< int, Real > > lapack, const std::vector< Real > &A, const std::vector< Real > &B, std::vector< Real > &X)