npstat is hosted by Hepforge, IPPP Durham
NPStat  5.10.0
npstat::AbsClassicalOrthoPoly1D Class Referenceabstract
Inheritance diagram for npstat::AbsClassicalOrthoPoly1D:
npstat::ChebyshevOrthoPoly1st npstat::ChebyshevOrthoPoly2nd npstat::ClassicalOrthoPoly1DFromWeight< Functor1D > npstat::DensityOrthoPoly1D npstat::HermiteProbOrthoPoly1D npstat::JacobiOrthoPoly1D npstat::JohnsonOrthoPoly1D npstat::LegendreOrthoPoly1D npstat::ShiftedLegendreOrthoPoly1D npstat::StorablePolySeries1D

Classes

class  MultiProdFcn
 
class  ProdFcn
 

Public Types

typedef OrthoPoly1DDeg degree_functor
 
typedef OrthoPoly1DWeight weight_functor
 

Public Member Functions

virtual AbsClassicalOrthoPoly1Dclone () const =0
 
virtual long double weight (long double x) const =0
 
virtual double xmin () const =0
 
virtual double xmax () const =0
 
virtual unsigned maxDegree () const
 
Matrix< long double > jacobiMatrix (unsigned n) const
 
long double poly (unsigned deg, long double x) const
 
long double monic (unsigned deg, long double x) const
 
std::pair< long double, long double > twopoly (unsigned deg1, unsigned deg2, long double x) const
 
void allpoly (long double x, long double *values, unsigned maxdeg) const
 
void allmonic (long double x, long double *values, unsigned maxdeg) const
 
double series (const double *coeffs, unsigned maxdeg, double x) const
 
double integrateSeries (const double *coeffs, unsigned maxdeg) const
 
template<class Functor , class Quadrature >
void calculateCoeffs (const Functor &fcn, const Quadrature &quad, double *coeffs, unsigned maxdeg) const
 
template<class Numeric >
void sampleCoeffs (const Numeric *coords, unsigned long lenCoords, double *coeffs, unsigned maxdeg) const
 
template<class Numeric >
void sampleCoeffVars (const Numeric *coords, unsigned long lenCoords, const double *coeffs, unsigned maxdeg, double *variances) const
 
template<class Numeric >
Matrix< double > sampleCoeffCovariance (const Numeric *coords, unsigned long lenCoords, const double *coeffs, unsigned maxdeg) const
 
template<class Numeric >
void sampleCoeffs (const Numeric *coords, unsigned long lenCoords, Numeric location, Numeric scale, double *coeffs, unsigned maxdeg) const
 
template<class Numeric >
void sampleCoeffVars (const Numeric *coords, unsigned long lenCoords, Numeric location, Numeric scale, const double *coeffs, unsigned maxdeg, double *variances) const
 
template<class Numeric >
Matrix< double > sampleCoeffCovariance (const Numeric *coords, unsigned long lenCoords, Numeric location, Numeric scale, const double *coeffs, unsigned maxdeg) const
 
template<class Pair >
void weightedSampleCoeffs (const Pair *points, unsigned long numPoints, double *coeffs, unsigned maxdeg) const
 
template<class Pair >
void weightedSampleCoeffVars (const Pair *points, unsigned long nPoints, const double *coeffs, unsigned maxdeg, double *variances) const
 
template<class Pair >
Matrix< double > weightedSampleCoeffCovariance (const Pair *points, unsigned long nPoints, const double *coeffs, unsigned maxdeg) const
 
template<class Pair >
void weightedSampleCoeffs (const Pair *points, unsigned long numPoints, typename Pair::first_type location, typename Pair::first_type scale, double *coeffs, unsigned maxdeg) const
 
template<class Pair >
void weightedSampleCoeffVars (const Pair *points, unsigned long nPoints, typename Pair::first_type location, typename Pair::first_type scale, const double *coeffs, unsigned maxdeg, double *variances) const
 
