]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Neue FEValues
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Thu, 3 Sep 1998 15:44:42 +0000 (15:44 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Thu, 3 Sep 1998 15:44:42 +0000 (15:44 +0000)
git-svn-id: https://svn.dealii.org/trunk@551 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/include/fe/fe_values.h
deal.II/deal.II/source/fe/fe_values.cc
deal.II/deal.II/source/numerics/assembler.cc
deal.II/deal.II/source/numerics/error_estimator.cc
deal.II/deal.II/source/numerics/matrices.cc
deal.II/deal.II/source/numerics/vectors.cc

index aa211dc7aa9dd755d512598b53e86511488e319c..cbaa74528df4327383865b43e3869fe8fa8e013d 100644 (file)
@@ -6,6 +6,7 @@
 /*----------------------------   fe.h     ---------------------------*/
 
 #include <base/exceptions.h>
+#include <base/subscriptor.h>
 #include <grid/point.h>
 #include <grid/dof.h>
 #include <grid/geometry_info.h>
@@ -212,7 +213,11 @@ struct FiniteElementData<2> {
  * @author Wolfgang Bangerth, 1998
  */
 template <int dim>
-struct FiniteElementBase : public FiniteElementData<dim> {
+struct FiniteElementBase :
+  public Subscriptor,
+  public FiniteElementData<dim>
+
+{
   public:
                                     /**
                                      * Construct an object of this type.
index 43fa7ca5219c4bb0ecab4442df9ca5db82c7eb3e..3b3f53c6e97a63f518fb51a4651d9b7cd1a82bfe 100644 (file)
@@ -6,9 +6,10 @@
 /*----------------------------   fe_values.h     ---------------------------*/
 
 
+#include <base/exceptions.h>
+#include <base/subscriptor.h>
 #include <lac/dfmatrix.h>
 #include <grid/point.h>
-#include <base/exceptions.h>
 #include <grid/tria.h>
 #include <fe/fe_update_flags.h>
 
@@ -617,10 +618,13 @@ class FEValues : public FEValuesBase<dim> {
                                      * may use other ways.)
                                      */
     void reinit (const typename DoFHandler<dim>::cell_iterator &,
-                const FiniteElement<dim> &,
                 const Boundary<dim> &);
 
   private:
+                                    /**
+                                     * Store the finite element for later use.
+                                     */
+  SmartPointer<const FiniteElement<dim> > fe;
                                     /**
                                      * Store the gradients of the shape
                                      * functions at the quadrature points on
@@ -889,8 +893,12 @@ class FEFaceValues : public FEFaceValuesBase<dim> {
                                      */
     void reinit (const typename DoFHandler<dim>::cell_iterator &cell,
                 const unsigned int                    face_no,
-                const FiniteElement<dim>             &fe,
                 const Boundary<dim>                  &boundary);
+private:
+                                    /**
+                                     * Store the finite element for later use.
+                                     */
+  SmartPointer<const FiniteElement<dim> > fe;
 };
 
 
@@ -1041,13 +1049,17 @@ class FESubfaceValues : public FEFaceValuesBase<dim> {
     void reinit (const typename DoFHandler<dim>::cell_iterator &cell,
                 const unsigned int                    face_no,
                 const unsigned int                    subface_no,
-                const FiniteElement<dim>             &fe,
                 const Boundary<dim>                  &boundary);
 
                                     /**
                                      * Exception
                                      */
     DeclException0 (ExcReinitCalledWithBoundaryFace);
+private:
+                                    /**
+                                     * Store the finite element for later use.
+                                     */
+  SmartPointer<const FiniteElement<dim> > fe;
 };
 
 
index 647351d8002345ec40dc895e73868f701ce87823..27b8320022c04774b5e6fa83c845904a8bf0c032 100644 (file)
@@ -164,6 +164,7 @@ FEValues<dim>::FEValues (const FiniteElement<dim> &fe,
                                   fe.n_transform_functions,
                                   1,
                                   update_flags),
+               fe(&fe),
                unit_shape_gradients(fe.total_dofs,
                                     vector<Point<dim> >(quadrature.n_quadrature_points)),
                unit_shape_gradients_transform(fe.n_transform_functions,
@@ -201,7 +202,6 @@ FEValues<dim>::FEValues (const FiniteElement<dim> &fe,
 
 template <int dim>
 void FEValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell,
-                           const FiniteElement<dim>                      &fe,
                            const Boundary<dim>                           &boundary) {
   present_cell = cell;
                                   // fill jacobi matrices and real
@@ -211,7 +211,7 @@ void FEValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell,
       (update_flags & update_q_points)  ||
       (update_flags & update_gradients) ||
       (update_flags & update_ansatz_points))
-    fe.fill_fe_values (cell,
+    fe->fill_fe_values (cell,
                       unit_quadrature_points,
                       jacobi_matrices,
                       update_flags & (update_jacobians  |
@@ -227,7 +227,7 @@ void FEValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell,
                                   // compute gradients on real element if
                                   // requested
   if (update_flags & update_gradients) 
-    for (unsigned int i=0; i<fe.total_dofs; ++i)
+    for (unsigned int i=0; i<fe->total_dofs; ++i)
       for (unsigned int j=0; j<n_quadrature_points; ++j)
        for (unsigned int s=0; s<dim; ++s)
          {
@@ -314,7 +314,8 @@ FEFaceValues<dim>::FEFaceValues (const FiniteElement<dim> &fe,
                                       fe.total_dofs,
                                       fe.n_transform_functions,
                                       GeometryInfo<dim>::faces_per_cell,
-                                      update_flags)
+                                      update_flags),
+               fe(&fe)
 {
   unit_face_quadrature_points = quadrature.get_quad_points();
   weights = quadrature.get_weights ();  
@@ -353,7 +354,6 @@ FEFaceValues<dim>::FEFaceValues (const FiniteElement<dim> &fe,
 template <int dim>
 void FEFaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell,
                                const unsigned int                             face_no,
-                               const FiniteElement<dim>                      &fe,
                                const Boundary<dim>                           &boundary) {
   present_cell  = cell;
   selected_dataset = face_no;
@@ -365,7 +365,7 @@ void FEFaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &c
       (update_flags & update_gradients) ||
       (update_flags & update_ansatz_points) ||
       (update_flags & update_JxW_values))
-    fe.fill_fe_face_values (cell,
+    fe->fill_fe_face_values (cell,
                            face_no,
                            unit_face_quadrature_points,
                            unit_quadrature_points[face_no],
@@ -388,7 +388,7 @@ void FEFaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &c
                                   // compute gradients on real element if
                                   // requested
   if (update_flags & update_gradients) 
-    for (unsigned int i=0; i<fe.total_dofs; ++i)
+    for (unsigned int i=0; i<fe->total_dofs; ++i)
       {
        fill_n (shape_gradients[i].begin(),
                n_quadrature_points,
@@ -432,7 +432,8 @@ FESubfaceValues<dim>::FESubfaceValues (const FiniteElement<dim> &fe,
                                       fe.total_dofs,
                                       fe.n_transform_functions,
                                       GeometryInfo<dim>::faces_per_cell * GeometryInfo<dim>::subfaces_per_face,
-                                      update_flags)
+                                      update_flags),
+               fe(&fe)
 {
   Assert ((update_flags & update_ansatz_points) == false,
          ExcInvalidUpdateFlag());
@@ -479,7 +480,6 @@ template <int dim>
 void FESubfaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell,
                                   const unsigned int         face_no,
                                   const unsigned int         subface_no,
-                                  const FiniteElement<dim>  &fe,
                                   const Boundary<dim>       &boundary) {
   Assert (cell->face(face_no)->at_boundary() == false,
          ExcReinitCalledWithBoundaryFace());
@@ -493,7 +493,7 @@ void FESubfaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator
       (update_flags & update_q_points)  ||
       (update_flags & update_gradients) ||
       (update_flags & update_JxW_values))
-    fe.fill_fe_subface_values (cell,
+    fe->fill_fe_subface_values (cell,
                               face_no,
                               subface_no,
                               unit_face_quadrature_points,
@@ -515,7 +515,7 @@ void FESubfaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator
                                   // compute gradients on real element if
                                   // requested
   if (update_flags & update_gradients) 
-    for (unsigned int i=0; i<fe.total_dofs; ++i) 
+    for (unsigned int i=0; i<fe->total_dofs; ++i) 
       {
        fill_n (shape_gradients[i].begin(),
                n_quadrature_points,
index 4881ac452fa7e14723d8da22d81e8f574a7c4b08..1cec0d6cbeb0553173a0b45d08a91af20edd37f8 100644 (file)
@@ -106,7 +106,6 @@ void Assembler<dim>::assemble (const Equation<dim> &equation) {
                                                    present_level,
                                                    present_index,
                                                    dof_handler),
-                   fe,
                    boundary);
   const unsigned int n_dofs = dof_handler->get_selected_fe().total_dofs;
 
index 647482a5017f00418e1f37a08b70753cb3ef425e..23b0002265b15b7a79ca5ee35e86e82089d445cf 100644 (file)
@@ -224,7 +224,7 @@ template <int dim>
 void KellyErrorEstimator<dim>::
 integrate_over_regular_face (const active_cell_iterator &cell,
                             const unsigned int          face_no,
-                            const FiniteElement<dim>   &fe,
+                            const FiniteElement<dim>   &,
                             const Boundary<dim>        &boundary,
                             const FunctionMap          &neumann_bc,
                             const unsigned int          n_q_points,
@@ -237,7 +237,7 @@ integrate_over_regular_face (const active_cell_iterator &cell,
   
                                   // initialize data of the restriction
                                   // of this cell to the present face
-  fe_face_values_cell.reinit (cell, face_no, fe, boundary);
+  fe_face_values_cell.reinit (cell, face_no, boundary);
   
                                   // set up a vector of the gradients
                                   // of the finite element function
@@ -275,8 +275,7 @@ integrate_over_regular_face (const active_cell_iterator &cell,
                                       // get restriction of finite element
                                       // function of #neighbor# to the
                                       // common face.
-      fe_face_values_neighbor.reinit (neighbor, neighbor_neighbor,
-                                     fe, boundary);
+      fe_face_values_neighbor.reinit (neighbor, neighbor_neighbor, boundary);
 
                                       // get a list of the gradients of
                                       // the finite element solution
@@ -380,7 +379,7 @@ template <int dim>
 void KellyErrorEstimator<dim>::
 integrate_over_irregular_face (const active_cell_iterator &cell,
                               const unsigned int          face_no,
-                              const FiniteElement<dim>   &fe,
+                              const FiniteElement<dim>   &,
                               const Boundary<dim>        &boundary,
                               const unsigned int          n_q_points,
                               FEFaceValues<dim>          &fe_face_values,
@@ -429,16 +428,14 @@ integrate_over_irregular_face (const active_cell_iterator &cell,
                                       // present cell to the subface and
                                       // store the gradient of the solution
                                       // in psi
-      fe_subface_values.reinit (cell, face_no, subface_no,
-                               fe, boundary);
+      fe_subface_values.reinit (cell, face_no, subface_no, boundary);
       fe_subface_values.get_function_grads (solution, psi);
 
                                       // restrict the finite element on the
                                       // neighbor cell to the common #subface#.
                                       // store the gradient in #neighbor_psi#
       vector<Point<dim> > neighbor_psi (n_q_points);
-      fe_face_values.reinit (neighbor_child, neighbor_neighbor,
-                            fe, boundary);
+      fe_face_values.reinit (neighbor_child, neighbor_neighbor, boundary);
       fe_face_values.get_function_grads (solution, neighbor_psi);
       
                                       // compute the jump in the gradients
index 938988f5024ed18698f6f986fd6b8e229f51b71b..2089de405f542ab8716bd662d8343ab299a3ed4a 100644 (file)
@@ -165,7 +165,7 @@ void MatrixCreator<dim>::create_boundary_mass_matrix (const DoFHandler<dim>    &
          cell_matrix.clear ();
          cell_vector.clear ();
          
-         fe_values.reinit (cell, face, fe, boundary);
+         fe_values.reinit (cell, face, boundary);
 
          const dFMatrix       &values    = fe_values.get_shape_values ();
          const vector<double> &weights   = fe_values.get_JxW_values ();
index 30e780f6a54e2901563e08093db88fb2d7ac090d..38f13f1a1e660ba3a076a94eccf4728a68fbc2cf 100644 (file)
@@ -370,7 +370,7 @@ void VectorTools<dim>::integrate_difference (const DoFHandler<dim>    &dof,
     {
       double diff=0;
                                       // initialize for this cell
-      fe_values.reinit (cell, fe, boundary);
+      fe_values.reinit (cell, boundary);
 
       switch (norm) 
        {

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.