nmrBernsteinPolynomialLineIntegral.h

Go to the documentation of this file.

/* -*- Mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-    */
/* ex: set filetype=cpp softtabstop=4 shiftwidth=4 tabstop=4 cindent expandtab: */

/*

  Author(s):    Ofri Sadowsky
  Created on:   2003-05-06

  (C) Copyright 2003-2007 Johns Hopkins University (JHU), All Rights
  Reserved.

--- begin cisst license - do not edit ---

This software is provided "as is" under an open source license, with
no warranty.  The complete license can be found in license.txt and
http://www.cisst.org/cisst/license.txt.

--- end cisst license ---
*/


#ifndef _nmrBernsteinPolynomialLineIntegral_h
#define _nmrBernsteinPolynomialLineIntegral_h

#include <cisstNumerical/nmrPolynomialBase.h>
#include <cisstNumerical/nmrBernsteinPolynomial.h>

#include <cisstNumerical/nmrExport.h>

#define CACHE_POWER_BASIS 1


#define USE_INTEGRAND_TABLEAU 1

class CISST_EXPORT nmrBernsteinPolynomialLineIntegral
{
public:
    // The value evaluated by the integral
    typedef nmrPolynomialBase::ValueType ValueType;

    // The type of value of a variable of the integrand function.
    typedef nmrPolynomialBase::VariableType VariableType;

    typedef nmrPolynomialBase::VariableIndexType VariableIndexType;

    // The type of a coefficient in the integrand function
    typedef nmrPolynomialBase::CoefficientType CoefficientType;

    typedef nmrPolynomialBase::PowerType PowerType;

    typedef nmrBernsteinPolynomial IntegrandType;

    typedef nmrPolynomialTermPowerIndex::MultinomialCoefficientType MultinomialType;


    // See the usage of roundoff tolerance to reduce the number of iterations.
    static const double DefaultRoundoffTolerance;

    // Profiling information that's accumulated during the integration.
    // Counts the invocations of different operations.
    // \note deprecated
    struct ProfilingInfo {
        unsigned long IntegrandTermCounter;
        unsigned long EvaluatePowerProductCounter;
        unsigned long TimeOfSingleTermIntegration;

        ProfilingInfo()
            : IntegrandTermCounter(0)
            , EvaluatePowerProductCounter(0)
            , TimeOfSingleTermIntegration(0)
        {}
    };

    // Initialize the integrator for a specific polynomial.
    nmrBernsteinPolynomialLineIntegral(const IntegrandType & integrand)
        : Integrand(integrand)
#if CACHE_POWER_BASIS
        , PowersAtP0(Integrand.GetNumVariables(), Integrand.GetMaxDegree())
        , PowersAtP1(Integrand.GetNumVariables(), Integrand.GetMaxDegree())
#endif  // CACHE_POWER_BASIS
    {
#if USE_INTEGRAND_TABLEAU
        InitializePolynomialIntegrationTableau(Integrand, IntegrationTableau);
#endif
    }

    void UpdateIntegrationTableau()
    {
        InitializePolynomialIntegrationTableau(Integrand, IntegrationTableau);
    }

    ValueType EvaluateForSegment(
        const VariableType p0[], const VariableType p1[], 
        const VariableType gradientNorm,
        const double roundoffTolerance = DefaultRoundoffTolerance,
        CoefficientType const * coefficients = NULL)
#if !CACHE_POWER_BASIS
        const
#endif
        ;

    const IntegrandType & GetIntegrand() const
    {
        return Integrand;
    }

protected:
    // The integrand polynomial
    const IntegrandType & Integrand;

#if CACHE_POWER_BASIS
    // Auxiliary objects to compute and contain the power basis of the variables
    nmrMultiVariablePowerBasis::StandardPowerBasis PowersAtP0;
    nmrMultiVariablePowerBasis::StandardPowerBasis PowersAtP1;
#endif  // CACHE_POWER_BASIS

    typedef std::vector<PowerType> PowerIndexType;
    typedef struct {
        PowerIndexType P0PowerIndex;
        MultinomialType P0MultinomialFactor;
        PowerIndexType P1PowerIndex;
        MultinomialType P1MultinomialFactor;
    } TermSummationElement;
    typedef std::vector<TermSummationElement> SingleTermIntegrationTableau;

    typedef std::vector<SingleTermIntegrationTableau> PolynomialIntegrationTableau;

    PolynomialIntegrationTableau IntegrationTableau;


    static ValueType IntegrateSingleTerm(
        const CoefficientType termCoefficient,
        const CoefficientType gradientNorm,
        const nmrPolynomialTermPowerIndex & termIndex,
        const nmrMultiVariablePowerBasis & powersAtP0,
        const nmrMultiVariablePowerBasis & powersAtP1,
        const double roundoffTolerance = DefaultRoundoffTolerance);



    typedef std::vector<VariableIndexType> IndexContainerForZeroVariables;

    static ValueType IntegrateSingleTerm(
        const CoefficientType termCoefficient,
        const CoefficientType gradientNorm,
        const PowerType termDegree,
        const SingleTermIntegrationTableau & termIntegrationTableau,
        const nmrMultiVariablePowerBasis & powersAtP0,
        const nmrMultiVariablePowerBasis & powersAtP1,
        const IndexContainerForZeroVariables & zerosOfP0,
        const IndexContainerForZeroVariables & zerosOfP1);

    static void InitializePolynomialIntegrationTableau( 
        const IntegrandType & integrand,
        PolynomialIntegrationTableau & tableau );

    static void InitializeSingleTermIntegrationTableau(
        const nmrPolynomialTermPowerIndex & termIndex,
        SingleTermIntegrationTableau & termTableau );

};

#endif