]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
making MeshWorker use data vectors
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 11 Jan 2010 08:51:37 +0000 (08:51 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 11 Jan 2010 08:51:37 +0000 (08:51 +0000)
git-svn-id: https://svn.dealii.org/trunk@20348 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/named_data.h
deal.II/deal.II/include/numerics/mesh_worker.h
deal.II/deal.II/include/numerics/mesh_worker_assembler.h
deal.II/deal.II/include/numerics/mesh_worker_info.h
deal.II/deal.II/include/numerics/mesh_worker_info.templates.h
deal.II/deal.II/include/numerics/mesh_worker_vector_selector.h
deal.II/deal.II/include/numerics/mesh_worker_vector_selector.templates.h
deal.II/deal.II/source/numerics/mesh_worker.cc

index f8c747dc468eb8abefc4324fdbfd9259b0361e59..b398600dee6fa7bd20145cd5d380d97abb4d2dd7 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2009 by the deal.II authors
+//    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -106,12 +106,26 @@ class NamedData : public Subscriptor
 /// Number of stored data objects.
     unsigned int size() const;
 
-/// Access to stored data object by index.
+                                    /**
+                                     * @brief Access to stored data
+                                     * object by index.
+                                     *
+                                     * @note This function throws an
+                                     * exception, if is_const()
+                                     * returns <tt>true</tt>. In such
+                                     * a case, either cast the
+                                     * NamedData object to a const
+                                     * reference, or use the function
+                                     * read() instead of this operator.
+                                     */
     DATA& operator() (unsigned int i);
       
 /// Read-only access to stored data object by index.
     const DATA& operator() (unsigned int i) const;
-      
+
+/// Read only access for a non-const object.
+    const DATA& read (unsigned int i) const;
+    
 /// Name of object at index.
     const std::string& name(unsigned int i) const;
       
@@ -260,10 +274,11 @@ void
 NamedData<DATA>::merge(NamedData<DATA2>& other)
 {
   Assert(!is_constant, ExcConstantObject());
+  
   for (unsigned int i=0;i<other.size();++i)
     {
       names.push_back(other.name(i));
-      data.push_back(other(i));
+      data.push_back(other.read(i));
     }
   is_constant = other.is_const();
 }
@@ -339,6 +354,16 @@ NamedData<DATA>::operator() (unsigned int i) const
 }
   
   
+template<typename DATA>
+inline
+const DATA&
+NamedData<DATA>::read (unsigned int i) const
+{
+  AssertIndexRange(i, size());
+  return data[i];
+}
+  
+  
 template<typename DATA>
 inline
 const std::string&
index e950c696372a415b2087d7274ab31e226a0b77da..5f956c2f3fe70085412eb1e1d04f3dd9a214a29e 100644 (file)
@@ -170,6 +170,32 @@ namespace MeshWorker
 /**
  * Worker object for integration of functionals, residuals or matrices.
  *
+ * In order to allow for sufficient generality, a few steps have to be
+ * undertaken to use this class. The constructor will only fill some
+ * default values.
+ *
+ * First, you should consider if you need values from any vectors in a
+ * NamedData object. If so, fill the VectorSelector objects
+ * #cell_selector, #bdry_selector and #face_selector with their names
+ * and the data type (value, gradient, Hessian) to be extracted.
+ *
+ * Afterwards, you will need to consider UpdateFlags for FEValues
+ * objects. A good start is initialize_update_flags(), which looks at
+ * the selectors filled before and adds all the flags needed to get
+ * the selection. Additional flags can be set with add_update_flags().
+ *
+ * Finally, we need to choose quadrature formulas. If you choose to
+ * use Gauss formulas only, use initialize_gauss_quadrature() with
+ * appropriate values. Otherwise, you can fill the variables
+ * #cell_quadrature, #bdry_quadrature and #face_quadrature directly.
+ *
+ * In order to save time, you can set the variables #boundary_fluxes
+ * and #interior_fluxes of the base class to false, thus telling the
+ * Meshworker::loop() not to loop over those faces.
+ *
+ * All the information in here is used to set up IntegrationInfo
+ * objects correctly, typically in an IntegrationInfoBox.
+ *
  * @ingroup MeshWorker
  * @author Guido Kanschat, 2009
  */
