From e48be39e8960a04fcd6aeb5a722d06d2afa860c2 Mon Sep 17 00:00:00 2001 From: Guido Kanschat Date: Mon, 11 Jan 2010 08:51:37 +0000 Subject: [PATCH] making MeshWorker use data vectors git-svn-id: https://svn.dealii.org/trunk@20348 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/base/include/base/named_data.h | 33 +++++++- .../deal.II/include/numerics/mesh_worker.h | 39 ++++++--- .../include/numerics/mesh_worker_assembler.h | 60 ++++++++++---- .../include/numerics/mesh_worker_info.h | 81 +++++++++++-------- .../numerics/mesh_worker_info.templates.h | 39 +++++---- .../numerics/mesh_worker_vector_selector.h | 81 +++++++++++++++++-- .../mesh_worker_vector_selector.templates.h | 31 ++++++- .../deal.II/source/numerics/mesh_worker.cc | 39 +-------- 8 files changed, 281 insertions(+), 122 deletions(-) diff --git a/deal.II/base/include/base/named_data.h b/deal.II/base/include/base/named_data.h index f8c747dc46..b398600dee 100644 --- a/deal.II/base/include/base/named_data.h +++ b/deal.II/base/include/base/named_data.h @@ -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 true. 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::merge(NamedData& other) { Assert(!is_constant, ExcConstantObject()); + for (unsigned int i=0;i::operator() (unsigned int i) const } +template +inline +const DATA& +NamedData::read (unsigned int i) const +{ + AssertIndexRange(i, size()); + return data[i]; +} + + template inline const std::string& diff --git a/deal.II/deal.II/include/numerics/mesh_worker.h b/deal.II/deal.II/include/numerics/mesh_worker.h index e950c69637..5f956c2f3f 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker.h +++ b/deal.II/deal.II/include/numerics/mesh_worker.h @@ -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. */ diff --git a/deal.II/deal.II/include/numerics/mesh_worker_assembler.h b/deal.II/deal.II/include/numerics/mesh_worker_assembler.h index c1cc2ce388..c822c195bb 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker_assembler.h +++ b/deal.II/deal.II/include/numerics/mesh_worker_assembler.h @@ -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& residuals); + void initialize(const BlockInfo* block_info, + NamedData& 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& local, - const BlockIndices& bi, const std::vector& dof); /** @@ -391,7 +391,12 @@ namespace MeshWorker * pointers. */ NamedData > > residuals; - }; + + /** + * A pointer to the object containing the block structure. + */ + SmartPointer block_info; + }; /** @@ -1134,8 +1139,10 @@ namespace MeshWorker template inline void - ResidualLocalBlocksToGlobalBlocks::initialize(NamedData& m) + ResidualLocalBlocksToGlobalBlocks::initialize(const BlockInfo* b, + NamedData& m) { + block_info = b; residuals = m; } @@ -1152,9 +1159,8 @@ namespace MeshWorker template inline void ResidualLocalBlocksToGlobalBlocks::assemble( - VECTOR global, + VECTOR& global, const BlockVector& local, - const BlockIndices& bi, const std::vector& dof) { for (unsigned int b=0;bbi.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::assemble( const DoFInfo& info) { - for (unsigned int i=0;i& info1, const DoFInfo& info2) { - for (unsigned int i=0;i + template + inline void + MatrixLocalBlocksToGlobalBlocks::initialize_info( + DoFInfo& info, + bool interior_face) const + { + info.initialize_matrices(matrices, interior_face); + } + + + template inline void MatrixLocalBlocksToGlobalBlocks::assemble( @@ -1504,6 +1522,18 @@ namespace MeshWorker } + template + template + inline void + MGMatrixLocalBlocksToGlobalBlocks::initialize_info( + DoFInfo& info, + bool interior_face) const + { + info.initialize_matrices(matrices, interior_face); + } + + + template inline void MGMatrixLocalBlocksToGlobalBlocks::initialize_edge_flux( diff --git a/deal.II/deal.II/include/numerics/mesh_worker_info.h b/deal.II/deal.II/include/numerics/mesh_worker_info.h index a2ceadb635..795410bece 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker_info.h +++ b/deal.II/deal.II/include/numerics/mesh_worker_info.h @@ -85,7 +85,7 @@ namespace MeshWorker * provide row and column info. */ template - void initialize_matrices(std::vector& matrices, + void initialize_matrices(const std::vector& matrices, bool both); /** @@ -526,7 +526,11 @@ namespace MeshWorker std::vector > >& 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 IntegrationInfoBox(const T&); - - template - void initialize_data(const VECTOR*, const std::string& name, - bool values, bool gradients, bool hessians); template void initialize(const WORKER&, const FiniteElement& el, const Mapping& mapping); + + template + void initialize(const WORKER&, + const FiniteElement& el, + const Mapping& mapping, + const NamedData& data); // private: - boost::shared_ptr > data; + boost::shared_ptr > cell_data; + boost::shared_ptr > bdry_data; + boost::shared_ptr > face_data; CellInfo cell_info; FaceInfo bdry_info; @@ -593,7 +601,7 @@ namespace MeshWorker template inline void LocalResults::initialize_matrices( - std::vector& matrices, + const std::vector& matrices, bool both) { M1.resize(matrices.size()); @@ -839,27 +847,6 @@ namespace MeshWorker {} - template - template - void - IntegrationInfoBox::initialize_data( - const VECTOR* v, const std::string& name, - bool values, bool gradients, bool hessians) - { - boost::shared_ptr > - p = boost::shared_ptr >(new VectorData ()); - 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 template void @@ -883,11 +870,41 @@ namespace MeshWorker subface_info.initialize >(el, mapping, integrator.face_quadrature, integrator.face_flags); neighbor_info.initialize >(el, mapping, integrator.face_quadrature, - integrator.ngbr_flags); - } -} + integrator.ngbr_flags); + } + template + template + void + IntegrationInfoBox::initialize( + const WORKER& integrator, + const FiniteElement& el, + const Mapping& mapping, + const NamedData& data) + { + initialize(integrator, el, mapping); + boost::shared_ptr > p; + + p = boost::shared_ptr >(new VectorData (integrator.cell_selector)); + p->initialize(data); + cell_data = p; + cell_info.initialize_data(p); + + p = boost::shared_ptr >(new VectorData (integrator.bdry_selector)); + p->initialize(data); + bdry_data = p; + bdry_info.initialize_data(p); + + p = boost::shared_ptr >(new VectorData (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 diff --git a/deal.II/deal.II/include/numerics/mesh_worker_info.templates.h b/deal.II/deal.II/include/numerics/mesh_worker_info.templates.h index 95c5274699..3e2da89c21 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker_info.templates.h +++ b/deal.II/deal.II/include/numerics/mesh_worker_info.templates.h @@ -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 + void + IntegrationInfo::initialize_data( + const boost::shared_ptr > 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;ilocal_indices().size()); + values[i].resize(n_components); +// deallog << ' ' << values[i].size() << " {"; // For all components for (unsigned int j=0;jn_gradients()); // For all selected finite // element functions for (unsigned int i=0;ilocal_indices().size()); + gradients[i].resize(n_components); // For all components for (unsigned int j=0;jlocal_indices().size()); + hessians[i].resize(n_components); // For all components for (unsigned int j=0;j - void - IntegrationInfo::initialize_data( - const boost::shared_ptr > data) - { - global_data = data; - } template 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 cab71251f3..fa4599bf89 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 @@ -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 void initialize(const NamedData&); @@ -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 { 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&); + + /** + * 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( 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 e691467949..6f1f855868 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 @@ -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::~VectorDataBase() {} + + template + VectorDataBase::VectorDataBase(const VectorSelector& v) + : + VectorSelector(v) + {} + + + template + VectorDataBase::VectorDataBase() + {} + + template void VectorDataBase::fill( @@ -38,6 +51,18 @@ namespace MeshWorker {} //----------------------------------------------------------------------// + template + VectorData::VectorData() + {} + + + template + VectorData::VectorData(const VectorSelector& s) + : + VectorDataBase(s) + {} + + template void VectorData::initialize(const NamedData& 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;in_values();++i) { const VECTOR& src = *data(this->value_index(i)); diff --git a/deal.II/deal.II/source/numerics/mesh_worker.cc b/deal.II/deal.II/source/numerics/mesh_worker.cc index bea9592615..a89e80acc9 100644 --- a/deal.II/deal.II/source/numerics/mesh_worker.cc +++ b/deal.II/deal.II/source/numerics/mesh_worker.cc @@ -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::IntegrationWorker () template void -IntegrationWorker::initialize_selectors( - const VectorSelector& cs, - const VectorSelector& bs, - const VectorSelector& fs) +IntegrationWorker::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::initialize_selectors( } -template -void -IntegrationWorker::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 void IntegrationWorker::add_update_flags( -- 2.39.5