template<class Pair >
Matrix< double > weightedSampleCoeffCovariance (const Pair *points, unsigned long nPoints, typename Pair::first_type location, typename Pair::first_type scale, const double *coeffs, unsigned maxdeg) const
 
template<class Quadrature >
double empiricalKroneckerDelta (const Quadrature &quad, unsigned deg1, unsigned deg2) const
 
Matrix< double > empiricalKroneckerMatrix (const AbsIntervalQuadrature1D &quad, unsigned maxdeg) const
 
template<class Quadrature >
double jointAverage (const Quadrature &q, unsigned deg1, unsigned deg2, unsigned deg3, unsigned deg4) const
 
template<class Quadrature >
double jointAverage (const Quadrature &q, const unsigned *degrees, unsigned nDegrees, bool checkedForZeros=false) const
 
template<class Quadrature >
double directQuadrature (const Quadrature &q, unsigned deg1, unsigned deg2, unsigned deg3, unsigned deg4) const
 
template<class Quadrature >
double directQuadrature (const Quadrature &q, const unsigned *degrees, unsigned nDegrees, bool checkedForZeros=false) const
 
double integratePoly (unsigned degree, unsigned power) const
 
double jointIntegral (const unsigned *degrees, unsigned nDegrees) const
 
template<typename Numeric , typename Real >
void sampleAverages (const Numeric *coords, unsigned long lenCoords, Real *averages, unsigned maxdeg) const
 
template<class Pair , typename Real >
void weightedPointsAverages (const Pair *points, unsigned long numPoints, Real *averages, unsigned maxdeg) const
 
template<class Functor , typename Real >
void extWeightAverages (const Functor &extWeight, const AbsIntervalQuadrature1D &quad, Real *averages, unsigned maxdeg) const
 
template<class Numeric >
Matrix< double > sampleProductAverages (const Numeric *coords, unsigned long lenCoords, unsigned maxdeg) const
 
template<class Functor >
Matrix< double > extWeightProductAverages (const Functor &extWeight, const AbsIntervalQuadrature1D &quad, unsigned maxdeg) const
 
template<class Numeric >
double jointSampleAverage (const Numeric *coords, unsigned long lenCoords, const unsigned *degrees, unsigned lenDegrees) const
 
template<class Numeric >
std::pair< double, double > twoJointSampleAverages (const Numeric *coords, unsigned long lenCoords, const unsigned *degrees1, unsigned lenDegrees1, const unsigned *degrees2, unsigned lenDegrees2) const
 
CPP11_auto_ptr< StorablePolySeries1DmakeStorablePolySeries (unsigned maxPolyDeg, const double *coeffs=0, unsigned maxdeg=0, double shift=0.0) const
 
long double location () const
 
long double scale () const
 
virtual std::pair< long double, long double > monicRecurrenceCoeffs (unsigned deg) const =0
 
virtual std::pair< long double, long double > recurrenceCoeffs (unsigned deg) const =0
 

Static Protected Member Functions

static int kdelta (const unsigned i, const unsigned j)
 

Friends

class MultiProdFcn
 

Member Function Documentation

◆ allmonic()

void npstat::AbsClassicalOrthoPoly1D::allmonic ( long double  x,
long double *  values,
unsigned  maxdeg 
) const

Values of all monic orthonormal polynomials up to some degree. The size of the "values" array should be at least maxdeg + 1.

◆ allpoly()

void npstat::AbsClassicalOrthoPoly1D::allpoly ( long double  x,
long double *  values,
unsigned  maxdeg 
) const

Values of all orthonormal polynomials up to some degree. Faster than calling "poly" multiple times. The size of the "values" array should be at least maxdeg + 1.

◆ calculateCoeffs()

template<class Functor , class Quadrature >
void npstat::AbsClassicalOrthoPoly1D::calculateCoeffs ( const Functor &  fcn,
const Quadrature &  quad,
double *  coeffs,
unsigned  maxdeg 
) const