@@ -211,17 +237,8 @@ namespace MeshWorker
                                        * add UpdateFlags to the
                                        * flags stored in this class.
                                        */
-      void initialize_selectors(const VectorSelector& cell_selector,
-                               const VectorSelector& bdry_selector,
-                               const VectorSelector& face_selector);
-
-                                      /**
-                                       * Add a vector to some or
-                                       * all selectors.
-                                       */
-      void add_selector(const std::string& name, bool values, bool gradients, bool hessians,
-                       bool cell, bool bdry, bool face);
-
+      void initialize_update_flags();
+      
                                       /**
                                        * Add additional values for update.
                                        */
index c1cc2ce38883871d551a8bc1ec6d8dec11ff0856..c822c195bbef535315eccdd56a371d8319bb46a7 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
 //
-//    Copyright (C) 2006, 2007, 2008, 2009 by Guido Kanschat
+//    Copyright (C) 2006, 2007, 2008, 2009, 2010 by Guido Kanschat
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -338,7 +338,8 @@ namespace MeshWorker
                                          * matrix pointers into local
                                          * variables.
                                          */
-       void initialize(NamedData<VECTOR*>& residuals);
+       void initialize(const BlockInfo* block_info,
+                       NamedData<VECTOR*>& residuals);
                                         /**
                                          * Initialize the local data
                                          * in the
@@ -380,9 +381,8 @@ namespace MeshWorker
                                          * Assemble a single local
                                          * residual into the global.
                                          */
-       void assemble(VECTOR global,
+       void assemble(VECTOR& global,
                      const BlockVector<double>& local,
-                     const BlockIndices& bi,
                      const std::vector<unsigned int>& dof);
        
                                         /**
@@ -391,7 +391,12 @@ namespace MeshWorker
                                          * pointers.
                                          */
        NamedData<SmartPointer<VECTOR,ResidualLocalBlocksToGlobalBlocks<VECTOR> > > residuals;
-    };
+       
+      /**
+       * A pointer to the object containing the block structure.
+       */
+      SmartPointer<const BlockInfo> block_info;
+   };
 
 
 /**
@@ -1134,8 +1139,10 @@ namespace MeshWorker
 
     template <class VECTOR>
     inline void
-    ResidualLocalBlocksToGlobalBlocks<VECTOR>::initialize(NamedData<VECTOR*>& m)
+    ResidualLocalBlocksToGlobalBlocks<VECTOR>::initialize(const BlockInfo* b,
+                                                         NamedData<VECTOR*>& m)
     {
+      block_info = b;
       residuals = m;      
     }
 
@@ -1152,9 +1159,8 @@ namespace MeshWorker
     template <class VECTOR>
     inline void
     ResidualLocalBlocksToGlobalBlocks<VECTOR>::assemble(
-      VECTOR global,
+      VECTOR& global,
       const BlockVector<double>& local,
-      const BlockIndices& bi,
       const std::vector<unsigned int>& dof)
     {
       for (unsigned int b=0;b<local.n_blocks();++b)
@@ -1168,8 +1174,8 @@ namespace MeshWorker
                                             // block-wise local
                                             // numbering we use in
                                             // our local vectors
-           const unsigned int jcell = this->bi.local_to_global(b, j);
-           (*global)(dof[jcell]) += local.block(b)(j);
+           const unsigned int jcell = this->block_info->local().local_to_global(b, j);
+           global(dof[jcell]) += local.block(b)(j);
          }
     }
 
@@ -1180,8 +1186,8 @@ namespace MeshWorker
     ResidualLocalBlocksToGlobalBlocks<VECTOR>::assemble(
       const DoFInfo<dim>& info)
     {
-      for (unsigned int i=0;i<residuals.n_vectors();++i)
-       assemble(residuals.vector(i), info.R[i], info.block_info.local(), info.indices);
+      for (unsigned int i=0;i<residuals.size();++i)
+       assemble(*residuals(i), info.R[i], info.indices);
     }
 
     
@@ -1192,10 +1198,10 @@ namespace MeshWorker
       const DoFInfo<dim>& info1,
       const DoFInfo<dim>& info2)
     {
-      for (unsigned int i=0;i<residuals.n_vectors();++i)
+      for (unsigned int i=0;i<residuals.size();++i)
        {
-         assemble(residuals.vector(i), info1.R[i], info1.block_info.local(), info1.indices);
-         assemble(residuals.vector(i), info2.R[i], info2.block_info.local(), info2.indices);
+         assemble(*residuals(i), info1.R[i], info1.indices);
+         assemble(*residuals(i), info2.R[i], info2.indices);
        }
     }
 
@@ -1402,6 +1408,18 @@ namespace MeshWorker
     
 
     
+    template <class MATRIX ,typename number>
+    template <int dim>
+    inline void
+    MatrixLocalBlocksToGlobalBlocks<MATRIX, number>::initialize_info(
+      DoFInfo<dim>& info,
+      bool interior_face) const
+    {
+      info.initialize_matrices(matrices, interior_face);
+    }
+
+
+
     template <class MATRIX, typename number>
     inline void
     MatrixLocalBlocksToGlobalBlocks<MATRIX, number>::assemble(
@@ -1504,6 +1522,18 @@ namespace MeshWorker
     }
     
 
+    template <class MATRIX ,typename number>
+    template <int dim>
+    inline void
+    MGMatrixLocalBlocksToGlobalBlocks<MATRIX, number>::initialize_info(
+      DoFInfo<dim>& info,
+      bool interior_face) const
+    {
+      info.initialize_matrices(matrices, interior_face);
+    }
+
+
+
     template <class MATRIX, typename number>
     inline void
     MGMatrixLocalBlocksToGlobalBlocks<MATRIX, number>::initialize_edge_flux(
index a2ceadb6351cef7d1d2591d0ec2fb6ec27558bf0..795410beceb21890ff6ad2814cfc1aaf750431fc 100644 (file)
@@ -85,7 +85,7 @@ namespace MeshWorker
                                        * provide row and column info.
                                        */
       template <class MatrixPtr>
