From: guido Date: Mon, 30 May 2005 13:17:12 +0000 (+0000) Subject: polynomial multiplication and Lagrange class started X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=fd1dda3cd8824555472f30beae3047361f295353;p=dealii-svn.git polynomial multiplication and Lagrange class started git-svn-id: https://svn.dealii.org/trunk@10776 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/base/include/base/polynomial.h b/deal.II/base/include/base/polynomial.h index 6c1365a1a2..9096b55076 100644 --- a/deal.II/base/include/base/polynomial.h +++ b/deal.II/base/include/base/polynomial.h @@ -21,6 +21,7 @@ #include +template class Point; /** * A namespace in which classes relating to the description of @@ -163,6 +164,11 @@ namespace Polynomials */ Polynomial& operator *= (const double s); + /** + * Multiply with another polynomial. + */ + Polynomial& operator *= (const Polynomial& p); + /** * Add a second polynomial. */ @@ -333,7 +339,31 @@ namespace Polynomials const unsigned int support_point); }; - +/** + * Lagrange polynomials for an arbistrary set of interpolation points. + * + * @author Guido Kanschat, 2005 + */ + class Lagrange + { + public: + /** + * Given a set of points, this + * function returns all + * Lagrange polynomials for + * interpolation of these + * points. The number of + * polynomials is equal to the + * number of points and the + * maximum degree is one less. + */ + static + std::vector > + generate_complete_basis (const std::vector >& points); + }; + + + /** * Legendre polynomials of arbitrary degree on [0,1]. * diff --git a/deal.II/base/source/polynomial.cc b/deal.II/base/source/polynomial.cc index a9d3f3faba..1e3e6589bd 100644 --- a/deal.II/base/source/polynomial.cc +++ b/deal.II/base/source/polynomial.cc @@ -13,6 +13,7 @@ #include +#include #include #include @@ -167,6 +168,24 @@ namespace Polynomials } + template + Polynomial& + Polynomial::operator *= (const Polynomial& p) + { + // Degree of the product + unsigned int new_degree = this->degree() + p.degree(); + + std::vector new_coefficients(new_degree+1, 0.); + + for (unsigned int i=0; icoefficients.size(); ++j) + new_coefficients[i+j] += this->coefficients[j]*p.coefficients[i]; + this->coefficients = new_coefficients; + + return *this; + } + + template Polynomial& Polynomial::operator += (const Polynomial& p) @@ -589,6 +608,38 @@ namespace Polynomials } +//----------------------------------------------------------------------// + + + std::vector > + Lagrange::generate_complete_basis (const std::vector >& points) + { + std::vector > p(points.size()); + // polynomials are built as + // products of linear + // factors. The coefficient in + // front of the linear term is + // always 1. + std::vector linear(2, 1.); + // We start with a constant polynomial + std::vector one(1, 1.); + + for (unsigned int i=0;i(one); + for (unsigned int k=0;k factor(linear); + factor *= 1./(points[i](0)-points[k](0)); + p[i] *= factor; + } + } + return p; + } + // ------------------ class Legendre --------------- // diff --git a/tests/base/polynomial1d.cc b/tests/base/polynomial1d.cc index dcd7647e8f..e2a5f78250 100644 --- a/tests/base/polynomial1d.cc +++ b/tests/base/polynomial1d.cc @@ -77,6 +77,9 @@ void polynomial_arithmetic () p1 *= 2.; std::cerr << "*2" << std::endl; p1.print(std::cerr); + std::cerr << "*P2" << std::endl; + p2 *= p1; + p2.print(std::cerr); for (unsigned int i=0;i<7;++i) { diff --git a/tests/results/i686-pc-linux-gnu+gcc3.2/base/polynomial1d.output b/tests/results/i686-pc-linux-gnu+gcc3.2/base/polynomial1d.output index d73e32dc27..fb8dcb3a55 100644 --- a/tests/results/i686-pc-linux-gnu+gcc3.2/base/polynomial1d.output +++ b/tests/results/i686-pc-linux-gnu+gcc3.2/base/polynomial1d.output @@ -32,6 +32,16 @@ P1+P2+x^5 7.40000 x^2 7.40000 x^1 4.80000 x^0 +*P2 +8.00000 x^8 +4.80000 x^7 +40.8000 x^6 +54.4000 x^5 +82.5600 x^4 +91.9200 x^3 +64.8000 x^2 +41.8400 x^1 +13.4400 x^0 derive 10.0000 x^4 0.00000 x^3