Build the coefficients of the orthogonal polynomial series for the given function. The length of the array "coeffs" should be at least maxdeg + 1. Note that the coefficients are returned in the order of increasing degree (same order is used by the "series" function).

◆ clone()

◆ directQuadrature() [1/2]

template<class Quadrature >
double npstat::AbsClassicalOrthoPoly1D::directQuadrature ( const Quadrature &  q,
const unsigned *  degrees,
unsigned  nDegrees,
bool  checkedForZeros = false 
) const

Direct unweighted integration of of a product of multiple orthonormal poly values. The integration function should support the "integrate" method without limits. Intended for use with GaussHermiteQuadrature and such. "checkedForZeros" argument should be set to "true" if we are sure that there are no zero elements in the "degrees" array.

◆ directQuadrature() [2/2]

template<class Quadrature >
double npstat::AbsClassicalOrthoPoly1D::directQuadrature ( const Quadrature &  q,
unsigned  deg1,
unsigned  deg2,
unsigned  deg3,
unsigned  deg4 
) const

Direct unweighted integration of of a product of four orthonormal poly values. The integration function should support the "integrate" method without limits. Intended for use with GaussHermiteQuadrature and such.

◆ empiricalKroneckerDelta()

template<class Quadrature >
double npstat::AbsClassicalOrthoPoly1D::empiricalKroneckerDelta ( const Quadrature &  quad,
unsigned  deg1,
unsigned  deg2 
) const

This method is useful for testing the numerical precision of the orthonormalization procedure. It returns the scalar products between various polynomials.

◆ empiricalKroneckerMatrix()

Matrix<double> npstat::AbsClassicalOrthoPoly1D::empiricalKroneckerMatrix ( const AbsIntervalQuadrature1D quad,
unsigned  maxdeg 
) const

A faster way to generate Kronecker deltas if the whole matrix of them is needed. The "maxdeg" argument must not exceed the return value of the "maxDegree" function. The resulting matrix will be dimensioned (maxdeg + 1) x (maxdeg + 1). Note that this function will produce meaningful results only for polynomials supported on compact intervals.

◆ extWeightAverages()

template<class Functor , typename Real >
void npstat::AbsClassicalOrthoPoly1D::extWeightAverages ( const Functor &  extWeight,
const AbsIntervalQuadrature1D quad,
Real *  averages,
unsigned  maxdeg 
) const

Averages of the polynomial values weighted by an external weight function. The length of array "averages" should be at least maxdeg + 1.

◆ extWeightProductAverages()

template<class Functor >
Matrix<double> npstat::AbsClassicalOrthoPoly1D::extWeightProductAverages ( const Functor &  extWeight,
const AbsIntervalQuadrature1D quad,
unsigned  maxdeg 
) const

Pairwise products of the polynomial values weighted by an external weight function. The returned matrix will be symmetric and will have the dimensions (maxdeg + 1) x (maxdeg + 1).

◆ integratePoly()

double npstat::AbsClassicalOrthoPoly1D::integratePoly ( unsigned  degree,
unsigned  power 
) const

An unweighted integral of the orthonormal polynomial with the given degree to some power. For this method, it is assumed that the polynomials are supported on a closed interval (without such an assumption unweighted integrals do not make much sense) and that Gauss-Legendre quadratures can be used.

◆ integrateSeries()

double npstat::AbsClassicalOrthoPoly1D::integrateSeries ( const double *  coeffs,
unsigned  maxdeg 
) const

Integral of the series from xmin() to xmax()

◆ jacobiMatrix()

Matrix<long double> npstat::AbsClassicalOrthoPoly1D::jacobiMatrix ( unsigned  n) const

Generate principal minor of order n of the Jacobi matrix. n must not exceed the output of "maxDegree()" method.

◆ jointAverage() [1/2]

