From: hartmann Date: Wed, 16 Dec 1998 15:09:30 +0000 (+0000) Subject: Implementation of the discontinuous elements, dg0 - dg4 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=424729e1b13adbc6bbf87bcd71b31fb00a327748;p=dealii-svn.git Implementation of the discontinuous elements, dg0 - dg4 git-svn-id: https://svn.dealii.org/trunk@715 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/source/fe/fe_lib.dg.cc b/deal.II/deal.II/source/fe/fe_lib.dg.cc new file mode 100644 index 0000000000..0c40e6599b --- /dev/null +++ b/deal.II/deal.II/source/fe/fe_lib.dg.cc @@ -0,0 +1,53 @@ +/* $Id$ */ +/* Ralf Hartmann, University of Heidelberg, Dez 98 */ + +#include + +template +FEDGLinear::FEDGLinear(): + FELinear(1) {}; + +template +FEDGQuadraticSub::FEDGQuadraticSub(): + FEQuadraticSub(1) {}; + +template +FEDGCubicSub::FEDGCubicSub(): + FECubicSub(1) {}; + +template +FEDGQuarticSub::FEDGQuarticSub(): + FEQuarticSub(1) {}; + + + +template +const dFMatrix & +FEDGLinear::restrict (const unsigned int child) const { + Assert (false, ExcNotImplemented()); + return restriction[child]; +}; + + +template +const dFMatrix & +FEDGQuadraticSub::restrict (const unsigned int child) const { + Assert (false, ExcNotImplemented()); + return restriction[child]; +}; + + +template +const dFMatrix & +FEDGCubicSub::restrict (const unsigned int child) const { + Assert (false, ExcNotImplemented()); + return restriction[child]; +}; + + +template +const dFMatrix & +FEDGQuarticSub::restrict (const unsigned int child) const { + Assert (false, ExcNotImplemented()); + return restriction[child]; +}; diff --git a/deal.II/deal.II/source/fe/fe_lib.dg.constant.cc b/deal.II/deal.II/source/fe/fe_lib.dg.constant.cc new file mode 100644 index 0000000000..2e8637edb0 --- /dev/null +++ b/deal.II/deal.II/source/fe/fe_lib.dg.constant.cc @@ -0,0 +1,219 @@ +/* Ralf Hartmann, University of Heidelberg, Dez 98 */ + + +#include +#include +#include +#include +#include + + + + +#if deal_II_dimension == 1 + +template <> +FEDGConstant<1>::FEDGConstant () : + FELinearMapping<1> (0, 1) +{ + // for restriction and prolongation matrices: + // note that we do not add up all the + // contributions since then we would get + // two summands per vertex in 1d (four + // in 2d, etc), but only one per line dof. + // We could accomplish for that by dividing + // the vertex dof values by 2 (4, etc), but + // would get into trouble at the boundary + // of the domain since there only one + // cell contributes to a vertex. Rather, + // we do not add up the contributions but + // set them right into the matrices! + + // The restriction matrizes got crazy values + // as it is yet not clear how they should work + // in the DG(0) case. In general + // the use of the restriction matrices + // is not yet finally decided about, too. + restriction[0](0,0) = 1e8; + restriction[1](0,0) = 1e8; + + prolongation[0](0,0) = 1.0; + + prolongation[1](0,0) = 1.0; +}; + + +template <> +void FEDGConstant<1>::get_face_support_points (const typename DoFHandler<1>::face_iterator &, + const Boundary<1> &, + vector > &) const { + Assert (false, ExcInternalError()); +}; + + +/* +template <> +void FEDGConstant<1>::get_local_mass_matrix (const DoFHandler<1>::cell_iterator &cell, + const Boundary<1> &, + dFMatrix &local_mass_matrix) const { + Assert (local_mass_matrix.n() == total_dofs, + ExcWrongFieldDimension(local_mass_matrix.n(),total_dofs)); + Assert (local_mass_matrix.m() == total_dofs, + ExcWrongFieldDimension(local_mass_matrix.m(),total_dofs)); + + const double h = cell->vertex(1)(0) - cell->vertex(0)(0); + Assert (h>0, ExcJacobiDeterminantHasWrongSign()); + + local_mass_matrix(0,0) = h; +}; +*/ +#endif + + + + +#if deal_II_dimension == 2 + +template <> +FEDGConstant<2>::FEDGConstant () : + FELinearMapping<2> (0, 0, 1) +{ + // The restriction matrizes got crazy values + // as it is yet not clear how they should work + // in the DG(0) case. In general + // the use of the restriction matrices + // is not yet finally decided about, too. + restriction[0](0,0) = 1e8; + restriction[1](0,0) = 1e8; + restriction[2](0,0) = 1e8; + restriction[3](0,0) = 1e8; + + prolongation[0](0,0) = 1.0; + + prolongation[1](0,0) = 1.0; + + prolongation[2](0,0) = 1.0; + + prolongation[3](0,0) = 1.0; +}; + +#endif + + + + + +template +inline +double +FEDGConstant::shape_value (const unsigned int i, + const Point&) const +{ + Assert((i +inline +Tensor<1,dim> +FEDGConstant::shape_grad (const unsigned int i, + const Point&) const +{ + Assert((i, so we + // still construct it as that. it should + // make no difference in practice, + // however + return Point (); +}; + + + +template +inline +Tensor<2,dim> +FEDGConstant::shape_grad_grad (const unsigned int i, + const Point &) const +{ + Assert((i(); +}; + + + +template +void FEDGConstant::get_local_mass_matrix (const DoFHandler::cell_iterator &cell, + const Boundary &, + dFMatrix &local_mass_matrix) const { + Assert (local_mass_matrix.n() == total_dofs, + ExcWrongFieldDimension(local_mass_matrix.n(),total_dofs)); + Assert (local_mass_matrix.m() == total_dofs, + ExcWrongFieldDimension(local_mass_matrix.m(),total_dofs)); + + local_mass_matrix(0,0) = cell->measure(); +}; + + + +template +void +FEDGConstant::get_unit_support_points (vector > &unit_points) const { + Assert (unit_points.size() == total_dofs, + ExcWrongFieldDimension (unit_points.size(), total_dofs)); + for (unsigned int d=0; d +void +FEDGConstant::get_support_points (const typename DoFHandler::cell_iterator &cell, + const Boundary &, + vector > &support_points) const { + Assert (support_points.size() == total_dofs, + ExcWrongFieldDimension (support_points.size(), total_dofs)); + + support_points[0] = cell->center(); +}; + + + +template +void +FEDGConstant::get_face_support_points (const typename DoFHandler::face_iterator &, + const Boundary &, + vector > &support_points) const { + Assert (false, ExcNotImplemented()); + + Assert ((support_points.size() == dofs_per_face) && + (support_points.size() == GeometryInfo::vertices_per_face), + ExcWrongFieldDimension (support_points.size(), + GeometryInfo::vertices_per_face)); +}; + + +template +const dFMatrix & +FEDGConstant::restrict (const unsigned int child) const { + Assert (false, ExcNotImplemented()); + return restriction[child]; +}; + + + +// explicit instantiations + +template class FEDGConstant; + + + + + + + + +