-      void initialize_matrices(std::vector<MatrixPtr>& matrices,
+      void initialize_matrices(const std::vector<MatrixPtr>& matrices,
                               bool both);
        
                                       /**
@@ -526,7 +526,11 @@ namespace MeshWorker
        std::vector<std::vector<std::vector<TYPE> > >& data,
        VectorSelector& selector,
        bool split_fevalues) const;
-       
+                                      /**
+                                       * Cache the number of
+                                       * components of the system element.
+                                       */
+      unsigned int n_components;
   };
 
 /**
@@ -549,19 +553,23 @@ namespace MeshWorker
                                        */
       template <typename T>
       IntegrationInfoBox(const T&);
-
-      template <class VECTOR>
-      void initialize_data(const VECTOR*, const std::string& name,
-                          bool values, bool gradients, bool hessians);
       
       template <class WORKER>
       void initialize(const WORKER&,
                      const FiniteElement<dim, spacedim>& el,
                      const Mapping<dim, spacedim>& mapping);
+
+      template <class WORKER, typename VECTOR>
+      void initialize(const WORKER&,
+                     const FiniteElement<dim, spacedim>& el,
+                     const Mapping<dim, spacedim>& mapping,
+                     const NamedData<VECTOR*>& data);
        
 //      private:
 
-      boost::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > data;
+      boost::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > cell_data;
+      boost::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > bdry_data;
+      boost::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > face_data;
       
       CellInfo cell_info;
       FaceInfo bdry_info;
@@ -593,7 +601,7 @@ namespace MeshWorker
   template <class MatrixPtr>
   inline void
   LocalResults<number>::initialize_matrices(
-    std::vector<MatrixPtr>& matrices,
+    const std::vector<MatrixPtr>& matrices,
     bool both)
   {
     M1.resize(matrices.size());
@@ -839,27 +847,6 @@ namespace MeshWorker
   {}
   
   
-  template <int dim, int sdim>
-  template <typename VECTOR>
-  void
-  IntegrationInfoBox<dim,sdim>::initialize_data(
-    const VECTOR* v, const std::string& name,
-    bool values, bool gradients, bool hessians)
-  {
-    boost::shared_ptr<VectorData<VECTOR, dim, sdim> >
-      p = boost::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> ());
-    p->add(name, values, gradients, hessians);
-    p->initialize(v, name);
-    data = p;
-    
-    cell_info.initialize_data(data);
-    bdry_info.initialize_data(data);
-    face_info.initialize_data(data);
-    subface_info.initialize_data(data);
-    neighbor_info.initialize_data(data);
-  }
-
-  
   template <int dim, int sdim>
   template <class WORKER>
   void