template<class Quadrature >
double npstat::AbsClassicalOrthoPoly1D::jointAverage ( const Quadrature &  q,
const unsigned *  degrees,
unsigned  nDegrees,
bool  checkedForZeros = false 
) const

A measure-weighted average of a product of multiple orthonormal poly values. "checkedForZeros" argument should be set to "true" if we are sure that there are no zero elements in the "degrees" array.

◆ jointAverage() [2/2]

template<class Quadrature >
double npstat::AbsClassicalOrthoPoly1D::jointAverage ( const Quadrature &  q,
unsigned  deg1,
unsigned  deg2,
unsigned  deg3,
unsigned  deg4 
) const

A measure-weighted average of a product of four orthonormal poly values

◆ jointIntegral()

double npstat::AbsClassicalOrthoPoly1D::jointIntegral ( const unsigned *  degrees,
unsigned  nDegrees 
) const

An unweighted integral of a product of multiple orthonormal polynomials. For this method, it is assumed that the polynomials are supported on a closed interval (without such an assumption unweighted integrals do not make much sense) and that Gauss-Legendre quadratures can be used.

◆ jointSampleAverage()

template<class Numeric >
double npstat::AbsClassicalOrthoPoly1D::jointSampleAverage ( const Numeric *  coords,
unsigned long  lenCoords,
const unsigned *  degrees,
unsigned  lenDegrees 
) const

Unweighted average of a product of polynomial values for the given sample. "degrees" is the collection of polynomial degrees. Polynomials of these degrees will be included in the product.

◆ maxDegree()

virtual unsigned npstat::AbsClassicalOrthoPoly1D::maxDegree ( ) const
inlinevirtual

◆ monic()

long double npstat::AbsClassicalOrthoPoly1D::monic ( unsigned  deg,
long double  x 
) const

Values of the corresponding monic polynomials

◆ poly()

long double npstat::AbsClassicalOrthoPoly1D::poly ( unsigned  deg,
long double  x 
) const

Polynomial values

◆ sampleAverages()

template<typename Numeric , typename Real >
void npstat::AbsClassicalOrthoPoly1D::sampleAverages ( const Numeric *  coords,
unsigned long  lenCoords,
Real *  averages,
unsigned  maxdeg 
) const

Unweighted averages of the polynomial values for the given sample. The length of array "averages" should be at least maxdeg + 1.

◆ sampleCoeffCovariance() [1/2]

template<class Numeric >
Matrix<double> npstat::AbsClassicalOrthoPoly1D::sampleCoeffCovariance ( const Numeric *  coords,
unsigned long  lenCoords,
const double *  coeffs,
unsigned  maxdeg 
) const

Estimate the covariances of the coefficients returned by the "sampleCoeffs" function. The "coeffs" array should be at least maxdeg + 1 elements long and should be filled by a previous call to "sampleCoeffs". The returned matrix will be dimensioned (maxdeg + 1) x (maxdeg + 1).

◆ sampleCoeffCovariance() [2/2]

template<class Numeric >
Matrix<double> npstat::AbsClassicalOrthoPoly1D::sampleCoeffCovariance ( const Numeric *  coords,
unsigned long  lenCoords,
Numeric  location,
Numeric  scale,
const double *  coeffs,
unsigned  maxdeg 
) const

Estimate the covariances of the coefficients returned by the "sampleCoeffs" function. The "coeffs" array should be at least maxdeg + 1 elements long and should be filled by a previous call to "sampleCoeffs" with the same location and scale. The returned matrix will be dimensioned (maxdeg + 1) x (maxdeg + 1).

◆ sampleCoeffs() [1/2]

template<class Numeric >
void npstat::AbsClassicalOrthoPoly1D::sampleCoeffs ( const Numeric *  coords,
unsigned long  lenCoords,
double *  coeffs,
unsigned  maxdeg 
) const

