ROL
NodalBasis.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Rapid Optimization Library (ROL) Package
4 //
5 // Copyright 2014 NTESS and the ROL contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 #include<iostream>
11 #include"ROL_StdVector.hpp"
13 #include"Lagrange.hpp"
14 #include"LinearAlgebra.hpp"
15 
16 #ifndef __NODAL_BASIS__
17 #define __NODAL_BASIS__
18 
19 
20 template<class Real>
21 struct NodalBasis{
22 
24  ROL::Ptr<Teuchos::LAPACK<int,Real> > lapack_;
25 
27  const int ni_;
28 
30  const int nq_;
31 
32  NodalBasis(ROL::Ptr<Teuchos::LAPACK<int,Real> > lapack, const int ni, const int nq);
33  ~NodalBasis();
34 
36  std::vector<Real> xi_;
37 
39  std::vector<Real> xq_;
40 
42  std::vector<Real> wq_;
43 
45  std::vector<Real> L_;
46  std::vector<Real> Lp_;
47 
49  ROL::Ptr<Lagrange<Real> > lagrange_;
50 
51 
52 };
53 
56 template<class Real>
57 NodalBasis<Real>::NodalBasis(ROL::Ptr<Teuchos::LAPACK<int,Real> > const lapack, const int ni, const int nq):
58  lapack_(lapack), ni_(ni), nq_(nq), xi_(ni_,0), xq_(nq_,0), wq_(nq_,0), L_(ni_*nq_,0), Lp_(ni_*nq_,0)
59  {
60 
61  // Generate Legendre-Gauss-Lobatto nodes and weights (unused)
62  std::vector<Real> ai(ni_,0);
63  std::vector<Real> bi(ni_,0);
64  std::vector<Real> wi(ni_,0);
65  rec_jacobi(this->lapack_,0,0,ai,bi);
66  rec_lobatto(this->lapack_,-1.0,1.0,ai,bi);
67  gauss(this->lapack_,ai,bi,xi_,wi);
68 
69  // Generate Legendre-Gauss nodes and weights
70  std::vector<Real> aq(nq_,0);
71  std::vector<Real> bq(nq_,0);
72  rec_jacobi(this->lapack_,0,0,aq,bq);
73  gauss(this->lapack_,aq,bq,xq_,wq_);
74 
75  std::vector<Real> e(ni_,0);
76  std::vector<Real> ell(nq,0);
77 
78  lagrange_ = ROL::makePtr<Lagrange<Real>>(xi_,xq_);
79 
80  // Loop over canonical vectors
81  for(int i=0;i<ni_;++i) {
82 
83  lagrange_->interpolant(i,ell);
84 
85  std::copy(ell.begin(),ell.end(),L_.begin()+i*nq_);
86 
87 
88  lagrange_->derivative(i,ell);
89  std::copy(ell.begin(),ell.end(),Lp_.begin()+i*nq_);
90 
91  // Rezero the vector
92  std::fill(e.begin(),e.end(),0);
93  }
94 
95 }
96 
97 template<class Real>
99 }
100 
101 #endif
std::vector< Real > Lp_
Definition: NodalBasis.hpp:46
void rec_lobatto(ROL::Ptr< Teuchos::LAPACK< int, Real > > const lapack, const double xl1, const double xl2, std::vector< Real > &a, std::vector< Real > &b)
Modify the given recurrence coefficients so that the set of zeros of the maximal order polynomial inc...
ROL::Ptr< Teuchos::LAPACK< int, Real > > lapack_
Definition: NodalBasis.hpp:24
NodalBasis(ROL::Ptr< Teuchos::LAPACK< int, Real > > lapack, const int ni, const int nq)
Set up quantities we will need repeatedly.
Definition: NodalBasis.hpp:57
const int ni_
Definition: NodalBasis.hpp:27
const int nq_
Definition: NodalBasis.hpp:30
std::vector< Real > xi_
Definition: NodalBasis.hpp:36
ROL::Ptr< Lagrange< Real > > lagrange_
Object for working with Lagrange polynomials and their derivatives.
Definition: NodalBasis.hpp:49
void gauss(ROL::Ptr< Teuchos::LAPACK< int, Real > > lapack, const std::vector< Real > &a, const std::vector< Real > &b, std::vector< Real > &x, std::vector< Real > &w)
Compute the Gauss quadrature nodes and weights for the polynomials generated by the recurrence coeffi...
std::vector< Real > xq_
Definition: NodalBasis.hpp:39
std::vector< Real > L_
Definition: NodalBasis.hpp:45
void rec_jacobi(ROL::Ptr< Teuchos::LAPACK< int, Real > > lapack, const double alpha, const double beta, std::vector< Real > &a, std::vector< Real > &b)
Generate the Jacobi polynomial recursion coeffcients .
std::vector< Real > wq_
Definition: NodalBasis.hpp:42