@@ -883,11 +870,41 @@ namespace MeshWorker
     subface_info.initialize<FESubfaceValues<dim,sdim> >(el, mapping, integrator.face_quadrature,
                            integrator.face_flags);
     neighbor_info.initialize<FEFaceValues<dim,sdim> >(el, mapping, integrator.face_quadrature,
-                            integrator.ngbr_flags);
-  }    
-}
+                                                     integrator.ngbr_flags);
+  }
 
 
+  template <int dim, int sdim>
+  template <class WORKER, typename VECTOR>
+  void
+  IntegrationInfoBox<dim,sdim>::initialize(
+    const WORKER& integrator,
+    const FiniteElement<dim,sdim>& el,
+    const Mapping<dim,sdim>& mapping,
+    const NamedData<VECTOR*>& data)
+  {
+    initialize(integrator, el, mapping);
+    boost::shared_ptr<VectorData<VECTOR, dim, sdim> > p;
+    
+    p = boost::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (integrator.cell_selector));
+    p->initialize(data);
+    cell_data = p;
+    cell_info.initialize_data(p);
+    
+    p = boost::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (integrator.bdry_selector));
+    p->initialize(data);
+    bdry_data = p;
+    bdry_info.initialize_data(p);
+        
+    p = boost::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (integrator.face_selector));
+    p->initialize(data);
+    face_data = p;
+    face_info.initialize_data(p);
+    subface_info.initialize_data(p);
+    neighbor_info.initialize_data(p);
+  }
+}
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif
index 95c5274699c24cc090934277026a460c42d7122f..3e2da89c2159b2d68fbd41743cd2412a6a3b7978 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2009 by the deal.II authors
+//    Copyright (C) 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -96,30 +96,46 @@ namespace MeshWorker
              new FEVALUES (mapping, el.base_element(i), quadrature, flags));
          }
       }
+    n_components = el.n_components();
+  }
+  
 
+  template<int dim, class FVB, int sdim>
+  void
+  IntegrationInfo<dim,FVB,sdim>::initialize_data(
+    const boost::shared_ptr<VectorDataBase<dim,sdim> > data)
+  {
+    global_data = data;
+    const unsigned int nqp = fevalv[0]->n_quadrature_points;
+    
     values.resize(global_data->n_values());
+//    deallog << "values: " << values.size() << " [";
                                     // For all selected finite
                                     // element functions
     for (unsigned int i=0;i<values.size();++i)
       {
-       values[i].resize(this->local_indices().size());
+       values[i].resize(n_components);
+//     deallog << ' ' << values[i].size() << " {";
                                         // For all components
        for (unsigned int j=0;j<values[i].size();++j)
          {
-           values[i][j].resize(quadrature.size());
+           values[i][j].resize(nqp);
+//         deallog << ' ' << values[i][j].size();
          }
+//     deallog << " }";
       }
+//    deallog << " ]" << std::endl;
     
     gradients.resize(global_data->n_gradients());
                                     // For all selected finite
                                     // element functions
     for (unsigned int i=0;i<gradients.size();++i)
       {
-       gradients[i].resize(this->local_indices().size());
+       gradients[i].resize(n_components);
                                         // For all components
        for (unsigned int j=0;j<gradients[i].size();++j)
          {
-           gradients[i][j].resize(quadrature.size());
+           gradients[i][j].resize(nqp);
          }
       }
     
@@ -128,23 +144,14 @@ namespace MeshWorker
                                     // element functions
     for (unsigned int i=0;i<hessians.size();++i)
       {
-       hessians[i].resize(this->local_indices().size());
+       hessians[i].resize(n_components);
                                         // For all components
        for (unsigned int j=0;j<hessians[i].size();++j)
          {
-           hessians[i][j].resize(quadrature.size());
+           hessians[i][j].resize(nqp);
          }
       }
   }
