]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Get rid of the weird SelectFEValues class since we can now always use hp::FEValues...
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 24 Feb 2006 16:43:23 +0000 (16:43 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 24 Feb 2006 16:43:23 +0000 (16:43 +0000)
git-svn-id: https://svn.dealii.org/trunk@12497 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/select_fe_values.h [deleted file]
deal.II/deal.II/source/numerics/data_out.cc
deal.II/deal.II/source/numerics/data_out_faces.cc
deal.II/deal.II/source/numerics/data_out_rotation.cc
deal.II/deal.II/source/numerics/data_out_stack.cc

diff --git a/deal.II/deal.II/include/fe/select_fe_values.h b/deal.II/deal.II/include/fe/select_fe_values.h
deleted file mode 100644 (file)
index 42c48fe..0000000
+++ /dev/null
@@ -1,69 +0,0 @@
-//----------------------------  select_fe_values.h  ---------------------------
-//    $Id$
-//    Version: $Name$
-//
-//    Copyright (C) 2003, 2006 by the deal.II authors
-//
-//    This file is subject to QPL and may not be  distributed
-//    without copyright and license information. Please refer
-//    to the file deal.II/doc/license.html for the  text  and
-//    further information on this license.
-//
-//----------------------------  select_fe_values.h  ---------------------------
-#ifndef __deal2__select_fe_values_h
-#define __deal2__select_fe_values_h
-
-#include <base/config.h>
-
-template <int> class DoFHandler;
-namespace hp
-{
-  template <int> class DoFHandler;
-}
-
-template <int> class FEValues;
-template <int> class FEFaceValues;
-template <int> class FESubfaceValues;
-
-namespace hp
-{
-  template <int> class FEValues;
-  template <int> class FEFaceValues;
-  template <int> class FESubfaceValues;
-}
-
-
-
-/**
- * 
- * @ingroup hp
- */  
-template <typename> struct SelectFEValues;
-
-/**
- * 
- * @ingroup hp
- */  
-template <int dim> struct SelectFEValues<DoFHandler<dim> >
-{
-    typedef FEValues<dim>        FEValues;
-    typedef FEFaceValues<dim>    FEFaceValues;
-    typedef FESubfaceValues<dim> FESubfaceValues;
-};
-
-
-/**
- * 
- * @ingroup hp
- */  
-template <int dim> struct SelectFEValues<hp::DoFHandler<dim> >
-{
-    typedef hp::FEValues<dim>        FEValues;
-    typedef hp::FEFaceValues<dim>    FEFaceValues;
-    typedef hp::FESubfaceValues<dim> FESubfaceValues;
-};
-
-
-
-
-#endif
index e8e1ced268b76adc97cf0c298570a36e1a205b8c..40b1ff5033a82b2fd48a2f1bbcd73caf6f048b07 100644 (file)
@@ -29,7 +29,6 @@
 #include <fe/fe_values.h>
 #include <fe/hp_fe_values.h>
 #include <fe/mapping_q1.h>
-#include <fe/select_fe_values.h>
 
 #include <sstream>
 
@@ -390,10 +389,24 @@ void DataOut<dim,DH>::build_some_patches (Data &data)
 
 //TODO[?]: This is strange -- Data has a member 'mapping' that should
 //be used here, but it isn't. Rather, up until version 1.94, we were
-//actually initializing a local mapping object and used that...
-  typename SelectFEValues<DH<dim> >::FEValues
-    x_fe_patch_values (this->dofs->get_fe(),
-                       patch_points, update_values);
+//actually initializing a local mapping object and used that... While
+//we use the mapping to transform the vertex coordinates, we do not
+//use the mapping to transform the shape functions (necessary for
+//example for Raviart-Thomas elements). This could lead to trouble
+//when someone tries to use MappingEulerian with such elements
+
+                                  // create collection objects from
+                                  // single quadratures,
+                                  // and finite elements. if we have
+                                  // an hp DoFHandler,
+                                  // dof_handler.get_fe() returns a
+                                  // collection of which we do a
+                                  // shallow copy instead
+  const hp::QCollection<dim>       q_collection (patch_points);
+  const hp::FECollection<dim>      fe_collection(this->dofs->get_fe());
+  
+  hp::FEValues<dim> x_fe_patch_values (fe_collection, q_collection,
+                                       update_values);
 
   const unsigned int n_q_points = patch_points.n_quadrature_points;
   
index 7b59a4f3c594a551882f400867eeb9cfe02c7266..31b2424b956c76ac53a4873e83d3e4dc473781f8 100644 (file)
@@ -24,7 +24,6 @@
 #include <fe/fe.h>
 #include <fe/fe_values.h>
 #include <fe/hp_fe_values.h>
-#include <fe/select_fe_values.h>
 #include <fe/mapping_q1.h>
 
 
@@ -35,9 +34,26 @@ void DataOutFaces<dim,DH>::build_some_patches (Data &data)
   QTrapez<1>        q_trapez;
   QIterated<dim-1>  patch_points (q_trapez, data.n_subdivisions);
   
-  typename SelectFEValues<DH<dim> >::FEFaceValues
-    x_fe_patch_values (this->dofs->get_fe(),
-                       patch_points, update_values);
+//TODO[?]: This is strange -- Data has a member 'mapping' that should
+//be used here, but it isn't. Rather, up until version 1.94, we were
+//actually initializing a local mapping object and used that... While
+//we use the mapping to transform the vertex coordinates, we do not
+//use the mapping to transform the shape functions (necessary for
+//example for Raviart-Thomas elements). This could lead to trouble
+//when someone tries to use MappingEulerian with such elements
+
+                                  // create collection objects from
+                                  // single quadratures,
+                                  // and finite elements. if we have
+                                  // an hp DoFHandler,
+                                  // dof_handler.get_fe() returns a
+                                  // collection of which we do a
+                                  // shallow copy instead
+  const hp::QCollection<dim-1>     q_collection (patch_points);
+  const hp::FECollection<dim>      fe_collection(this->dofs->get_fe());
+  
+  hp::FEFaceValues<dim> x_fe_patch_values (fe_collection, q_collection,
+                                           update_values);
 
   const unsigned int n_q_points = patch_points.n_quadrature_points;
   
index 5caa9afd7b17f1989939a762936edbda6715b0f6..679f237a451bc1411925c3f0c80155c573de50b2 100644 (file)
@@ -23,7 +23,6 @@
 #include <fe/fe.h>
 #include <fe/fe_values.h>
 #include <fe/hp_fe_values.h>
-#include <fe/select_fe_values.h>
 #include <fe/mapping_q1.h>
 #include <numerics/data_out_rotation.h>
 
@@ -38,14 +37,24 @@ void DataOutRotation<dim,DH>::build_some_patches (Data &data)
   QTrapez<1>     q_trapez;
   QIterated<dim> patch_points (q_trapez, data.n_subdivisions);
 
+                                  // create collection objects from
+                                  // single quadratures,
+                                  // and finite elements. if we have
+                                  // an hp DoFHandler,
+                                  // dof_handler.get_fe() returns a
+                                  // collection of which we do a
+                                  // shallow copy instead
+                                   //
                                   // since most output formats can't
                                   // handle cells that are not
                                   // transformed using a Q1 mapping,
                                   // we don't support anything else
                                   // as well
-  typename SelectFEValues<DH<dim> >::FEValues
-    x_fe_patch_values (StaticMappingQ1<dim>::mapping, this->dofs->get_fe(),
-                       patch_points, update_values);
+  const hp::QCollection<dim>       q_collection (patch_points);
+  const hp::FECollection<dim>      fe_collection(this->dofs->get_fe());
+  
+  hp::FEValues<dim> x_fe_patch_values (fe_collection, q_collection,
+                                       update_values);
 
   const unsigned int n_patches_per_circle = data.n_patches_per_circle;
 
index 41fb42845acf3cd8bcb5f09f5e79beecbe5bcfa1..4a800adc67ebaa578f96dd28d9d047283b09da18 100644 (file)
@@ -23,7 +23,6 @@
 #include <fe/fe.h>
 #include <fe/fe_values.h>
 #include <fe/hp_fe_values.h>
-#include <fe/select_fe_values.h>
 #include <fe/mapping_q1.h>
 
 #include <sstream>
@@ -257,9 +256,19 @@ void DataOutStack<dim,DH>::build_patches (const unsigned int n_subdivisions)
   QTrapez<1>     q_trapez;
   QIterated<dim> patch_points (q_trapez, n_subdivisions);
   
-  typename SelectFEValues<DH<dim> >::FEValues
-    x_fe_patch_values (dof_handler->get_fe(),
-                       patch_points, update_values);
+                                  // create collection objects from
+                                  // single quadratures,
+                                  // and finite elements. if we have
+                                  // an hp DoFHandler,
+                                  // dof_handler.get_fe() returns a
+                                  // collection of which we do a
+                                  // shallow copy instead
+  const hp::QCollection<dim>       q_collection (patch_points);
+  const hp::FECollection<dim>      fe_collection(dof_handler->get_fe());
+  
+  hp::FEValues<dim> x_fe_patch_values (fe_collection, q_collection,
+                                       update_values);
+
   const unsigned int n_q_points = patch_points.n_quadrature_points;
   std::vector<double>          patch_values (n_q_points);
   std::vector<Vector<double> > patch_values_system (n_q_points,

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.