Build the coefficients of the orthogonal polynomial series for the given sample of points (empirical density function). The length of the array "coeffs" should be at least maxdeg + 1. Note that the coefficients are returned in the order of increasing degree (same order is used by the "series" function).

◆ sampleCoeffs() [2/2]

template<class Numeric >
void npstat::AbsClassicalOrthoPoly1D::sampleCoeffs ( const Numeric *  coords,
unsigned long  lenCoords,
Numeric  location,
Numeric  scale,
double *  coeffs,
unsigned  maxdeg 
) const

Build the coefficients of the orthogonal polynomial series for the given sample of points (empirical density function). The length of the array "coeffs" should be at least maxdeg + 1. Note that the coefficients are returned in the order of increasing degree (same order is used by the "series" function). Before calculating the coefficients, the coordinates are shifted and scaled according to x_new = (x_original - location)/scale. The resulting coefficients are also divided by scale.

◆ sampleCoeffVars() [1/2]

template<class Numeric >
void npstat::AbsClassicalOrthoPoly1D::sampleCoeffVars ( const Numeric *  coords,
unsigned long  lenCoords,
const double *  coeffs,
unsigned  maxdeg,
double *  variances 
) const

Estimate the variances of the coefficients returned by the "sampleCoeffs" function. The "coeffs" array should be at least maxdeg + 1 elements long and should be filled by a previous call to "sampleCoeffs". The "variances" array should have at least maxdeg + 1 elements. It will contain the respective variances upon return.

◆ sampleCoeffVars() [2/2]

template<class Numeric >
void npstat::AbsClassicalOrthoPoly1D::sampleCoeffVars ( const Numeric *  coords,
unsigned long  lenCoords,
Numeric  location,
Numeric  scale,
const double *  coeffs,
unsigned  maxdeg,
double *  variances 
) const

Estimate the variances of the coefficients returned by the "sampleCoeffs" function. The "coeffs" array should be at least maxdeg + 1 elements long and should be filled by a previous call to "sampleCoeffs" with the same location and scale. The "variances" array should have at least maxdeg + 1 elements. It will contain the respective variances upon return.

◆ sampleProductAverages()

template<class Numeric >
Matrix<double> npstat::AbsClassicalOrthoPoly1D::sampleProductAverages ( const Numeric *  coords,
unsigned long  lenCoords,
unsigned  maxdeg 
) const

Unweighted averages of the pairwise products of the polynomial values for the given sample. The returned matrix will be symmetric and will have the dimensions (maxdeg + 1) x (maxdeg + 1).

◆ series()

double npstat::AbsClassicalOrthoPoly1D::series ( const double *  coeffs,
unsigned  maxdeg,
double  x 
) const

Polynomial series

◆ twopoly()

std::pair<long double,long double> npstat::AbsClassicalOrthoPoly1D::twopoly ( unsigned  deg1,
unsigned  deg2,
long double  x 
) const

Values of two orthonormal polynomials. Faster than calling "poly" two times.

◆ weight()

◆ weightedPointsAverages()

template<class Pair , typename Real >
void npstat::AbsClassicalOrthoPoly1D::weightedPointsAverages ( const Pair *  points,
unsigned long  numPoints,
Real *  averages,
unsigned  maxdeg 
) const

Averages of the polynomial values for the given sample of weighted points (not weighting by the polynomial weight function). The first element of the pair will be treated as the coordinate and the second element will be treated as weight. Weights must not be negative. The length of array "averages" should be at least maxdeg + 1.

◆ weightedSampleCoeffCovariance() [1/2]

template<class Pair >
Matrix<double> npstat::AbsClassicalOrthoPoly1D::weightedSampleCoeffCovariance ( const Pair *  points,
unsigned long  nPoints,
const double *  coeffs,
unsigned  maxdeg 
) const