-  
-
-  template<int dim, class FVB, int sdim>
-  void
-  IntegrationInfo<dim,FVB,sdim>::initialize_data(
-    const boost::shared_ptr<VectorDataBase<dim,sdim> > data)
-  {
-    global_data = data;
-  }
 
 
   template<int dim, class FVB, int sdim>
index cab71251f3d83558a077e55d58c060dd50cc6ebe..fa4599bf89d1deb203c77bde54c64768403cfd2a 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2009 by the deal.II authors
+//    Copyright (C) 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -27,7 +27,7 @@ namespace MeshWorker
 /**
  * A class that selects vectors from a list of named vectors.
  *
- * Since the number of vectors in FEVectors may grow with every
+ * Since the number of vectors in a NamedData object may grow with every
  * nesting of applications or loops, it is important to be able to
  * select those, which are actually used in computing residuals etc.
  * This class organizes the selection.
@@ -44,11 +44,22 @@ namespace MeshWorker
     public:
                                       /**
                                        * Add a vector to the
-                                       * selection. The arguments are
+                                       * selection of finite element
+                                       * functions. The arguments are
                                        * the name of the vector and
                                        * indicators, which
                                        * information is to be
-                                       * extracted from the vector.
+                                       * extracted from the
+                                       * vector. The name refers to
+                                       * an entry in a NamedData
+                                       * object, which will be
+                                       * identified by initialize().
+                                       * The three bool parameters
+                                       * indicate, whether values,
+                                       * greadients and Hessians of
+                                       * the finite element function
+                                       * are to be computed on each
+                                       * cell or face.
                                        */
       void add(const std::string& name,
               bool values = true,
@@ -57,7 +68,20 @@ namespace MeshWorker
 
                                       /**
                                        * Initialize the selection
-                                       * field with a data vector.
+                                       * field with a data
+                                       * vector. While add() only
+                                       * enters the names of vectors,
+                                       * which will be used in the integration
+                                       * loop over cells and faces,
+                                       * this function links the
+                                       * names to actual vectos in a
+                                       * NamedData object.
+                                       *
+                                       * @note This function caches
+                                       * the index associated with a
+                                       * name. Therefore, it must be
+                                       * called everytime after the
+                                       * NamedData object has changed.
                                        */
       template <class DATA>
       void initialize(const NamedData<DATA>&);
@@ -159,12 +183,22 @@ namespace MeshWorker
       public VectorSelector
   {
     public:
+                                      /**
+                                       * Constructor
+                                       */
+      VectorDataBase();
+      
+                                      /**
+                                       * Constructor from a base
+                                       * class object
+                                       */
+      VectorDataBase(const VectorSelector&);
+      
                                       /**
                                        * Virtual, but empty
                                        * destructor.
                                        */
       virtual ~VectorDataBase();
-      
                                       /**
                                        * The only function added to
                                        * VectorSelector is an
@@ -257,7 +291,42 @@ namespace MeshWorker
       public VectorDataBase<dim, spacedim>
   {
     public:
+                                      /**
+                                       * Constructor.
+                                       */
+      VectorData();
+                                      /**
+                                       * Constructor using a
+                                       * prefilled VectorSelector
+                                       */
+      VectorData(const VectorSelector&);
+
+                                      /**
+                                       * Initialize with a NamedData
+                                       * object and cache the indices
+                                       * in the VectorSelector base
+                                       * class.
+                                       *
+                                       * @note Make sure the
+                                       * VectorSelector base class
+                                       * was filled with reasonable
+                                       * data before calling this
+                                       * function.
+                                       */
       void initialize(const NamedData<VECTOR*>&);
+
+                                      /**
+                                       * Initialize with a single
+                                       * vector and cache the indices
+                                       * in the VectorSelector base
+                                       * class.
+                                       *
+                                       * @note Make sure the
+                                       * VectorSelector base class
+                                       * was filled with reasonable
+                                       * data before calling this
+                                       * function.
+                                       */
       void initialize(const VECTOR*, const std::string& name);
       
       virtual void fill(
index e691467949996e712d72000756944903fe45f4a6..6f1f855868539f9e9beb4bd88a79c86413d5ff50 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2009 by the deal.II authors
+//    Copyright (C) 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -23,6 +23,19 @@ namespace MeshWorker
   VectorDataBase<dim, spacedim>::~VectorDataBase()
   {}
   
+
+  template <int dim, int spacedim>
+  VectorDataBase<dim, spacedim>::VectorDataBase(const VectorSelector& v)
+                 :
+                 VectorSelector(v)
+  {}  
+
+
+  template <int dim, int spacedim>
+  VectorDataBase<dim, spacedim>::VectorDataBase()
+  {}  
+
+
   template <int dim, int spacedim>
   void
   VectorDataBase<dim, spacedim>::fill(
@@ -38,6 +51,18 @@ namespace MeshWorker
   {}
 //----------------------------------------------------------------------//
 
+  template <class VECTOR, int dim, int spacedim>
+  VectorData<VECTOR, dim, spacedim>::VectorData()
+  {}
+  
+  
+  template <class VECTOR, int dim, int spacedim>
+  VectorData<VECTOR, dim, spacedim>::VectorData(const VectorSelector& s)
+                 :
+                 VectorDataBase<dim, spacedim>(s)
+  {}
+  
+  
   template <class VECTOR, int dim, int spacedim>
   void
   VectorData<VECTOR, dim, spacedim>::initialize(const NamedData<VECTOR*>& d)
@@ -70,6 +95,10 @@ namespace MeshWorker
        unsigned int start,
        unsigned int size) const
   {
+    AssertDimension(values.size(), this->n_values());
+    AssertDimension(gradients.size(), this->n_gradients());
+    AssertDimension(hessians.size(), this->n_hessians());
+    
     for (unsigned int i=0;i<this->n_values();++i)
       {
        const VECTOR& src = *data(this->value_index(i));
index bea95926157ca1aff40eb88216baf0ee88b654c6..a89e80acc964cc3952b86c44117ee90f3f424218 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
 //
-//    Copyright (C) 2006, 2007, 2008, 2009 by Guido Kanschat
+//    Copyright (C) 2006, 2007, 2008, 2009, 2010 by Guido Kanschat
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -31,15 +31,8 @@ IntegrationWorker<dim>::IntegrationWorker ()
 
 template<int dim>
 void
-IntegrationWorker<dim>::initialize_selectors(
-  const VectorSelector& cs,
-  const VectorSelector& bs,
-  const VectorSelector& fs)
+IntegrationWorker<dim>::initialize_update_flags ()
 {
-  cell_selector = cs;
-  bdry_selector = bs;
-  face_selector = fs;
-  
   if (cell_selector.has_values() != 0) cell_flags |= update_values;
   if (cell_selector.has_gradients() != 0) cell_flags |= update_gradients;
   if (cell_selector.has_hessians() != 0) cell_flags |= update_hessians;
@@ -58,34 +51,6 @@ IntegrationWorker<dim>::initialize_selectors(
 }
 
 
-template<int dim>
-void
-IntegrationWorker<dim>::add_selector(
-  const std::string& name, bool values, bool gradients, bool hessians,
-  bool cell, bool bdry, bool face)
-{
-  if (cell) cell_selector.add(name, values, gradients, hessians);
-  if (bdry) bdry_selector.add(name, values, gradients, hessians);
-  if (face) face_selector.add(name, values, gradients, hessians);  
-  
-  if (cell_selector.has_values() != 0) cell_flags |= update_values;
-  if (cell_selector.has_gradients() != 0) cell_flags |= update_gradients;
-  if (cell_selector.has_hessians() != 0) cell_flags |= update_hessians;
-
-  if (bdry_selector.has_values() != 0) bdry_flags |= update_values;
-  if (bdry_selector.has_gradients() != 0) bdry_flags |= update_gradients;
-  if (bdry_selector.has_hessians() != 0) bdry_flags |= update_hessians;
-  
-  if (face_selector.has_values() != 0) face_flags |= update_values;
-  if (face_selector.has_gradients() != 0) face_flags |= update_gradients;
-  if (face_selector.has_hessians() != 0) face_flags |= update_hessians;
-  
-  if (face_selector.has_values() != 0) ngbr_flags |= update_values;
-  if (face_selector.has_gradients() != 0) ngbr_flags |= update_gradients;
-  if (face_selector.has_hessians() != 0) ngbr_flags |= update_hessians;  
-}
-
-
 template<int dim>
 void
 IntegrationWorker<dim>::add_update_flags(

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.