]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add some functionality to the finite element class and its siblings.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 20 Jul 1998 03:37:26 +0000 (03:37 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 20 Jul 1998 03:37:26 +0000 (03:37 +0000)
git-svn-id: https://svn.dealii.org/trunk@449 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/include/fe/fe_lib.lagrange.h
deal.II/deal.II/source/fe/fe.cc
deal.II/deal.II/source/fe/fe_lib.cubic.cc
deal.II/deal.II/source/fe/fe_lib.linear.cc
deal.II/deal.II/source/fe/fe_lib.quadratic.cc
deal.II/deal.II/source/fe/fe_lib.quartic.cc

index 87c4b073e38e4c83ba37068cffeb51087db902f7..0dcc0001f46c7c60f79c77363700dbfe6de69e85 100644 (file)
@@ -777,6 +777,24 @@ class FiniteElement : public FiniteElementBase<dim> {
                                         const vector<vector<Point<dim> > > &shape_grads_transform,
                                         const Boundary<dim> &boundary) const;
 
+                                    /**
+                                     * Return the ansatz points of the
+                                     * ansatz functions on the unit cell.
+                                     *
+                                     * The function assumes that the
+                                     * #unit_points# array already has the
+                                     * right size. The order of points in
+                                     * the array matches that returned by
+                                     * the #cell->get_dof_indices# function.
+                                     *
+                                     * For one space dimension there is a
+                                     * standard implementation assuming
+                                     * equidistant off-points on the unit
+                                     * line. For all other dimensions, an
+                                     * overwritten function has to be provided.
+                                     */
+    virtual void get_unit_ansatz_points (vector<Point<dim> > &unit_points) const;
+    
                                     /**
                                      * Compute the off-points of the finite
                                      * element basis functions on the given
@@ -797,7 +815,7 @@ class FiniteElement : public FiniteElementBase<dim> {
                                      * #ansatz_points# array already has the
                                      * right size. The order of points in
                                      * the array matches that returned by
-                                     * the #face->get_dof_indices# function.
+                                     * the #cell->get_dof_indices# function.
                                      *
                                      * For one space dimension there is a
                                      * standard implementation assuming
index 45777de5c8deecda254fec71aef55fddb8326d9a..608dd3345f852c535e9b632ca41232c3028d1507 100644 (file)
@@ -46,6 +46,12 @@ class FELinear : public FELinearMapping<dim> {
     virtual Point<dim> shape_grad(const unsigned int i,
                                  const Point<dim>& p) const;
 
+                                    /**
+                                     * Refer to the base class for detailed
+                                     * information on this function.
+                                     */
+    virtual void get_unit_ansatz_points (vector<Point<dim> > &ansatz_points) const;
+
                                     /**
                                      * Refer to the base class for detailed
                                      * information on this function.
@@ -101,6 +107,12 @@ class FEQuadraticSub : public FELinearMapping<dim> {
     virtual Point<dim> shape_grad(const unsigned int i,
                                  const Point<dim>& p) const;
 
+                                    /**
+                                     * Refer to the base class for detailed
+                                     * information on this function.
+                                     */
+    virtual void get_unit_ansatz_points (vector<Point<dim> > &ansatz_points) const;
+
                                     /**
                                      * Refer to the base class for detailed
                                      * information on this function.
@@ -175,6 +187,12 @@ class FECubicSub : public FELinearMapping<dim> {
     virtual Point<dim> shape_grad(const unsigned int i,
                                  const Point<dim>& p) const;
 
+                                    /**
+                                     * Refer to the base class for detailed
+                                     * information on this function.
+                                     */
+    virtual void get_unit_ansatz_points (vector<Point<dim> > &ansatz_points) const;
+
                                     /**
                                      * Refer to the base class for detailed
                                      * information on this function.
@@ -250,6 +268,12 @@ class FEQuarticSub : public FELinearMapping<dim> {
     virtual Point<dim> shape_grad(const unsigned int i,
                                  const Point<dim>& p) const;
 
+                                    /**
+                                     * Refer to the base class for detailed
+                                     * information on this function.
+                                     */
+    virtual void get_unit_ansatz_points (vector<Point<dim> > &ansatz_points) const;
+
                                     /**
                                      * Refer to the base class for detailed
                                      * information on this function.
index cc8e4a672edd8e08c003e3958a72362a1a6e7cfb..ece34091c5f15c95bf360a826a233bf57b2520e9 100644 (file)
@@ -257,6 +257,28 @@ void FiniteElement<1>::fill_fe_subface_values (const DoFHandler<1>::cell_iterato
 
 
 
+template <>
+void FiniteElement<1>::get_unit_ansatz_points (vector<Point<1> > &unit_points) const {
+  Assert (unit_points.size() == total_dofs,
+         ExcWrongFieldDimension(unit_points.size(), total_dofs));
+                                  // compute ansatz points. The first ones
+                                  // belong to vertex one, the second ones
+                                  // to vertex two, all following are
+                                  // equally spaced along the line
+  unsigned int next = 0;
+                                  // first the dofs in the vertices
+  for (unsigned int i=0; i<dofs_per_vertex; ++i)
+    ansatz_points[next++] = Point<1>(0);
+  for (unsigned int i=0; i<dofs_per_vertex; ++i)
+    ansatz_points[next++] = Point<1>(1);
+  
+                                  // now dofs on line
+  for (unsigned int i=0; i<dofs_per_line; ++i) 
+    ansatz_points[next++] = Point<1>((i+1.0)/(dofs_per_line+1.0));
+};
+
+
+
 template <>
 void FiniteElement<1>::get_ansatz_points (const DoFHandler<1>::cell_iterator &cell,
                                          const Boundary<1> &,
@@ -278,7 +300,7 @@ void FiniteElement<1>::get_ansatz_points (const DoFHandler<1>::cell_iterator &ce
                                   // now dofs on line
   for (unsigned int i=0; i<dofs_per_line; ++i) 
     ansatz_points[next++] = cell->vertex(0) +
-                           Point<1>((i+1.0)/(total_dofs+1.0)*h);
+                           Point<1>((i+1.0)/(dofs_per_line+1.0)*h);
 };
 
 #endif
@@ -329,7 +351,7 @@ void FiniteElement<dim>::fill_fe_face_values (const DoFHandler<dim>::cell_iterat
   Assert (ansatz_points.size() == dofs_per_face,
          ExcWrongFieldDimension(ansatz_points.size(), dofs_per_face));
   
-  static vector<Point<dim> > dummy(total_dofs);
+  vector<Point<dim> > dummy(total_dofs);
   fill_fe_values (cell, global_unit_points,
                  jacobians, compute_jacobians,
                  dummy, false,
@@ -395,6 +417,13 @@ void FiniteElement<dim>::fill_fe_subface_values (const DoFHandler<dim>::cell_ite
 
 
 
+template <int dim>
+void FiniteElement<dim>::get_unit_ansatz_points (vector<Point<dim> > &) const {
+  Assert (false, ExcPureFunctionCalled());
+};
+
+
+
 template <int dim>
 void FiniteElement<dim>::get_ansatz_points (const DoFHandler<dim>::cell_iterator &,
                                            const Boundary<dim> &,
index 2263f03d14d3c101e42d563aefbd5a64849488a9..b9179eaa76c414ac77f64fc6f8350e64d001fd19 100644 (file)
@@ -378,6 +378,13 @@ FECubicSub<1>::shape_grad(const unsigned int i,
 
 
 
+template <>
+void FECubicSub<1>::get_unit_ansatz_points (vector<Point<1> > &unit_points) const {
+  FiniteElement<1>::get_unit_ansatz_points (unit_points);
+};
+
+
+
 template <>
 void FECubicSub<1>::get_ansatz_points (const typename DoFHandler<1>::cell_iterator &cell,
                                           const Boundary<1>  &boundary,
@@ -1539,6 +1546,31 @@ void FECubicSub<2>::get_local_mass_matrix (const DoFHandler<2>::cell_iterator &c
 
 
 
+template <>
+void FECubicSub<2>::get_unit_ansatz_points (vector<Point<2> > &unit_points) const {
+  Assert (unit_points.size() == total_dofs,
+         ExcWrongFieldDimension (unit_points.size(), total_dofs));
+
+  unit_points[0] = Point<2>(0,0);
+  unit_points[1] = Point<2>(1,0);
+  unit_points[2] = Point<2>(1,1);
+  unit_points[3] = Point<2>(0,1);
+  unit_points[4] = Point<2>(1./3,0);
+  unit_points[5] = Point<2>(2./3,0);
+  unit_points[6] = Point<2>(1,1./3);
+  unit_points[7] = Point<2>(1,2./3);
+  unit_points[8] = Point<2>(1./3,1);
+  unit_points[9] = Point<2>(2./3,1);
+  unit_points[10]= Point<2>(0,1./3);
+  unit_points[11]= Point<2>(0,2./3);
+  unit_points[12]= Point<2>(1./3,1./3);
+  unit_points[13]= Point<2>(2./3,1./3);
+  unit_points[14]= Point<2>(2./3,2./3);
+  unit_points[15]= Point<2>(1./3,2./3);
+};
+
+
+
 template <>
 void FECubicSub<2>::get_ansatz_points (const typename DoFHandler<2>::cell_iterator &cell,
                                           const Boundary<2>&,
index f293ae8dc9b899772466b529f68b3d4a3c138243..3eb0afa921834c79fa3b9299d77816ff0bf38213 100644 (file)
@@ -79,6 +79,13 @@ FELinear<1>::shape_grad(const unsigned int i,
 
 
 
+template <>
+void FELinear<1>::get_unit_ansatz_points (vector<Point<1> >  &ansatz_points) const {
+  FiniteElement<1>::get_unit_ansatz_points (ansatz_points);
+};
+
+
+
 template <>
 void FELinear<1>::get_ansatz_points (const typename DoFHandler<1>::cell_iterator &cell,
                                     const Boundary<1>  &boundary,
@@ -352,6 +359,18 @@ void FELinear<2>::get_local_mass_matrix (const DoFHandler<2>::cell_iterator &cel
 
 
 
+template <>
+void FELinear<2>::get_unit_ansatz_points (vector<Point<2> > &unit_points) const {
+  Assert (unit_points.size() == total_dofs,
+         ExcWrongFieldDimension (unit_points.size(), total_dofs));
+
+  unit_points[0] = Point<2> (0,0);
+  unit_points[1] = Point<2> (1,0);
+  unit_points[2] = Point<2> (1,1);
+  unit_points[3] = Point<2> (0,1);
+};
+  
+
 #endif
 
 
@@ -361,7 +380,7 @@ void FELinear<dim>::get_ansatz_points (const typename DoFHandler<dim>::cell_iter
                                       const Boundary<dim>  &,
                                       vector<Point<dim> >  &ansatz_points) const {
   Assert (ansatz_points.size() == total_dofs,
-         ExcWrongFieldDimension (ansatz_points.size(), 1<<(dim-1)));
+         ExcWrongFieldDimension (ansatz_points.size(), total_dofs));
   
   for (unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell; ++vertex)
     ansatz_points[vertex] = cell->vertex(vertex);
index 2f10e9d1a8b7a968c17bce6be83e62224b35f090..260711105c9a64d5fe3d4c8789cb0daead37585e 100644 (file)
@@ -98,6 +98,13 @@ FEQuadraticSub<1>::shape_grad(const unsigned int i,
 
 
 
+template <>
+void FEQuadraticSub<1>::get_unit_ansatz_points (vector<Point<1> > &unit_points) const {
+  FiniteElement<1>::get_unit_ansatz_points (unit_points);
+};
+
+
+
 template <>
 void FEQuadraticSub<1>::get_ansatz_points (const typename DoFHandler<1>::cell_iterator &cell,
                                           const Boundary<1>  &boundary,
@@ -856,6 +863,24 @@ void FEQuadraticSub<2>::get_local_mass_matrix (const DoFHandler<2>::cell_iterato
 
 
 
+template <>
+void FEQuadraticSub<2>::get_unit_ansatz_points (vector<Point<2> > &unit_points) const {
+  Assert (unit_points.size() == total_dofs,
+         ExcWrongFieldDimension (unit_points.size(), total_dofs));
+  
+  unit_points[0] = Point<2> (0,0);
+  unit_points[1] = Point<2> (1,0);
+  unit_points[2] = Point<2> (1,1);
+  unit_points[3] = Point<2> (0,1);
+  unit_points[4] = Point<2> (0.5,0);
+  unit_points[5] = Point<2> (1,0.5);
+  unit_points[6] = Point<2> (0.5,1);
+  unit_points[7] = Point<2> (0,0.5);
+  unit_points[8] = Point<2> (0.5,0.5);
+};
+
+  
+
 template <>
 void FEQuadraticSub<2>::get_ansatz_points (const typename DoFHandler<2>::cell_iterator &cell,
                                           const Boundary<2>&,
index 71ad977592173ecc3e30c4ba2532ce9f52b810db..9c55c2f8e742ba2a09d5f5adfe9ae7489b9d1f9b 100644 (file)
@@ -391,6 +391,13 @@ FEQuarticSub<1>::shape_grad(const unsigned int i,
 
 
 
+template <>
+void FEQuarticSub<1>::get_unit_ansatz_points (vector<Point<1> > &unit_points) const {
+  FiniteElement<1>::get_unit_ansatz_points (unit_points);
+};
+
+
+
 template <>
 void FEQuarticSub<1>::get_ansatz_points (const typename DoFHandler<1>::cell_iterator &cell,
                                         const Boundary<1>  &boundary,
@@ -2652,6 +2659,40 @@ void FEQuarticSub<2>::get_local_mass_matrix (const DoFHandler<2>::cell_iterator
 
 
 
+template <>
+void FEQuarticSub<2>::get_unit_ansatz_points (vector<Point<2> > &unit_points) const {
+  Assert (unit_points.size() == total_dofs,
+         ExcWrongFieldDimension (unit_points.size(), total_dofs));
+
+  unit_points[0] = Point<2>(0,0);
+  unit_points[1] = Point<2>(1,0);
+  unit_points[2] = Point<2>(1,1);
+  unit_points[3] = Point<2>(0,1);
+  unit_points[4] = Point<2>(1./4,0);
+  unit_points[5] = Point<2>(2./4,0);
+  unit_points[6] = Point<2>(3./4,0);
+  unit_points[7] = Point<2>(1,1./4);
+  unit_points[8] = Point<2>(1,2./4);
+  unit_points[9] = Point<2>(1,3./4);
+  unit_points[10] = Point<2>(1./4,1);
+  unit_points[11] = Point<2>(2./4,1);
+  unit_points[12] = Point<2>(3./4,1);
+  unit_points[13] = Point<2>(0,1./4);
+  unit_points[14] = Point<2>(0,2./4);
+  unit_points[15] = Point<2>(0,3./4);
+  unit_points[16] = Point<2>(1./4,1./4);
+  unit_points[17] = Point<2>(3./4,1./4);
+  unit_points[18] = Point<2>(3./4,3./4);
+  unit_points[19] = Point<2>(1./4,3./4);
+  unit_points[20] = Point<2>(1./2,1./4);
+  unit_points[21] = Point<2>(3./4,1./2);
+  unit_points[22] = Point<2>(1./2,3./4);
+  unit_points[23] = Point<2>(1./4,1./2);
+  unit_points[24] = Point<2>(1./2,1./2);
+};
+
+
+
 template <>
 void FEQuarticSub<2>::get_ansatz_points (const typename DoFHandler<2>::cell_iterator &cell,
                                         const Boundary<2>&,

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.