From: Wolfgang Bangerth Date: Mon, 20 Jul 1998 03:37:26 +0000 (+0000) Subject: Add some functionality to the finite element class and its siblings. X-Git-Tag: v8.0.0~22805 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=34ceab90910a160166ce7272e4cbd8fb1eba491b;p=dealii.git Add some functionality to the finite element class and its siblings. git-svn-id: https://svn.dealii.org/trunk@449 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe.h b/deal.II/deal.II/include/fe/fe.h index 87c4b073e3..0dcc0001f4 100644 --- a/deal.II/deal.II/include/fe/fe.h +++ b/deal.II/deal.II/include/fe/fe.h @@ -777,6 +777,24 @@ class FiniteElement : public FiniteElementBase { const vector > > &shape_grads_transform, const Boundary &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 > &unit_points) const; + /** * Compute the off-points of the finite * element basis functions on the given @@ -797,7 +815,7 @@ class FiniteElement : public FiniteElementBase { * #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 diff --git a/deal.II/deal.II/include/fe/fe_lib.lagrange.h b/deal.II/deal.II/include/fe/fe_lib.lagrange.h index 45777de5c8..608dd3345f 100644 --- a/deal.II/deal.II/include/fe/fe_lib.lagrange.h +++ b/deal.II/deal.II/include/fe/fe_lib.lagrange.h @@ -46,6 +46,12 @@ class FELinear : public FELinearMapping { virtual Point shape_grad(const unsigned int i, const Point& p) const; + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_unit_ansatz_points (vector > &ansatz_points) const; + /** * Refer to the base class for detailed * information on this function. @@ -101,6 +107,12 @@ class FEQuadraticSub : public FELinearMapping { virtual Point shape_grad(const unsigned int i, const Point& p) const; + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_unit_ansatz_points (vector > &ansatz_points) const; + /** * Refer to the base class for detailed * information on this function. @@ -175,6 +187,12 @@ class FECubicSub : public FELinearMapping { virtual Point shape_grad(const unsigned int i, const Point& p) const; + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_unit_ansatz_points (vector > &ansatz_points) const; + /** * Refer to the base class for detailed * information on this function. @@ -250,6 +268,12 @@ class FEQuarticSub : public FELinearMapping { virtual Point shape_grad(const unsigned int i, const Point& p) const; + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_unit_ansatz_points (vector > &ansatz_points) const; + /** * Refer to the base class for detailed * information on this function. diff --git a/deal.II/deal.II/source/fe/fe.cc b/deal.II/deal.II/source/fe/fe.cc index cc8e4a672e..ece34091c5 100644 --- a/deal.II/deal.II/source/fe/fe.cc +++ b/deal.II/deal.II/source/fe/fe.cc @@ -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 > &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(0); + for (unsigned int i=0; i(1); + + // now dofs on line + for (unsigned int i=0; i((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; ivertex(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::fill_fe_face_values (const DoFHandler::cell_iterat Assert (ansatz_points.size() == dofs_per_face, ExcWrongFieldDimension(ansatz_points.size(), dofs_per_face)); - static vector > dummy(total_dofs); + vector > dummy(total_dofs); fill_fe_values (cell, global_unit_points, jacobians, compute_jacobians, dummy, false, @@ -395,6 +417,13 @@ void FiniteElement::fill_fe_subface_values (const DoFHandler::cell_ite +template +void FiniteElement::get_unit_ansatz_points (vector > &) const { + Assert (false, ExcPureFunctionCalled()); +}; + + + template void FiniteElement::get_ansatz_points (const DoFHandler::cell_iterator &, const Boundary &, diff --git a/deal.II/deal.II/source/fe/fe_lib.cubic.cc b/deal.II/deal.II/source/fe/fe_lib.cubic.cc index 2263f03d14..b9179eaa76 100644 --- a/deal.II/deal.II/source/fe/fe_lib.cubic.cc +++ b/deal.II/deal.II/source/fe/fe_lib.cubic.cc @@ -378,6 +378,13 @@ FECubicSub<1>::shape_grad(const unsigned int i, +template <> +void FECubicSub<1>::get_unit_ansatz_points (vector > &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 > &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>&, diff --git a/deal.II/deal.II/source/fe/fe_lib.linear.cc b/deal.II/deal.II/source/fe/fe_lib.linear.cc index f293ae8dc9..3eb0afa921 100644 --- a/deal.II/deal.II/source/fe/fe_lib.linear.cc +++ b/deal.II/deal.II/source/fe/fe_lib.linear.cc @@ -79,6 +79,13 @@ FELinear<1>::shape_grad(const unsigned int i, +template <> +void FELinear<1>::get_unit_ansatz_points (vector > &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 > &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::get_ansatz_points (const typename DoFHandler::cell_iter const Boundary &, vector > &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::vertices_per_cell; ++vertex) ansatz_points[vertex] = cell->vertex(vertex); diff --git a/deal.II/deal.II/source/fe/fe_lib.quadratic.cc b/deal.II/deal.II/source/fe/fe_lib.quadratic.cc index 2f10e9d1a8..260711105c 100644 --- a/deal.II/deal.II/source/fe/fe_lib.quadratic.cc +++ b/deal.II/deal.II/source/fe/fe_lib.quadratic.cc @@ -98,6 +98,13 @@ FEQuadraticSub<1>::shape_grad(const unsigned int i, +template <> +void FEQuadraticSub<1>::get_unit_ansatz_points (vector > &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 > &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>&, diff --git a/deal.II/deal.II/source/fe/fe_lib.quartic.cc b/deal.II/deal.II/source/fe/fe_lib.quartic.cc index 71ad977592..9c55c2f8e7 100644 --- a/deal.II/deal.II/source/fe/fe_lib.quartic.cc +++ b/deal.II/deal.II/source/fe/fe_lib.quartic.cc @@ -391,6 +391,13 @@ FEQuarticSub<1>::shape_grad(const unsigned int i, +template <> +void FEQuarticSub<1>::get_unit_ansatz_points (vector > &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 > &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>&,