#include <base/named_data.h>
#include <base/tensor.h>
#include <base/smartpointer.h>
+#include <multigrid/mg_level_object.h>
DEAL_II_NAMESPACE_OPEN
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;
};
};
+/**
+ * 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
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>
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