From: Ralf Hartmann Date: Fri, 15 Dec 2000 06:59:08 +0000 (+0000) Subject: Initial implementation of new classes Polynomial (includes Horner scheme) and Lagrang... X-Git-Tag: v8.0.0~19899 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=19e3520c0dc6af8395ecb83482d141c5d4ec8aef;p=dealii.git Initial implementation of new classes Polynomial (includes Horner scheme) and LagrangeEquidistant for Q0...Q4. git-svn-id: https://svn.dealii.org/trunk@3534 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/base/include/base/polynomial.h b/deal.II/base/include/base/polynomial.h new file mode 100644 index 0000000000..2bc65b764c --- /dev/null +++ b/deal.II/base/include/base/polynomial.h @@ -0,0 +1,95 @@ +/*---------------------------- polynomial.h ---------------------------*/ +/* $Id$ */ +/* Ralf Hartmann, University of Heidelberg */ +#ifndef __polynomial_H +#define __polynomial_H +/*---------------------------- polynomial.h ---------------------------*/ + + +#include + +#include + +/** + * Base class for all 1D polynomials. + * + * @author Ralf Hartmann, 2000 + */ +class Polynomial +{ + public: + /** + * Constructor. + */ + Polynomial(const vector &a); + + /** + * Returns the values and the + * derivatives of the @p{Polynomial} + * at point @p{x}. @p{values[i], + * i=0,...,values.size()} + * includes the @p{i}th + * derivative. + * + * This function uses the Horner + * scheme. + */ + void value(double x, vector &values) const; + + protected: + + /** + * Coefficients of the polynomial + * $\sum_ia_ix^i$. This vector is + * filled by the constructor of + * derived classes. + */ + const vector coefficients; +}; + + + +/** + * Class of Lagrange polynomials with equidistant interpolation + * points. The polynomial of order @p{n} has got @p{n+1} interpolation + * points. The interpolation points are x=0, x=1 and x=intermediate + * points in ]0,1[ in ascending order. This order gives an index to + * each interpolation point. A Lagrangian polynomial equals 1 at one + * interpolation point that is called `support point', and 0 at all other + * interpolation points. + * + * @author Ralf Hartmann, 2000 + */ +class LagrangeEquidistant: public Polynomial +{ + public: + /** + * Constructor. Takes the order + * @p{n} of the Lagrangian + * polynom and the index + * @p{support_point} of the + * support point. Fills the + * @p{coefficients} of the base + * class @p{Polynomial}. + */ + LagrangeEquidistant(unsigned int n, unsigned int support_point); + + private: + + /** + * Computes the @p{coefficients} + * of the base class + * @p{Polynomial}. This function + * is static to allow the + * @p{coefficients} to be a + * @p{const} vector. + */ + static vector compute_coefficients(unsigned int n, unsigned int support_point); +}; + + + +/*---------------------------- polynomial.h ---------------------------*/ +/* end of #ifndef __polynomial_H */ +#endif +/*---------------------------- polynomial.h ---------------------------*/ diff --git a/deal.II/base/source/polynomial.cc b/deal.II/base/source/polynomial.cc new file mode 100644 index 0000000000..62401eb3ce --- /dev/null +++ b/deal.II/base/source/polynomial.cc @@ -0,0 +1,180 @@ +//---------------------------- polynomial.cc --------------------------- +// $Id$ */ +// Version: $Name$ +// +// Copyright (C) 2000 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- polynomial.cc ----------------------- + + +#include + + +Polynomial::Polynomial(const vector &a): + coefficients(a) +{} + + +void Polynomial::value(double x, vector &values) const +{ + const unsigned int m=coefficients.size(); + vector a(coefficients); + // Horner scheme + unsigned int j_faculty=1; + for (unsigned int j=0; j=static_cast(j); --k) + a[k]+=x*a[k+1]; + values[j]=j_faculty*a[j]; + + j_faculty*=j+1; + } +} + + +// ------------------------------------------------------------ // + + + +LagrangeEquidistant::LagrangeEquidistant(unsigned int n, unsigned int support_point): + Polynomial(compute_coefficients(n,support_point)) +{} + + +vector LagrangeEquidistant::compute_coefficients(unsigned int n, unsigned int support_point) +{ + vector a; + a.resize(n+1); + Assert(support_point