60 numIterations_(numIterations),
61 weight_(numIterations + 1),
62 root_(numIterations + 1)
64 calculateWeightAndRoot();
73 const std::vector<double>&
getWeight() const noexcept
84 const std::vector<double>&
getRoot() const noexcept
97 double derivative{ 0.0 };
106 Result(
const double val,
const double deriv) noexcept :
116 void calculateWeightAndRoot() noexcept
118 const auto numIterationsDouble =
static_cast<double>(numIterations_);
119 for (
uint32 step = 0; step <= numIterations_; ++step)
121 double root =
std::cos(
constants::pi * (
static_cast<double>(step) - 0.25) / (numIterationsDouble + 0.5));
122 Result result = calculatePolynomialValueAndDerivative(root);
124 double newtonRaphsonRatio;
127 newtonRaphsonRatio = result.value / result.derivative;
128 root -= newtonRaphsonRatio;
129 result = calculatePolynomialValueAndDerivative(root);
130 }
while (std::fabs(newtonRaphsonRatio) > EPSILON);
133 weight_[step] = 2.0 / ((1.0 -
utils::sqr(root)) * result.derivative * result.derivative);
144 Result calculatePolynomialValueAndDerivative(
const double x) noexcept
146 Result result(x, 0.0);
148 double value_minus_1 = 1.0;
150 for (
uint32 step = 2; step <= numIterations_; ++step)
152 const auto stepDouble =
static_cast<double>(step);
153 const double value = ((2.0 * stepDouble - 1.0) * x * result.value - (stepDouble - 1.0) * value_minus_1) / stepDouble;
154 result.derivative = stepDouble *
f * (x * value - result.value);
156 value_minus_1 = result.value;
157 result.value = value;
164 const double EPSILON{ 1
e-15 };
166 const uint32 numIterations_;
167 std::vector<double> weight_;
168 std::vector<double> root_;
183 const std::function<
double(
double)>&
f)
186 const std::vector<double>& weight = legendrePolynomial.
getWeight();
187 const std::vector<double>& root = legendrePolynomial.
getRoot();
189 const double width = 0.5 * (high - low);
190 const double mean = 0.5 * (low + high);
192 double gaussLegendre = 0.0;
193 for (
uint32 step = 1; step <= n; ++step)
195 gaussLegendre += weight[step] *
f(width * root[step] +
mean);
198 return gaussLegendre * width;