From: Wolfgang Bangerth Date: Fri, 6 Mar 1998 17:33:13 +0000 (+0000) Subject: Small changes and interface extension. X-Git-Tag: v8.0.0~23208 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=fcf91aa0daa58b371d1f606234338617aca1c971;p=dealii.git Small changes and interface extension. git-svn-id: https://svn.dealii.org/trunk@43 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/base/include/base/quadrature.h b/deal.II/base/include/base/quadrature.h index fbe09067ed..7f06f7c545 100644 --- a/deal.II/base/include/base/quadrature.h +++ b/deal.II/base/include/base/quadrature.h @@ -11,7 +11,11 @@ /** - Base class for quadrature formulae in arbitrary dimensions. + Base class for quadrature formulae in arbitrary dimensions. This class + stores quadrature points and weights on the unit line [0,1], unit + square [0,1]x[0,1], etc. This information is used together with + objects of the \Ref{FiniteElement} class to compute the values stored + in the \Ref{FEValues} objects. */ template class Quadrature { @@ -19,12 +23,12 @@ class Quadrature { /** * Number of quadrature points. */ - const unsigned int n_quad_points; + const unsigned int n_quadrature_points; /** * Constructor. */ - Quadrature (const unsigned int n_quad_points); + Quadrature (const unsigned int n_quadrature_points); /** * Return the #i#th quadrature point. diff --git a/deal.II/base/source/quadrature_lib.cc b/deal.II/base/source/quadrature_lib.cc index 92dd78e099..0f21f502fb 100644 --- a/deal.II/base/source/quadrature_lib.cc +++ b/deal.II/base/source/quadrature_lib.cc @@ -10,7 +10,7 @@ QGauss2<1>::QGauss2 () : static const double xpts[] = { 0.288675135, 0.71132486 }; static const double wts[] = { 0.5, 0.5 }; - for (unsigned int i=0; i (xpts[i])); weights.push_back (wts[i]); @@ -34,7 +34,7 @@ QGauss2x4<1>::QGauss2x4 () : static const double wts[] = { 0.5*W0, 0.5*W1, 0.5*W1, 0.5*W0, 0.5*W0, 0.5*W1, 0.5*W1, 0.5*W0 }; - for (unsigned int i=0; i (xpts[i])); weights.push_back (wts[i]); @@ -56,7 +56,7 @@ QGauss4<1>::QGauss4 () : static const double xpts[] = { G0, G1, G2, G3 }; static const double wts[] = { W0, W1, W1, W0 }; - for (unsigned int i=0; i (xpts[i])); weights.push_back (wts[i]); @@ -84,7 +84,7 @@ QGauss8<1>::QGauss8 () : static const double xpts[] = { G0, G1, G2, G3, G4, G5, G6, G7 }; static const double wts[] = { W0, W1, W2, W3, W3, W2, W1, W0 }; - for (unsigned int i=0; i (xpts[i])); weights.push_back (wts[i]); @@ -108,7 +108,7 @@ QSimpson<1>::QSimpson () : static const double xpts[] = { 0.0, 0.5, 1.0 }; static const double wts[] = { 1./6., 2./3., 1./6. }; - for (unsigned int i=0; i (xpts[i])); weights.push_back (wts[i]); @@ -123,7 +123,7 @@ QTrapez<1>::QTrapez () : static const double xpts[] = { 0.0, 1.0 }; static const double wts[] = { 0.5, 0.5 }; - for (unsigned int i=0; i (xpts[i])); weights.push_back (wts[i]); diff --git a/deal.II/deal.II/include/fe/fe.h b/deal.II/deal.II/include/fe/fe.h index 5006d86763..ca758eb7e6 100644 --- a/deal.II/deal.II/include/fe/fe.h +++ b/deal.II/deal.II/include/fe/fe.h @@ -21,8 +21,10 @@ template class Quadrature; This class is an optimization which avoids evaluating the shape functions at the quadrature points each time a quadrature takes place. Rather, the values and gradients (and possibly higher order derivatives in future - versions of this library) are evaluated once and for all before doing the - quadrature itself. + versions of this library) are evaluated once and for all on the unit + cell before doing the quadrature itself. Only the Jacobian matrix of + the transformation from the unit cell to the real cell and the integration + points in real space are calculated each time we move on to a new cell. Objects of this class store a multitude of different values needed to do the assemblage steps on real cells rather than on the unit cell. Among @@ -43,6 +45,16 @@ template class Quadrature; template class FEValues { public: + /** + * Number of quadrature points. + */ + const unsigned int n_quadrature_points; + + /** + * Total number of shape functions. + */ + const unsigned int total_dofs; + /** * Constructor. Fill all arrays with the * values of the shape functions of the @@ -76,6 +88,19 @@ class FEValues { const Point & shape_grad (const unsigned int i, const unsigned int j) const; + /** + * Return the position of the #i#th + * quadrature point in real space. + */ + const Point & quadrature_point (const unsigned int i) const; + + /** + * Return the Jacobi determinant times + * the weight of the #i#th quadrature + * point. + */ + double JxW (const unsigned int i) const; + /** * Reinitialize the gradients, Jacobi * determinants, etc for the given cell @@ -394,7 +419,7 @@ class FiniteElementBase { /** - Define a finite element type. + Define a finite element class. */ template class FiniteElement; diff --git a/deal.II/deal.II/source/fe/fe.cc b/deal.II/deal.II/source/fe/fe.cc index 79db61b61c..6d70b92e3b 100644 --- a/deal.II/deal.II/source/fe/fe.cc +++ b/deal.II/deal.II/source/fe/fe.cc @@ -15,26 +15,28 @@ extern TriaIterator<1,CellAccessor<1> > __dummy2687; // for gcc2.8 template FEValues::FEValues (const FiniteElement &fe, const Quadrature &quadrature) : - shape_values(fe.total_dofs, quadrature.n_quad_points), + n_quadrature_points(quadrature.n_quadrature_points), + total_dofs(fe.total_dofs), + shape_values(fe.total_dofs, quadrature.n_quadrature_points), shape_gradients(fe.total_dofs, - vector >(quadrature.n_quad_points)), + vector >(quadrature.n_quadrature_points)), unit_shape_gradients(fe.total_dofs, - vector >(quadrature.n_quad_points)), - weights(quadrature.n_quad_points, 0), - JxW_values(quadrature.n_quad_points, 0), - quadrature_points(quadrature.n_quad_points, Point()), - unit_quadrature_points(quadrature.n_quad_points, Point()), - jacobi_matrices (quadrature.n_quad_points) + vector >(quadrature.n_quadrature_points)), + weights(quadrature.n_quadrature_points, 0), + JxW_values(quadrature.n_quadrature_points, 0), + quadrature_points(quadrature.n_quadrature_points, Point()), + unit_quadrature_points(quadrature.n_quadrature_points, Point()), + jacobi_matrices (quadrature.n_quadrature_points) { for (unsigned int i=0; i::shape_grad (const unsigned int i, +template +const Point & FEValues::quadrature_point (const unsigned int i) const { + Assert (i +double FEValues::JxW (const unsigned int i) const { + Assert (i::reinit (const Triangulation<1>::cell_iterator &cell, const FiniteElement<1> &fe) { const unsigned int dim=1; @@ -77,7 +97,7 @@ void FEValues<1>::reinit (const Triangulation<1>::cell_iterator &cell, quadrature_points); // compute gradients on real element for (unsigned int i=0; i::reinit (const Triangulation<1>::cell_iterator &cell, // refer to the general doc for // why we take the inverse of the // determinant - for (unsigned int i=0; i::reinit (const Triangulation<2>::cell_iterator &cell, quadrature_points); // compute gradients on real element for (unsigned int i=0; i::reinit (const Triangulation<2>::cell_iterator &cell, // refer to the general doc for // why we take the inverse of the // determinant - for (unsigned int i=0; i::fill_fe_values (const Triangulation<1>::cell_iterator &ce // local mesh width double h=(cell->vertex(1)(0) - cell->vertex(0)(0)); - unsigned int n_points = unit_points.size(); - for (unsigned int i=0; ivertex(0) + h*unit_points[i];