00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012 #ifndef __deal2__polynomial_space_h
00013 #define __deal2__polynomial_space_h
00014
00015
00016 #include <deal.II/base/config.h>
00017 #include <deal.II/base/exceptions.h>
00018 #include <deal.II/base/tensor.h>
00019 #include <deal.II/base/point.h>
00020 #include <deal.II/base/polynomial.h>
00021 #include <deal.II/base/smartpointer.h>
00022
00023 #include <vector>
00024
00025 DEAL_II_NAMESPACE_OPEN
00026
00083 template <int dim>
00084 class PolynomialSpace
00085 {
00086 public:
00093 static const unsigned int dimension = dim;
00094
00109 template <class Pol>
00110 PolynomialSpace (const std::vector<Pol> &pols);
00111
00116 template <class STREAM>
00117 void output_indices(STREAM &out) const;
00118
00126 void set_numbering(const std::vector<unsigned int> &renumber);
00127
00153 void compute (const Point<dim> &unit_point,
00154 std::vector<double> &values,
00155 std::vector<Tensor<1,dim> > &grads,
00156 std::vector<Tensor<2,dim> > &grad_grads) const;
00157
00165 double compute_value (const unsigned int i,
00166 const Point<dim> &p) const;
00167
00175 Tensor<1,dim> compute_grad (const unsigned int i,
00176 const Point<dim> &p) const;
00177
00186 Tensor<2,dim> compute_grad_grad (const unsigned int i,
00187 const Point<dim> &p) const;
00188
00201 unsigned int n () const;
00202
00213 unsigned int degree () const;
00214
00220 static unsigned int compute_n_pols (const unsigned int n);
00221
00222 protected:
00223
00233 void compute_index (const unsigned int n,
00234 unsigned int (&index)[dim]) const;
00235
00236 private:
00242 const std::vector<Polynomials::Polynomial<double> > polynomials;
00243
00249 const unsigned int n_pols;
00250
00255 std::vector<unsigned int> index_map;
00256
00261 std::vector<unsigned int> index_map_inverse;
00262 };
00263
00264
00265
00266
00267 template <>
00268 void PolynomialSpace<1>::compute_index(const unsigned int n,
00269 unsigned int (&index)[1]) const;
00270 template <>
00271 void PolynomialSpace<2>::compute_index(const unsigned int n,
00272 unsigned int (&index)[2]) const;
00273 template <>
00274 void PolynomialSpace<3>::compute_index(const unsigned int n,
00275 unsigned int (&index)[3]) const;
00276
00277
00278
00279
00280
00281 template <int dim>
00282 template <class Pol>
00283 PolynomialSpace<dim>::PolynomialSpace (const std::vector<Pol> &pols)
00284 :
00285 polynomials (pols.begin(), pols.end()),
00286 n_pols (compute_n_pols(polynomials.size())),
00287 index_map(n_pols),
00288 index_map_inverse(n_pols)
00289 {
00290
00291
00292
00293
00294 for (unsigned int i=0; i<n_pols; ++i)
00295 {
00296 index_map[i]=i;
00297 index_map_inverse[i]=i;
00298 }
00299 }
00300
00301
00302 template<int dim>
00303 inline
00304 unsigned int
00305 PolynomialSpace<dim>::n() const
00306 {
00307 return n_pols;
00308 }
00309
00310
00311
00312 template<int dim>
00313 inline
00314 unsigned int
00315 PolynomialSpace<dim>::degree() const
00316 {
00317 return polynomials.size();
00318 }
00319
00320
00321 template <int dim>
00322 template <class STREAM>
00323 void
00324 PolynomialSpace<dim>::output_indices(STREAM &out) const
00325 {
00326 unsigned int ix[dim];
00327 for (unsigned int i=0; i<n_pols; ++i)
00328 {
00329 compute_index(i,ix);
00330 out << i << "\t";
00331 for (unsigned int d=0; d<dim; ++d)
00332 out << ix[d] << " ";
00333 out << std::endl;
00334 }
00335 }
00336
00337
00338 DEAL_II_NAMESPACE_CLOSE
00339
00340 #endif
00341