Estimate the covariances of the coefficients returned by the "weightedSampleCoeffs" function. The "coeffs" array should be at least maxdeg + 1 elements long and should be filled by a previous call to "weightedSampleCoeffs". The returned matrix will be dimensioned (maxdeg + 1) x (maxdeg + 1). This code assumes that weights and coordinates are statistically independent from each other.

◆ weightedSampleCoeffCovariance() [2/2]

template<class Pair >
Matrix<double> npstat::AbsClassicalOrthoPoly1D::weightedSampleCoeffCovariance ( const Pair *  points,
unsigned long  nPoints,
typename Pair::first_type  location,
typename Pair::first_type  scale,
const double *  coeffs,
unsigned  maxdeg 
) const

Estimate the covariances of the coefficients returned by the "weightedSampleCoeffs" function. The "coeffs" array should be at least maxdeg + 1 elements long and should be filled by a previous call to "weightedSampleCoeffs" with the same location and scale. The returned matrix will be dimensioned (maxdeg + 1) x (maxdeg + 1). This code assumes that weights and coordinates are statistically independent from each other.

◆ weightedSampleCoeffs() [1/2]

template<class Pair >
void npstat::AbsClassicalOrthoPoly1D::weightedSampleCoeffs ( const Pair *  points,
unsigned long  numPoints,
double *  coeffs,
unsigned  maxdeg 
) const

Build the coefficients of the orthogonal polynomial series for the given sample of weighted points (empirical density function). The first element of the pair will be treated as the coordinate and the second element will be treated as weight. Weights must not be negative. The length of the array "coeffs" should be at least maxdeg + 1. Note that the coefficients are returned in the order of increasing degree (same order is used by the "series" function).

◆ weightedSampleCoeffs() [2/2]

template<class Pair >
void npstat::AbsClassicalOrthoPoly1D::weightedSampleCoeffs ( const Pair *  points,
unsigned long  numPoints,
typename Pair::first_type  location,
typename Pair::first_type  scale,
double *  coeffs,
unsigned  maxdeg 
) const

Build the coefficients of the orthogonal polynomial series for the given sample of weighted points (empirical density function). The first element of the pair will be treated as the coordinate and the second element will be treated as weight. Weights must not be negative. The length of the array "coeffs" should be at least maxdeg + 1. Note that the coefficients are returned in the order of increasing degree (same order is used by the "series" function). Before calculating the coefficients, the coordinates are shifted and scaled according to x_new = (x_original - location)/scale. The resulting coefficients are also divided by scale.

◆ weightedSampleCoeffVars() [1/2]

template<class Pair >
void npstat::AbsClassicalOrthoPoly1D::weightedSampleCoeffVars ( const Pair *  points,
unsigned long  nPoints,
const double *  coeffs,
unsigned  maxdeg,
double *  variances 
) const

Estimate the variances of the coefficients returned by the "weightedSampleCoeffs" function. The "coeffs" array should be at least maxdeg + 1 elements long and should be filled by a previous call to "weightedSampleCoeffs". The "variances" array should have at least maxdeg + 1 elements. It will contain the respective variances upon return. This code assumes that weights and coordinates are statistically independent from each other.

◆ weightedSampleCoeffVars() [2/2]

template<class Pair >
void npstat::AbsClassicalOrthoPoly1D::weightedSampleCoeffVars ( const Pair *  points,
unsigned long  nPoints,
typename Pair::first_type  location,
typename Pair::first_type  scale,
const double *  coeffs,
unsigned  maxdeg,
double *  variances 
) const

Estimate the variances of the coefficients returned by the "weightedSampleCoeffs" function. The "coeffs" array should be at least maxdeg + 1 elements long and should be filled by a previous call to "weightedSampleCoeffs" with the same location and scale. The "variances" array should have at least maxdeg + 1 elements. It will contain the respective variances upon return. This code assumes that weights and coordinates are statistically independent from each other.

◆ xmin()


The documentation for this class was generated from the following file: