]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix the computation of the number of polynomials the ABF element has.
authorWolfgang Bangerth <bangerth@colostate.edu>
Mon, 9 Jan 2017 04:28:25 +0000 (21:28 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Mon, 9 Jan 2017 04:28:25 +0000 (21:28 -0700)
source/base/polynomials_abf.cc

index f2718a00a2175ef345c58d32bb8824535f46ac96..ae13827fa4cad300c80b61586829a343c264a45a 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2004 - 2015 by the deal.II authors
+// Copyright (C) 2004 - 2017 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -31,13 +31,21 @@ PolynomialsABF<dim>::PolynomialsABF (const unsigned int k)
 {
   std::vector<std::vector< Polynomials::Polynomial< double > > > pols(dim);
   pols[0] = Polynomials::LagrangeEquidistant::generate_complete_basis(k+2);
+
   if (k == 0)
     for (unsigned int d=1; d<dim; ++d)
       pols[d] = Polynomials::Legendre::generate_complete_basis(0);
   else
     for (unsigned int d=1; d<dim; ++d)
       pols[d] = Polynomials::LagrangeEquidistant::generate_complete_basis(k);
+
   polynomial_space = new AnisotropicPolynomials<dim>(pols);
+
+  // check that the dimensions match. we only store one of the 'dim'
+  // anisotropic polynomials that make up the vector-valued space, so
+  // multiply by 'dim'
+  Assert (dim * polynomial_space->n() == compute_n_pols(k),
+          ExcInternalError());
 }
 
 
@@ -135,14 +143,28 @@ PolynomialsABF<dim>::compute (const Point<dim>            &unit_point,
 
 template <int dim>
 unsigned int
-PolynomialsABF<dim>::compute_n_pols(unsigned int k)
+PolynomialsABF<dim>::compute_n_pols(const unsigned int k)
 {
-  if (dim == 1) return k+1;
-  if (dim == 2) return 2*(k+1)*(k+3);
-  //TODO:Check what are the correct numbers ...
-  if (dim == 3) return 3*(k+1)*(k+1)*(k+2);
+  switch (dim)
+    {
+    case 1:
+      // in 1d, we simply have Q_{k+2}, which has dimension k+3
+      return k+3;
+
+    case 2:
+      // the polynomial space is Q_{k+2,k} \times Q_{k,k+2}, which has
+      // 2(k+3)(k+1) DoFs
+      return 2*(k+3)*(k+1);
+
+    case 3:
+      // the polynomial space is Q_{k+2,k,k} \times Q_{k,k+2,k} \times Q_{k,k,k+2},
+      // which has 3(k+3)(k+1)(k+1) DoFs
+      return 3*(k+3)*(k+1)*(k+1);
+
+    default:
+      Assert(false, ExcNotImplemented());
+    }
 
-  Assert(false, ExcNotImplemented());
   return 0;
 }
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.