Reference documentation for deal.II version GIT 3779fa9eb4 2023-09-28 13:00:02+00:00
polynomials_integrated_legendre_sz.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2015 - 2021 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15
16
17
18
20
22  : Polynomials::Polynomial<double>(get_coefficients(k))
23 {}
24
25
26
27 std::vector<double>
29 {
30  std::vector<double> coefficients(k + 1);
31
32  // first two polynomials are hard-coded:
33  if (k == 0)
34  {
35  coefficients[0] = -1.;
36  return coefficients;
37  }
38  else if (k == 1)
39  {
40  coefficients[0] = 0.;
41  coefficients[1] = 1.;
42  return coefficients;
43  }
44
45  // General formula is:
46  // k*L_{k}(x) = (2*k-3)*x*L_{k-1} - (k-3)*L_{k-2}.
47  std::vector<double> coefficients_km2 = get_coefficients(k - 2);
48  std::vector<double> coefficients_km1 = get_coefficients(k - 1);
49
50  const double a = 1.0 / k;
51  const double b = 2.0 * k - 3.0;
52  const double c = k - 3.0;
53
54  // To maintain stability, delay the division (multiplication by a) until the
55  // end.
56  for (unsigned int i = 1; i <= k - 2; ++i)
57  {
58  coefficients[i] = b * coefficients_km1[i - 1] - c * coefficients_km2[i];
59  }
60
61  coefficients[0] = -c * coefficients_km2[0];
62  coefficients[k] = b * coefficients_km1[k - 1];
63  coefficients[k - 1] = b * coefficients_km1[k - 2];
64
65  for (double &coefficient : coefficients)
66  {
67  coefficient *= a;
68  }
69
70  return coefficients;
71 }
72
73
74
75 std::vector<Polynomials::Polynomial<double>>
77 {
78  std::vector<Polynomials::Polynomial<double>> v;
79  v.reserve(degree + 1);
80  for (unsigned int i = 0; i <= degree; ++i)
81  {
82  v.push_back(IntegratedLegendreSZ(i));
83  }
84  return v;
85 }
86
87
88
static std::vector< Polynomials::Polynomial< double > > generate_complete_basis(const unsigned int degree)
static std::vector< double > get_coefficients(const unsigned int k)
std::vector< double > coefficients
Definition: polynomial.h:302
unsigned int degree() const
Definition: polynomial.h:800
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)