]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Throw an exception in 3d if degree>3 since then the mapping is not
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 30 Oct 2003 16:30:48 +0000 (16:30 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 30 Oct 2003 16:30:48 +0000 (16:30 +0000)
implemented. Previously, we would allow the creation of such an
object, but that didn't make much sense since the object couldn't be
used.

git-svn-id: https://svn.dealii.org/trunk@8176 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/source/fe/mapping_q.cc

index 848af86bd5b2b70fbffda399be83a8abb1cd1baf..39b945af7785f06ff70c6fd663ae3f23ed8b2183 100644 (file)
@@ -475,15 +475,33 @@ MappingQ<dim>::set_laplace_on_quad_vector(Table<2,double> &loqvs) const
                                       // not precomputed, then do so now
       if (dim==2)
        compute_laplace_vector(loqvs);
-      
-                                      // for dim==3 don't throw an
-                                      // ExcNotImplemented here to
-                                      // allow the creating of that
-                                      // MappingQ<3> object. But an
-                                      // ExcLaplaceVectorNotSet
-                                      // assertion is thrown when the
-                                      // apply_laplace_vector
-                                      // function is called.
+      else
+                                        // computing the Laplace
+                                        // vector for faces is not
+                                        // supported in 3d at
+                                        // present. presumably, doing
+                                        // so would not be so hard:
+                                        // we would only have to call
+                                        // the function in 2d,
+                                        // i.e. the quad values in 3d
+                                        // are equal to the cell
+                                        // values in 2d. however,
+                                        // that would require us to
+                                        // link in the 2d library,
+                                        // which is kind of awkward
+                                        // (note that
+                                        // compute_laplace_vector
+                                        // really makes use of a lot
+                                        // of 2d stuff, such as
+                                        // FEValues etc). an
+                                        // alternative would be to
+                                        // precompute the values of
+                                        // this array for a couple of
+                                        // higher mapping orders, pin
+                                        // down their values and
+                                        // insert them into the array
+                                        // above.
+       Assert (false, ExcNotImplemented());
     }
 
                                   // the sum of weights of the points
@@ -592,9 +610,9 @@ MappingQ<dim>::compute_laplace_vector(Table<2,double> &lvs) const
   for (unsigned int point=0; point<n_q_points; ++point)
     for (unsigned int i=0; i<n_inner; ++i)
       for (unsigned int j=0; j<n_inner; ++j)
-       S(i,j)+=contract(quadrature_data.derivative(point, n_outer+i),
-                        quadrature_data.derivative(point, n_outer+j))
-               *quadrature.weight(point);
+       S(i,j) += contract(quadrature_data.derivative(point, n_outer+i),
+                          quadrature_data.derivative(point, n_outer+j))
+                 * quadrature.weight(point);
   
                                   // Compute the components of T to be the
                                   // product of gradients of inner and
@@ -603,9 +621,9 @@ MappingQ<dim>::compute_laplace_vector(Table<2,double> &lvs) const
   for (unsigned int point=0; point<n_q_points; ++point)
     for (unsigned int i=0; i<n_inner; ++i)
       for (unsigned int k=0; k<n_outer; ++k)
-       T(i,k)+=contract(quadrature_data.derivative(point, n_outer+i),
-                        quadrature_data.derivative(point, k))
-               *quadrature.weight(point);
+       T(i,k) += contract(quadrature_data.derivative(point, n_outer+i),
+                          quadrature_data.derivative(point, k))
+                 *quadrature.weight(point);
   
   FullMatrix<double> S_1(n_inner);
   S_1.invert(S);

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.