Add another likelihood marginalization test

This one is now testing where we treat the marginal parameter pdf
as the weighting function in the quadrature rule. Here, we use
a zero mean/unit variance Gaussian for the 2D marginal parameter space
and then use Gauss-Hermite quadrature.

Author: pbauman

test/unit/marg_likelihood.C

@@ -35,6 +35,7 @@
 #include <queso/LikelihoodBase.h>
 #include <queso/UniformVectorRV.h>
 #include <queso/1DQuadrature.h>
+#include <queso/GaussianVectorRV.h>
 
 #include <queso/GslVector.h>
 #include <queso/GslMatrix.h>
@@ -120,6 +121,46 @@ namespace QUESOTesting
 
   };
 
+  template <class V, class M>
+  class GaussModelLikelihood : public TestlingLikelihoodBase<V,M>
+  {
+  public:
+
+    GaussModelLikelihood(const char * prefix,
+                         const QUESO::VectorSet<V, M> & domainSet,
+                         const V & observations,
+                         typename QUESO::SharedPtr<QUESO::BaseVectorRV<V,M> >::Type & marg_param_pdf,
+                         typename QUESO::SharedPtr<QUESO::MultiDQuadratureBase<V,M> >::Type & marg_integration,
+                         bool marg_pdf_is_weight_func)
+      : TestlingLikelihoodBase<V,M>(prefix,domainSet,observations,marg_param_pdf,marg_integration,marg_pdf_is_weight_func)
+    {}
+
+    virtual void evaluateModel(const V & domainVector,
+                               const V & marginalVector,
+                               V & modelOutput) const
+    {
+      queso_require_equal_to(1,domainVector.sizeGlobal());
+      queso_require_equal_to(2,marginalVector.sizeGlobal());
+      queso_require_equal_to(1,modelOutput.sizeGlobal());
+
+      double q1 = marginalVector[0];
+      double q2 = marginalVector[1];
+      double m = domainVector[0];
+
+      modelOutput[0] = 3.14*m + 1.1*q1 + 2.2*q2*q2;
+    }
+
+    virtual double lnMargLikelihood( const V & domainVector ) const
+    {
+      queso_require_equal_to(1,domainVector.sizeGlobal());
+
+      double m = domainVector[0];
+
+      return std::log(3.14*m*M_PI + 2.2*M_PI/2.0);
+    }
+
+  };
+
   template <class V, class M>
   class MarginalLikelihoodTestBase : public CppUnit::TestCase
   {
@@ -184,6 +225,64 @@ namespace QUESOTesting
                                           param_space, tol, likelihood );
     }
 
+    void test_linear_func_gaussian_marg_space()
+    {
+      // Instantiate the parameter space
+      unsigned int param_dim = 1;
+      unsigned int marg_dim = 2;
+      QUESO::VectorSpace<V,M> param_space( (*this->_env), "param_", param_dim, NULL);
+      QUESO::VectorSpace<V,M> marg_space( (*this->_env), "marg_", marg_dim, NULL);
+
+      double param_min_domain_value = 0.0;
+      double param_max_domain_value = 1.0;
+
+      double marg_min_domain_value = -INFINITY;
+      double marg_max_domain_value = INFINITY;
+
+      typename QUESO::ScopedPtr<V>::Type param_min_values( param_space.newVector(param_min_domain_value) );
+      typename QUESO::ScopedPtr<V>::Type param_max_values( param_space.newVector(param_max_domain_value) );
+
+      typename QUESO::ScopedPtr<V>::Type marg_min_values( marg_space.newVector(marg_min_domain_value) );
+      typename QUESO::ScopedPtr<V>::Type marg_max_values( marg_space.newVector(marg_max_domain_value) );
+
+      QUESO::BoxSubset<V,M> param_domain( "param_domain_", param_space, (*param_min_values), (*param_max_values) );
+      QUESO::BoxSubset<V,M> marg_domain( "marg_domain_", marg_space, (*marg_min_values), (*marg_max_values) );
+
+      typename QUESO::ScopedPtr<V>::Type var_vec( param_space.newVector(1.0) );
+
+      // Zero mean, unit variance Gaussian
+      typename QUESO::SharedPtr<QUESO::BaseVectorRV<V,M> >::Type
+        marg_param_rv( new QUESO::GaussianVectorRV<V,M>("marg_param_rv_", marg_domain,
+                                                        marg_space.zeroVector(),
+                                                        *var_vec) );
+
+      const V & data = param_space.zeroVector();
+
+      unsigned int int_order = 1;
+
+      QUESO::SharedPtr<QUESO::Base1DQuadrature>::Type
+        qrule_1d( new QUESO::GaussianHermite1DQuadrature(0.0,1.0,int_order) );
+
+      std::vector<QUESO::SharedPtr<QUESO::Base1DQuadrature>::Type> all_1d_qrules(marg_dim,qrule_1d);
+
+      typename QUESO::SharedPtr<QUESO::MultiDQuadratureBase<V,M> >::Type
+        marg_integration( new QUESO::TensorProductQuadrature<V,M>(marg_domain,all_1d_qrules) );
+
+      GaussModelLikelihood<V,M> likelihood( "likelihood_test_",
+                                             param_domain,
+                                             data,
+                                             marg_param_rv,
+                                             marg_integration,
+                                             true );
+
+      // Test marginalized likelihood over a range of points in the parameter domain
+      unsigned int n_intervals = 20;
+      double tol = std::numeric_limits<double>::epsilon()*50;
+
+      this->test_likelihood_values_range( n_intervals, param_min_domain_value, param_max_domain_value,
+                                          param_space, tol, likelihood );
+    }
+
   protected:
 
     QUESO::EnvOptionsValues _options;
@@ -216,6 +315,7 @@ namespace QUESOTesting
     CPPUNIT_TEST_SUITE( MarginalLikelihoodGslTest );
 
     CPPUNIT_TEST( test_linear_func_uniform_marg_space );
+    CPPUNIT_TEST( test_linear_func_gaussian_marg_space );
 
     CPPUNIT_TEST_SUITE_END();
   };