Page MenuHomeHEPForge
Files F10881448

No OneTemporary
Actions

View Options

Index: trunk/tests/test_ContOrthoPoly1D.cc
===================================================================
--- trunk/tests/test_ContOrthoPoly1D.cc (revision 562)
+++ trunk/tests/test_ContOrthoPoly1D.cc (revision 563)
@@ -1,138 +1,179 @@
#include "UnitTest++.h"
#include "test_utils.hh"
#include "npstat/nm/ContOrthoPoly1D.hh"
#include "npstat/rng/MersenneTwister.hh"
using namespace npstat;
using namespace std;
inline static int kdelta(const unsigned i, const unsigned j)
{
return i == j ? 1 : 0;
}
namespace {
TEST(ContOrthoPoly1D_orthonormalization)
{
const double eps = 1.0e-15;
const OrthoPolyMethod method[] = {OPOLY_STIELTJES, OPOLY_LANCZOS};
const unsigned nMethods = sizeof(method)/sizeof(method[0]);
const unsigned npoints = 64;
const unsigned maxdeg = 10;
const unsigned ntries = 10;
std::vector<double> points(2U*npoints);
std::vector<ContOrthoPoly1D::MeasurePoint> measure;
measure.reserve(npoints);
std::vector<double> coeffs(maxdeg + 1);
MersenneTwister rng;
for (unsigned imeth=0; imeth<nMethods; ++imeth)
for (unsigned itry=0; itry<ntries; ++itry)
{
const double xmin = rng()*5.0 - 6.0;
const double xmax = rng()*2.0 + 1.0;
const double range = xmax - xmin;
rng.run(&points[0], 2U*npoints, 2U*npoints);
measure.clear();
for (unsigned i=0; i<npoints; ++i)
{
ContOrthoPoly1D::MeasurePoint p(points[2*i]*range + xmin, points[2*i+1]);
measure.push_back(p);
}
ContOrthoPoly1D poly(maxdeg, measure, method[imeth]);
for (unsigned i=0; i<=maxdeg; ++i)
for (unsigned j=0; j<=maxdeg; ++j)
{
const double d = poly.empiricalKroneckerDelta(i, j);
CHECK_CLOSE(static_cast<double>(kdelta(i, j)), d, eps);
}
measure.clear();
for (unsigned i=0; i<npoints; ++i)
{
ContOrthoPoly1D::MeasurePoint p(points[2*i]*range + xmin, 1.0);
measure.push_back(p);
}
ContOrthoPoly1D poly2(maxdeg, measure, method[imeth]);
poly2.weightCoeffs(&coeffs[0], maxdeg);
for (unsigned i=0; i<=maxdeg; ++i)
CHECK_CLOSE(static_cast<double>(kdelta(i, 0)), coeffs[i], eps);
}
}
+ TEST(ContOrthoPoly1D_weightCoeffs)
+ {
+ const double eps = 1.0e-8;
+
+ const unsigned maxdeg = 10;
+ const unsigned ntries = 10;
+ const unsigned npoints = maxdeg + 1U;
+
+ std::vector<double> points(2U*npoints);
+ std::vector<ContOrthoPoly1D::MeasurePoint> measure;
+ measure.reserve(npoints);
+ std::vector<double> coeffs(maxdeg + 1);
+
+ MersenneTwister rng;
+
+ for (unsigned itry=0; itry<ntries; ++itry)
+ {
+ const double xmin = rng()*5.0 - 6.0;
+ const double xmax = rng()*2.0 + 1.0;
+ const double range = xmax - xmin;
+ rng.run(&points[0], 2U*npoints, 2U*npoints);
+
+ measure.clear();
+ for (unsigned i=0; i<npoints; ++i)
+ {
+ ContOrthoPoly1D::MeasurePoint p(points[2*i]*range + xmin, points[2*i+1]);
+ measure.push_back(p);
+ }
+ ContOrthoPoly1D poly(maxdeg, measure);
+ poly.weightCoeffs(&coeffs[0], maxdeg);
+
+ for (unsigned i=0; i<npoints; ++i)
+ {
+ const double x = measure[i].first;
+ const double w = measure[i].second;
+ const double fvalue = poly.series(&coeffs[0], maxdeg, x);
+ CHECK_CLOSE(w, fvalue, eps);
+ }
+ }
+ }
+
TEST(ContOrthoPoly1D_cov8)
{
const double eps = 1.0e-12;
MersenneTwister rng;
const unsigned npoints = 64;
const unsigned maxdeg = 6;
const unsigned ntries = 5;
const unsigned degtries = 100;
std::vector<double> points(npoints);
for (unsigned itry=0; itry<ntries; ++itry)
{
rng.run(&points[0], npoints, npoints);
ContOrthoPoly1D poly(maxdeg, points);
for (unsigned ideg=0; ideg<degtries; ++ideg)
{
unsigned d[8];
bool has0 = true;
while (has0)
{
for (unsigned i=0; i<8; ++i)
d[i] = (maxdeg + 0.99)*rng();
has0 = (d[0] == 0 && d[1] == 0) ||
(d[2] == 0 && d[3] == 0) ||
(d[4] == 0 && d[5] == 0) ||
(d[6] == 0 && d[7] == 0);
}
const double cov8 = poly.cov8(d[0],d[1],d[2],d[3],d[4],d[5],d[6],d[7]);
const double check = poly.slowCov8(d[0],d[1],d[2],d[3],d[4],d[5],d[6],d[7]);
CHECK_CLOSE(check, cov8, eps);
}
}
}
TEST(ContOrthoPoly1D_epsCovariance)
{
const double eps = 1.0e-15;
const unsigned npoints = 64;
const unsigned maxdeg = 6;
const unsigned ntries = 5;
MersenneTwister rng;
std::vector<double> points(npoints);
for (unsigned itry=0; itry<ntries; ++itry)
{
rng.run(&points[0], npoints, npoints);
ContOrthoPoly1D poly(maxdeg, points);
for (unsigned i=0; i<=maxdeg; ++i)
for (unsigned j=0; j<=i; ++j)
{
CHECK(poly.epsExpectation(i, j, false) == 0.0);
for (unsigned k=0; k<=maxdeg; ++k)
for (unsigned m=0; m<=k; ++m)
{
const double eps1 = poly.empiricalKroneckerCovariance(i, j, k, m);
const double eps2 = poly.epsCovariance(i, j, k, m, false);
CHECK_CLOSE(eps1, eps2, eps);
}
}
}
}
}

File Metadata

Mime Type
text/x-diff
Expires
Sat, May 3, 6:24 AM (1 d, 9 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
4962065
Default Alt Text
(6 KB)

Event Timeline