]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add MGVectorData
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Fri, 12 Mar 2010 15:48:10 +0000 (15:48 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Fri, 12 Mar 2010 15:48:10 +0000 (15:48 +0000)
git-svn-id: https://svn.dealii.org/trunk@20805 0785d39b-7218-0410-832d-ea1e28bc413d

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_vector_selector.inst.in

index fa4599bf89d1deb203c77bde54c64768403cfd2a..125568eff8f109ff336280fb90d85f6e673356e6 100644 (file)
@@ -16,6 +16,7 @@
 #include <base/named_data.h>
 #include <base/tensor.h>
 #include <base/smartpointer.h>
+#include <multigrid/mg_level_object.h>
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -276,6 +277,29 @@ namespace MeshWorker
        unsigned int n_comp,
        unsigned int start,
        unsigned int size) const;
+
+                                      /**
+                                       * Fill the local data vector
+                                       * from level vectors. Performs
+                                       * exactly what the other
+                                       * fill() does, but uses the
+                                       * cell level to access a
+                                       * single level out of a
+                                       * hierarchy of level vectors,
+                                       * instead of a global data
+                                       * vector on the active cells.
+                                       */
+       virtual void mg_fill(
+       std::vector<std::vector<std::vector<double> > >& values,
+       std::vector<std::vector<std::vector<Tensor<1,dim> > > >& gradients,
+       std::vector<std::vector<std::vector<Tensor<2,dim> > > >& hessians,
+       const FEValuesBase<dim,spacedim>& fe,
+       unsigned int level,
+       const std::vector<unsigned int>& index,
+       unsigned int component,
+       unsigned int n_comp,
+       unsigned int start,
+       unsigned int size) const;     
   };
 
   
@@ -344,6 +368,72 @@ namespace MeshWorker
   };
 
   
+/**
+ * Based on VectorSelector, this is the class that implements the
+ * function VectorDataBase::fill() for a certain type of multilevel vectors, using
+ * NamedData to identify vectors by name.
+ *
+ * @author Guido Kanschat, 2009
+ */
+  template <class VECTOR, int dim, int spacedim = dim>
+  class MGVectorData :
+      public VectorDataBase<dim, spacedim>
+  {
+    public:
+                                      /**
+                                       * Constructor.
+                                       */
+      MGVectorData();
+                                      /**
+                                       * Constructor using a
+                                       * prefilled VectorSelector
+                                       */
+      MGVectorData(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<MGLevelObject<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 MGLevelObject<VECTOR>*, const std::string& name);
+      
+      virtual void mg_fill(
+       std::vector<std::vector<std::vector<double> > >& values,
+       std::vector<std::vector<std::vector<Tensor<1,dim> > > >& gradients,
+       std::vector<std::vector<std::vector<Tensor<2,dim> > > >& hessians,
+       const FEValuesBase<dim,spacedim>& fe,
+       unsigned int level,
+       const std::vector<unsigned int>& index,
+       unsigned int component,
+       unsigned int n_comp,
+       unsigned int start,
+       unsigned int size) const;
+    private:
+      NamedData<SmartPointer<const MGLevelObject<VECTOR>,MGVectorData<VECTOR,dim,spacedim> > > data;
+  };
+
+  
 //----------------------------------------------------------------------//
 
   inline void
index 6f1f855868539f9e9beb4bd88a79c86413d5ff50..b54dbe429494bfb58aa3b15e9dff574e8a4a43c8 100644 (file)
@@ -49,6 +49,22 @@ namespace MeshWorker
        unsigned int,
        unsigned int) const
   {}
+
+
+  template <int dim, int spacedim>
+  void
+  VectorDataBase<dim, spacedim>::mg_fill(
+       std::vector<std::vector<std::vector<double> > >&,
+       std::vector<std::vector<std::vector<Tensor<1,dim> > > >&,
+       std::vector<std::vector<std::vector<Tensor<2,dim> > > >&,
+       const FEValuesBase<dim,spacedim>&,
+       unsigned int,
+       const std::vector<unsigned int>&,
+       unsigned int,
+       unsigned int,
+       unsigned int,
+       unsigned int) const
+  {}
 //----------------------------------------------------------------------//
 
   template <class VECTOR, int dim, int spacedim>
@@ -120,6 +136,80 @@ namespace MeshWorker
        fe.get_function_hessians(src, make_slice(index, start, size), dst, true);
       }
   }
+
+//----------------------------------------------------------------------//
+
+  template <class VECTOR, int dim, int spacedim>
+  MGVectorData<VECTOR, dim, spacedim>::MGVectorData()
+  {}
+  
+  
+  template <class VECTOR, int dim, int spacedim>
+  MGVectorData<VECTOR, dim, spacedim>::MGVectorData(const VectorSelector& s)
+                 :
+                 VectorDataBase<dim, spacedim>(s)
+  {}
+  
+  
+  template <class VECTOR, int dim, int spacedim>
+  void
+  MGVectorData<VECTOR, dim, spacedim>::initialize(const NamedData<MGLevelObject<VECTOR>*>& d)
+  {
+    data = d;
+    VectorSelector::initialize(d);
+  }
+  
+  
+  template <class VECTOR, int dim, int spacedim>
+  void
+  MGVectorData<VECTOR, dim, spacedim>::initialize(const MGLevelObject<VECTOR>* v, const std::string& name)
+  {
+    SmartPointer<const MGLevelObject<VECTOR>, MGVectorData<VECTOR, dim, spacedim> >
+      p = v;
+    data.add(p, name);
+    VectorSelector::initialize(data);
+  }
+  
+  
+  template <class VECTOR, int dim, int spacedim>
+  void
+  MGVectorData<VECTOR, dim, spacedim>::mg_fill(
+       std::vector<std::vector<std::vector<double> > >& values,
+       std::vector<std::vector<std::vector<Tensor<1,dim> > > >& gradients,
+       std::vector<std::vector<std::vector<Tensor<2,dim> > > >& hessians,
+       const FEValuesBase<dim,spacedim>& fe,
+       unsigned int level,
+       const std::vector<unsigned int>& index,
+       unsigned int component,
+       unsigned int n_comp,
+       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)))[level];
+       VectorSlice<std::vector<std::vector<double> > > dst(values[i], component, n_comp);
+       fe.get_function_values(src, make_slice(index, start, size), dst, true);
+      }
+    
+    for (unsigned int i=0;i<this->n_gradients();++i)
+      {
+       const VECTOR& src = (*data(this->value_index(i)))[level];
+       VectorSlice<std::vector<std::vector<Tensor<1,dim> > > > dst(gradients[i], component, n_comp);
+       fe.get_function_gradients(src, make_slice(index, start, size), dst, true);
+      }
+    
+    for (unsigned int i=0;i<this->n_hessians();++i)
+      {
+       const VECTOR& src = (*data(this->value_index(i)))[level];
+       VectorSlice<std::vector<std::vector<Tensor<2,dim> > > > dst(hessians[i], component, n_comp);
+       fe.get_function_hessians(src, make_slice(index, start, size), dst, true);
+      }
+  }
 }
 
 DEAL_II_NAMESPACE_CLOSE
index 17ba91a9e52a93ee019532c5efdc6c148d42e971..743ef908ceb5b4d157f6a4d19b7e0f6f056dc195 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id: fe_field_function.inst.in 19046 2009-07-08 19:30:23Z bangerth $
 //
-//    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
@@ -14,4 +14,5 @@
 for (VECTOR : SERIAL_VECTORS)
 {
   template class VectorData<VECTOR,deal_II_dimension>;
+  template class MGVectorData<VECTOR,deal_II_dimension>;
 }

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.