00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012 #ifndef __deal2__tensor_product_polynomials_h
00013 #define __deal2__tensor_product_polynomials_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
00022 #include <vector>
00023
00024 DEAL_II_NAMESPACE_OPEN
00025
00060 template <int dim>
00061 class TensorProductPolynomials
00062 {
00063 public:
00070 static const unsigned int dimension = dim;
00071
00081 template <class Pol>
00082 TensorProductPolynomials (const std::vector<Pol> &pols);
00083
00088 void output_indices(std::ostream &out) const;
00089
00097 void set_numbering(const std::vector<unsigned int> &renumber);
00098
00103 const std::vector<unsigned int> &get_numbering() const;
00104
00109 const std::vector<unsigned int> &get_numbering_inverse() const;
00110
00135 void compute (const Point<dim> &unit_point,
00136 std::vector<double> &values,
00137 std::vector<Tensor<1,dim> > &grads,
00138 std::vector<Tensor<2,dim> > &grad_grads) const;
00139
00163 double compute_value (const unsigned int i,
00164 const Point<dim> &p) const;
00165
00190 Tensor<1,dim> compute_grad (const unsigned int i,
00191 const Point<dim> &p) const;
00192
00218 Tensor<2,dim> compute_grad_grad (const unsigned int i,
00219 const Point<dim> &p) const;
00220
00226 unsigned int n () const;
00227
00228
00229 private:
00235 std::vector<Polynomials::Polynomial<double> > polynomials;
00236
00241 unsigned int n_tensor_pols;
00242
00247 std::vector<unsigned int> index_map;
00248
00253 std::vector<unsigned int> index_map_inverse;
00254
00265
00266
00267 void compute_index (const unsigned int i,
00268 unsigned int (&indices)[(dim>0?dim:1)]) const;
00269
00276 static
00277 unsigned int x_to_the_dim (const unsigned int x);
00278 };
00279
00280 #ifndef DOXYGEN
00281
00282 template <int dim>
00283 inline
00284 const std::vector<unsigned int> &
00285 TensorProductPolynomials<dim>::get_numbering() const
00286 {
00287 return index_map;
00288 }
00289
00290
00291 template <int dim>
00292 inline
00293 const std::vector<unsigned int> &
00294 TensorProductPolynomials<dim>::get_numbering_inverse() const
00295 {
00296 return index_map_inverse;
00297 }
00298
00299
00300 #endif // DOXYGEN
00301
00302
00322 template <int dim>
00323 class AnisotropicPolynomials
00324 {
00325 public:
00340 AnisotropicPolynomials (const std::vector<std::vector<Polynomials::Polynomial<double> > > &pols);
00341
00366 void compute (const Point<dim> &unit_point,
00367 std::vector<double> &values,
00368 std::vector<Tensor<1,dim> > &grads,
00369 std::vector<Tensor<2,dim> > &grad_grads) const;
00370
00395 double compute_value (const unsigned int i,
00396 const Point<dim> &p) const;
00397
00422 Tensor<1,dim> compute_grad (const unsigned int i,
00423 const Point<dim> &p) const;
00424
00450 Tensor<2,dim> compute_grad_grad (const unsigned int i,
00451 const Point<dim> &p) const;
00452
00460 unsigned int n () const;
00461
00462 private:
00468 std::vector<std::vector<Polynomials::Polynomial<double> > > polynomials;
00469
00477 unsigned int n_tensor_pols;
00478
00489 void compute_index (const unsigned int i,
00490 unsigned int (&indices)[dim]) const;
00491
00497 static
00498 unsigned int
00499 get_n_tensor_pols (const std::vector<std::vector<Polynomials::Polynomial<double> > > &pols);
00500 };
00501
00504 #ifndef DOXYGEN
00505
00506
00507
00508 template <>
00509 void
00510 TensorProductPolynomials<1>::compute_index(const unsigned int n,
00511 unsigned int (&index)[1]) const;
00512 template <>
00513 void
00514 TensorProductPolynomials<2>::compute_index(const unsigned int n,
00515 unsigned int (&index)[2]) const;
00516 template <>
00517 void
00518 TensorProductPolynomials<3>::compute_index(const unsigned int n,
00519 unsigned int (&index)[3]) const;
00520
00521
00522
00523
00524 template <int dim>
00525 inline
00526 unsigned int
00527 TensorProductPolynomials<dim>::
00528 x_to_the_dim (const unsigned int x)
00529 {
00530 unsigned int y = 1;
00531 for (int d=0; d<dim; ++d)
00532 y *= x;
00533 return y;
00534 }
00535
00536
00537
00538 template <int dim>
00539 template <class Pol>
00540 TensorProductPolynomials<dim>::
00541 TensorProductPolynomials(const std::vector<Pol> &pols)
00542 :
00543 polynomials (pols.begin(), pols.end()),
00544 n_tensor_pols(x_to_the_dim(pols.size())),
00545 index_map(n_tensor_pols),
00546 index_map_inverse(n_tensor_pols)
00547 {
00548
00549
00550
00551
00552 for (unsigned int i=0; i<n_tensor_pols; ++i)
00553 {
00554 index_map[i]=i;
00555 index_map_inverse[i]=i;
00556 }
00557 }
00558
00559
00560
00561 template <>
00562 void
00563 AnisotropicPolynomials<1>::compute_index(const unsigned int n,
00564 unsigned int (&index)[1]) const;
00565 template <>
00566 void
00567 AnisotropicPolynomials<2>::compute_index(const unsigned int n,
00568 unsigned int (&index)[2]) const;
00569 template <>
00570 void
00571 AnisotropicPolynomials<3>::compute_index(const unsigned int n,
00572 unsigned int (&index)[3]) const;
00573
00574 #endif // DOXYGEN
00575 DEAL_II_NAMESPACE_CLOSE
00576
00577 #endif
00578