diff --git a/GPU/TPCFastTransformation/MultivariatePolynomial.h b/GPU/TPCFastTransformation/MultivariatePolynomial.h index 4fd2157409133..1454194f9e3b4 100644 --- a/GPU/TPCFastTransformation/MultivariatePolynomial.h +++ b/GPU/TPCFastTransformation/MultivariatePolynomial.h @@ -56,7 +56,7 @@ class MultivariatePolynomial : public FlatObject, public MultivariatePolynomialH /// constructor for compile time evaluation of polynomial formula template ::type = 0> - MultivariatePolynomial() : mNParams{this->getNParameters(Degree, Dim, InteractionOnly)} + MultivariatePolynomial() : mNParams{this->template getNParameters()} { construct(); } diff --git a/GPU/TPCFastTransformation/MultivariatePolynomialHelper.h b/GPU/TPCFastTransformation/MultivariatePolynomialHelper.h index 6cefab5d1c573..2dd186a859ab0 100644 --- a/GPU/TPCFastTransformation/MultivariatePolynomialHelper.h +++ b/GPU/TPCFastTransformation/MultivariatePolynomialHelper.h @@ -57,7 +57,7 @@ struct MultivariatePolynomialContainer { class MultivariatePolynomialParametersHelper { public: - /// \returns number of parameters for given dimension and degree of polynomials + /// \returns number of parameters for given dimension and degree of polynomials at compile time /// calculates the number of parameters for a multivariate polynomial for given degree: nParameters = (n+d-1 d) -> binomial coefficient /// see: https://mathoverflow.net/questions/225953/number-of-polynomial-terms-for-certain-degree-and-certain-number-of-variables template @@ -70,7 +70,19 @@ class MultivariatePolynomialParametersHelper } } - /// \returns the number of parameters for interaction terms only (see: https://en.wikipedia.org/wiki/Combination) + /// \returns number of parameters for given dimension and degree of polynomials + /// calculates the number of parameters for a multivariate polynomial for given degree: nParameters = (n+d-1 d) -> binomial coefficient + /// see: https://mathoverflow.net/questions/225953/number-of-polynomial-terms-for-certain-degree-and-certain-number-of-variables + GPUd() static constexpr uint32_t getNParametersAllTerms(uint32_t degree, uint32_t dim) + { + if (degree == 0) { + return binomialCoeff(dim - 1, 0); + } else { + return binomialCoeff(dim - 1 + degree, degree) + getNParametersAllTerms(degree - 1, dim); + } + } + + /// \returns the number of parameters at compile time for interaction terms only (see: https://en.wikipedia.org/wiki/Combination) template GPUd() static constexpr uint32_t getNParametersInteractionOnly() { @@ -81,6 +93,16 @@ class MultivariatePolynomialParametersHelper } } + /// \returns the number of parameters for interaction terms only (see: https://en.wikipedia.org/wiki/Combination) + GPUd() static constexpr uint32_t getNParametersInteractionOnly(uint32_t degree, uint32_t dim) + { + if (degree == 0) { + return binomialCoeff(dim - 1, 0); + } else { + return binomialCoeff(dim, degree) + getNParametersInteractionOnly(degree - 1, dim); + } + } + template GPUd() static constexpr uint32_t getNParameters() { @@ -91,8 +113,17 @@ class MultivariatePolynomialParametersHelper } } + GPUd() static constexpr uint32_t getNParameters(uint32_t degree, uint32_t dim, bool interactionOnly) + { + if (interactionOnly) { + return getNParametersInteractionOnly(degree, dim); + } else { + return getNParametersAllTerms(degree, dim); + } + } + private: - /// calculate factorial of n + /// calculate factorial of n at compile time /// \return returns n! template GPUd() static constexpr uint32_t factorial() @@ -104,13 +135,24 @@ class MultivariatePolynomialParametersHelper } } - /// calculates binomial coefficient + /// calculate factorial of n + /// \return returns n! + GPUd() static constexpr uint32_t factorial(uint32_t n) { return n == 0 || n == 1 ? 1 : n * factorial(n - 1); } + + /// calculates binomial coefficient at compile time /// \return returns (n k) template GPUd() static constexpr uint32_t binomialCoeff() { return factorial() / (factorial() * factorial()); } + + /// calculates binomial coefficient + /// \return returns (n k) + GPUd() static constexpr uint32_t binomialCoeff(uint32_t n, uint32_t k) + { + return factorial(n) / (factorial(k) * factorial(n - k)); + } }; /// Helper struct for evaluating a multidimensional polynomial using compile time evaluated formula diff --git a/GPU/TPCFastTransformation/NDPiecewisePolynomials.inc b/GPU/TPCFastTransformation/NDPiecewisePolynomials.inc index 2538e30056448..1cbcf9cd8e23e 100644 --- a/GPU/TPCFastTransformation/NDPiecewisePolynomials.inc +++ b/GPU/TPCFastTransformation/NDPiecewisePolynomials.inc @@ -165,7 +165,7 @@ void NDPiecewisePolynomials::performFits(const std // check if data points are in the grid if (index == indexClamped) { // index of the polyniomial - const uint32_t idx = getDataIndex(index.data()) / MultivariatePolynomialParametersHelper::getNParameters(Degree, Dim, InteractionOnly); + const uint32_t idx = getDataIndex(index.data()) / MultivariatePolynomialParametersHelper::getNParameters(); // store index to data point dataPointsIndices[idx].emplace_back(i); @@ -216,7 +216,7 @@ void NDPiecewisePolynomials::performFits(const std const auto params = MultivariatePolynomialHelper<0, 0, false>::fit(fitter, xCords, response, error, true); // store parameters - std::copy(params.begin(), params.end(), &mParams[i * MultivariatePolynomialParametersHelper::getNParameters(Degree, Dim, InteractionOnly)]); + std::copy(params.begin(), params.end(), &mParams[i * MultivariatePolynomialParametersHelper::getNParameters()]); } }