42template <
typename ordinal_type,
typename value_type>
45 const ordinal_type& p,
46 value_type (*weightFn)(
const value_type&),
47 const value_type& leftEndPt,
48 const value_type& rightEndPt,
51 RecurrenceBasis<ordinal_type,value_type>(
std::string(
"DiscretizedStieltjes -- ") + label, p, normalize, growth),
53 leftEndPt_(leftEndPt),
54 rightEndPt_(rightEndPt),
58 Teuchos::RCP<const Stokhos::LegendreBasis<ordinal_type,value_type> > quadBasis =
66template <
typename ordinal_type,
typename value_type>
71 scaleFactor(basis.scaleFactor),
72 leftEndPt_(basis.leftEndPt_),
73 rightEndPt_(basis.rightEndPt_),
74 weightFn_(basis.weightFn_)
77 Teuchos::RCP<const Stokhos::LegendreBasis<ordinal_type,value_type> > quadBasis =
89template <
typename ordinal_type,
typename value_type>
95template <
typename ordinal_type,
typename value_type>
99 Teuchos::Array<value_type>& alpha,
100 Teuchos::Array<value_type>& beta,
101 Teuchos::Array<value_type>& delta,
102 Teuchos::Array<value_type>& gamma)
const
114 value_type oneNorm = expectedValue_J_nsquared(0, alpha, beta);
116 scaleFactor = 1/oneNorm;
119 value_type integral2;
124 value_type past_integral = expectedValue_J_nsquared(0, alpha, beta);
125 alpha[0] = expectedValue_tJ_nsquared(0, alpha, beta)/past_integral;
131 for (ordinal_type n = 1; n<k; n++){
132 integral2 = expectedValue_J_nsquared(n, alpha, beta);
133 alpha[n] = expectedValue_tJ_nsquared(n, alpha, beta)/integral2;
134 beta[n] = integral2/past_integral;
135 past_integral = integral2;
143template <
typename ordinal_type,
typename value_type>
148 return (x < leftEndPt_ || x > rightEndPt_) ? 0: scaleFactor*weightFn_(x);
151template <
typename ordinal_type,
typename value_type>
155 const Teuchos::Array<value_type>& alpha,
156 const Teuchos::Array<value_type>& beta)
const
160 value_type integral = 0;
161 for(ordinal_type quadIdx = 0;
162 quadIdx < static_cast<ordinal_type>(quad_points.size()); quadIdx++) {
163 value_type x = (rightEndPt_ - leftEndPt_)*.5*quad_points[quadIdx] +
164 (rightEndPt_ + leftEndPt_)*.5;
165 value_type
val = evaluateRecurrence(x,order,alpha,beta);
166 integral += x*
val*
val*evaluateWeight(x)*quad_weights[quadIdx];
169 return integral*(rightEndPt_ - leftEndPt_);
172template <
typename ordinal_type,
typename value_type>
176 const Teuchos::Array<value_type>& alpha,
177 const Teuchos::Array<value_type>& beta)
const
181 value_type integral = 0;
182 for(ordinal_type quadIdx = 0;
183 quadIdx < static_cast<ordinal_type>(quad_points.size()); quadIdx++){
184 value_type x = (rightEndPt_ - leftEndPt_)*.5*quad_points[quadIdx] +
185 (rightEndPt_ + leftEndPt_)*.5;
186 value_type
val = evaluateRecurrence(x,order,alpha,beta);
187 integral +=
val*
val*evaluateWeight(x)*quad_weights[quadIdx];
190 return integral*(rightEndPt_ - leftEndPt_);
193template <
typename ordinal_type,
typename value_type>
201 value_type integral = 0;
202 for(ordinal_type quadIdx = 0;
203 quadIdx < static_cast<ordinal_type>(quad_points.size()); quadIdx++){
204 value_type x = (rightEndPt_ - leftEndPt_)*.5*quad_points[quadIdx] +
205 (rightEndPt_ + leftEndPt_)*.5;
206 integral += this->evaluate(x,order1)*this->evaluate(x,order2)*evaluateWeight(x)*quad_weights[quadIdx];
209 return integral*(rightEndPt_ - leftEndPt_);
212template <
typename ordinal_type,
typename value_type>
217 const Teuchos::Array<value_type>& alpha,
218 const Teuchos::Array<value_type>& beta)
const
221 return value_type(1.0);
225 value_type v0 = value_type(1.0);
226 value_type v1 = x-alpha[0]*v0;
227 value_type v2 = value_type(0.0);
228 for (ordinal_type i=2; i<=k; i++) {
229 v2 = (x-alpha[i-1])*v1 - beta[i-1]*v0;
237template <
typename ordinal_type,
typename value_type>
238Teuchos::RCP<Stokhos::OneDOrthogPolyBasis<ordinal_type,value_type> >
Generates three-term recurrence using the Discretized Stieltjes procedure.
value_type expectedValue_tJ_nsquared(const ordinal_type &order, const Teuchos::Array< value_type > &alpha, const Teuchos::Array< value_type > &beta) const
Approximates where = order.
Teuchos::Array< value_type > quad_weights
Quadrature points for discretized stieltjes procedure.
~DiscretizedStieltjesBasis()
Destructor.
virtual bool computeRecurrenceCoefficients(ordinal_type n, Teuchos::Array< value_type > &alpha, Teuchos::Array< value_type > &beta, Teuchos::Array< value_type > &delta, Teuchos::Array< value_type > &gamma) const
Compute recurrence coefficients.
Teuchos::Array< Teuchos::Array< value_type > > quad_values
Quadrature values for discretized stieltjes procedure.
Teuchos::Array< value_type > quad_points
Quadrature points for discretized stieltjes procedure.
DiscretizedStieltjesBasis(const std::string &name, const ordinal_type &p, value_type(*weightFn)(const value_type &), const value_type &leftEndPt, const value_type &rightEndPt, bool normalize, GrowthPolicy growth=SLOW_GROWTH)
Constructor.
value_type eval_inner_product(const ordinal_type &order1, const ordinal_type &order2) const
Evaluate inner product of two basis functions to test orthogonality.
value_type evaluateWeight(const value_type &x) const
Evaluates the scaled weight function.
value_type evaluateRecurrence(const value_type &point, ordinal_type order, const Teuchos::Array< value_type > &alpha, const Teuchos::Array< value_type > &beta) const
Evaluate 3-term recurrence formula using supplied coefficients.
value_type expectedValue_J_nsquared(const ordinal_type &order, const Teuchos::Array< value_type > &alpha, const Teuchos::Array< value_type > &beta) const
Approximates where = order.
virtual Teuchos::RCP< OneDOrthogPolyBasis< ordinal_type, value_type > > cloneWithOrder(ordinal_type p) const
Clone this object with the option of building a higher order basis.
Legendre polynomial basis.
Implementation of OneDOrthogPolyBasis based on the general three-term recurrence relationship:
Teuchos::Array< value_type > alpha
Recurrence coefficients.
Teuchos::Array< value_type > beta
Recurrence coefficients.
ordinal_type p
Order of basis.
Teuchos::Array< value_type > gamma
Recurrence coefficients.
Teuchos::Array< value_type > delta
Recurrence coefficients.
virtual void setup()
Setup basis after computing recurrence coefficients.
GrowthPolicy
Enumerated type for determining Smolyak growth policies.