Utilities/test/testFunctions.cc

0001 #include <cppunit/extensions/HelperMacros.h>
0002 #include "PhysicsTools/Utilities/interface/Operations.h"
0003 #include "PhysicsTools/Utilities/interface/Gaussian.h"
0004 #include "PhysicsTools/Utilities/interface/Exponential.h"
0005 #include "PhysicsTools/Utilities/interface/Identity.h"
0006 #include "PhysicsTools/Utilities/interface/Composition.h"
0007 #include "PhysicsTools/Utilities/interface/Convolution.h"
0008 #include "PhysicsTools/Utilities/interface/Functions.h"
0009 #include "PhysicsTools/Utilities/interface/Variables.h"
0010 #include "PhysicsTools/Utilities/interface/Numerical.h"
0011 #include "PhysicsTools/Utilities/interface/Fraction.h"
0012 #include "PhysicsTools/Utilities/interface/Expression.h"
0013 #include "PhysicsTools/Utilities/interface/FunctClone.h"
0014 #include <iostream>
0015
0016 class testFunctions : public CppUnit::TestFixture {
0017   CPPUNIT_TEST_SUITE(testFunctions);
0018   CPPUNIT_TEST(checkAll);
0019   CPPUNIT_TEST_SUITE_END();
0020
0021 public:
0022   void setUp() {}
0023   void tearDown() {}
0024   void checkAll();
0025 };
0026
0027 CPPUNIT_TEST_SUITE_REGISTRATION(testFunctions);
0028
0029 struct TestFun {
0030   TestFun() : gauss_(0, 1) {}
0031   double operator()(double x) {
0032     ++counter_;
0033     return gauss_(x);
0034   }
0035   void reset() { counter_ = 0; }
0036   static size_t counter_;
0037
0038 private:
0039   funct::Gaussian gauss_;
0040 };
0041
0042 size_t TestFun::counter_ = 0;
0043
0044 void testFunctions::checkAll() {
0045   using namespace funct;
0046   {
0047     Gaussian g1(0, 1);
0048     Gaussian g2(1, 1);
0049     Identity i;
0050     const double epsilon = 1.e-6;
0051     Sum<Gaussian, Gaussian>::type g1plus2 = g1 + g2;
0052     Product<Gaussian, Gaussian>::type g1times2 = g1 * g2;
0053     Difference<Gaussian, Gaussian>::type g1minus2 = g1 - g2;
0054     Ratio<Gaussian, Gaussian>::type g1over2 = g1 / g2;
0055     Minus<Gaussian>::type gm1 = -g1;
0056     Composition<Identity, Gaussian>::type gg1 = compose(i, g1);
0057     double x = 0.5;
0058     CPPUNIT_ASSERT(std::abs(g1plus2(x) - (g1(x) + g2(x))) < epsilon);
0059     CPPUNIT_ASSERT(std::abs(g1times2(x) - (g1(x) * g2(x))) < epsilon);
0060     CPPUNIT_ASSERT(std::abs(g1minus2(x) - (g1(x) - g2(x))) < epsilon);
0061     CPPUNIT_ASSERT(std::abs(g1over2(x) - (g1(x) / g2(x))) < epsilon);
0062     CPPUNIT_ASSERT(std::abs(gm1(x) - (-g1(x))) < epsilon);
0063     Convolution<Gaussian, Gaussian, TrapezoidIntegrator>::type ggt(g1, g1, -5, 5, TrapezoidIntegrator(1000));
0064     CPPUNIT_ASSERT(std::abs(ggt(0) - g1(0) / sqrt(2.0)) < epsilon);
0065     Convolution<Gaussian, Gaussian, GaussLegendreIntegrator>::type gggl(
0066         g1, g1, -5, 5, GaussLegendreIntegrator(1000, epsilon));
0067     CPPUNIT_ASSERT(std::abs(gggl(0) - g1(0) / sqrt(2.0)) < epsilon);
0068     CPPUNIT_ASSERT(gg1(0) == g1(0));
0069   }
0070   {
0071     double value = 0.123;
0072     double epsilon = 1.e-8;
0073     X x;
0074     x = value;
0075
0076     Exp<X>::type f_exp = exp(x);
0077     CPPUNIT_ASSERT(std::abs(f_exp() - exp(value)) < epsilon);
0078     Sin<X>::type f_sin = sin(x);
0079     CPPUNIT_ASSERT(std::abs(f_sin() - sin(value)) < epsilon);
0080     Cos<X>::type f_cos = cos(x);
0081     CPPUNIT_ASSERT(std::abs(f_cos() - cos(value)) < epsilon);
0082     Log<X>::type f_log = log(x);
0083     CPPUNIT_ASSERT(std::abs(f_log() - log(value)) < epsilon);
0084   }
0085   {
0086     Numerical<1> _1;
0087     CPPUNIT_ASSERT(_1 == 1);
0088     Numerical<2> _2;
0089     CPPUNIT_ASSERT(_2 == 2);
0090     Numerical<3> _3;
0091     CPPUNIT_ASSERT(_3 == 3);
0092     CPPUNIT_ASSERT(num<1>() == 1);
0093     CPPUNIT_ASSERT(num<2>() == 2);
0094     CPPUNIT_ASSERT(num<3>() == 3);
0095   }
0096   {
0097     Fraction<1, 2>::type _1_2;
0098     CPPUNIT_ASSERT(_1_2 == 0.5);
0099   }
0100   {
0101     X x = 3.141516;
0102     Expression f = sin(x) * cos(x);
0103     const double epsilon = 1.e-4;
0104     CPPUNIT_ASSERT(std::abs(f() - (sin(x) * cos(x))()) < epsilon);
0105   }
0106   {
0107     TestFun f;
0108     Master<TestFun> g(f);
0109     Slave<TestFun> g1(g), g2(g);
0110     const double epsilon = 1.e-5;
0111     double y, y1, y2, x;
0112     CPPUNIT_ASSERT(f.counter_ == 0);
0113     x = 0.5;
0114     y = g(x), y1 = g1(x), y2 = g2(x);
0115     CPPUNIT_ASSERT(y == y1 && y1 == y2);
0116     CPPUNIT_ASSERT(f.counter_ == 1);
0117     CPPUNIT_ASSERT(std::abs(f(x) - y) < epsilon);
0118     f.reset();
0119     x = 1.5;
0120     y1 = g1(x), y = g(x), y2 = g2(x);
0121     CPPUNIT_ASSERT(y == y1 && y1 == y2);
0122     CPPUNIT_ASSERT(f.counter_ == 1);
0123     CPPUNIT_ASSERT(std::abs(f(x) - y) < epsilon);
0124     f.reset();
0125     x = 0.765;
0126     y2 = g2(x), y1 = g1(x), y = g(x);
0127     CPPUNIT_ASSERT(y == y1 && y1 == y2);
0128     CPPUNIT_ASSERT(f.counter_ == 1);
0129     CPPUNIT_ASSERT(std::abs(f(x) - y) < epsilon);
0130     f.reset();
0131     g(0.5);
0132     CPPUNIT_ASSERT(f.counter_ == 1);
0133     g(0.5);
0134     CPPUNIT_ASSERT(f.counter_ == 2);
0135     g1(0.5);
0136     CPPUNIT_ASSERT(f.counter_ == 2);
0137     g1(0.5);
0138     CPPUNIT_ASSERT(f.counter_ == 3);
0139     g(0.5);
0140     CPPUNIT_ASSERT(f.counter_ == 3);
0141     f.reset();
0142     // odd case: slaves don't catch a change in values
0143     x = 0.123;
0144     y = g(x), y1 = g1(0.5), y2 = g2(0.7);
0145     CPPUNIT_ASSERT(y == y1 && y1 == y2);
0146     CPPUNIT_ASSERT(f.counter_ == 1);
0147     f.reset();
0148     CPPUNIT_ASSERT(std::abs(f(x) - y) < epsilon);
0149   }
0150 }