From 956862a01d59267c0d5faeccd0dc9a64c2f146d0 Mon Sep 17 00:00:00 2001 From: Guido Kanschat Date: Fri, 12 Mar 2010 15:48:10 +0000 Subject: [PATCH] add MGVectorData git-svn-id: https://svn.dealii.org/trunk@20805 0785d39b-7218-0410-832d-ea1e28bc413d --- .../numerics/mesh_worker_vector_selector.h | 90 +++++++++++++++++++ .../mesh_worker_vector_selector.templates.h | 90 +++++++++++++++++++ .../mesh_worker_vector_selector.inst.in | 3 +- 3 files changed, 182 insertions(+), 1 deletion(-) diff --git a/deal.II/deal.II/include/numerics/mesh_worker_vector_selector.h b/deal.II/deal.II/include/numerics/mesh_worker_vector_selector.h index fa4599bf89..125568eff8 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker_vector_selector.h +++ b/deal.II/deal.II/include/numerics/mesh_worker_vector_selector.h @@ -16,6 +16,7 @@ #include #include #include +#include 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 > >& values, + std::vector > > >& gradients, + std::vector > > >& hessians, + const FEValuesBase& fe, + unsigned int level, + const std::vector& 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 MGVectorData : + public VectorDataBase + { + 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*>&); + + /** + * 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*, const std::string& name); + + virtual void mg_fill( + std::vector > >& values, + std::vector > > >& gradients, + std::vector > > >& hessians, + const FEValuesBase& fe, + unsigned int level, + const std::vector& index, + unsigned int component, + unsigned int n_comp, + unsigned int start, + unsigned int size) const; + private: + NamedData,MGVectorData > > data; + }; + + //----------------------------------------------------------------------// inline void diff --git a/deal.II/deal.II/include/numerics/mesh_worker_vector_selector.templates.h b/deal.II/deal.II/include/numerics/mesh_worker_vector_selector.templates.h index 6f1f855868..b54dbe4294 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker_vector_selector.templates.h +++ b/deal.II/deal.II/include/numerics/mesh_worker_vector_selector.templates.h @@ -49,6 +49,22 @@ namespace MeshWorker unsigned int, unsigned int) const {} + + + template + void + VectorDataBase::mg_fill( + std::vector > >&, + std::vector > > >&, + std::vector > > >&, + const FEValuesBase&, + unsigned int, + const std::vector&, + unsigned int, + unsigned int, + unsigned int, + unsigned int) const + {} //----------------------------------------------------------------------// template @@ -120,6 +136,80 @@ namespace MeshWorker fe.get_function_hessians(src, make_slice(index, start, size), dst, true); } } + +//----------------------------------------------------------------------// + + template + MGVectorData::MGVectorData() + {} + + + template + MGVectorData::MGVectorData(const VectorSelector& s) + : + VectorDataBase(s) + {} + + + template + void + MGVectorData::initialize(const NamedData*>& d) + { + data = d; + VectorSelector::initialize(d); + } + + + template + void + MGVectorData::initialize(const MGLevelObject* v, const std::string& name) + { + SmartPointer, MGVectorData > + p = v; + data.add(p, name); + VectorSelector::initialize(data); + } + + + template + void + MGVectorData::mg_fill( + std::vector > >& values, + std::vector > > >& gradients, + std::vector > > >& hessians, + const FEValuesBase& fe, + unsigned int level, + const std::vector& 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;in_values();++i) + { + const VECTOR& src = (*data(this->value_index(i)))[level]; + VectorSlice > > dst(values[i], component, n_comp); + fe.get_function_values(src, make_slice(index, start, size), dst, true); + } + + for (unsigned int i=0;in_gradients();++i) + { + const VECTOR& src = (*data(this->value_index(i)))[level]; + VectorSlice > > > dst(gradients[i], component, n_comp); + fe.get_function_gradients(src, make_slice(index, start, size), dst, true); + } + + for (unsigned int i=0;in_hessians();++i) + { + const VECTOR& src = (*data(this->value_index(i)))[level]; + VectorSlice > > > dst(hessians[i], component, n_comp); + fe.get_function_hessians(src, make_slice(index, start, size), dst, true); + } + } } DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/deal.II/source/numerics/mesh_worker_vector_selector.inst.in b/deal.II/deal.II/source/numerics/mesh_worker_vector_selector.inst.in index 17ba91a9e5..743ef908ce 100644 --- a/deal.II/deal.II/source/numerics/mesh_worker_vector_selector.inst.in +++ b/deal.II/deal.II/source/numerics/mesh_worker_vector_selector.inst.in @@ -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; + template class MGVectorData; } -- 2.39.5