From d9da6b6d8b09165893cd1f1599b7f7ee045884f4 Mon Sep 17 00:00:00 2001 From: bangerth Date: Wed, 20 Jan 2010 19:59:05 +0000 Subject: [PATCH] Get rid of the IntegratingAssembler class git-svn-id: https://svn.dealii.org/trunk@20411 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/include/numerics/mesh_worker.h | 150 +--------------- .../include/numerics/mesh_worker_info.h | 168 +++++++++--------- .../include/numerics/mesh_worker_loop.h | 13 +- deal.II/examples/step-38/step-38.cc | 31 ++-- 4 files changed, 109 insertions(+), 253 deletions(-) diff --git a/deal.II/deal.II/include/numerics/mesh_worker.h b/deal.II/deal.II/include/numerics/mesh_worker.h index 5a29fa3b73..c0ad600857 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker.h +++ b/deal.II/deal.II/include/numerics/mesh_worker.h @@ -345,113 +345,6 @@ namespace MeshWorker }; -/** - * A worker class for integration_loop(), separating assembling and - * local integration. Objects of this class rely on the base class - * provided by the template argument ASSEMBLER for assembling local - * data into global. The integration of local data is delegated to the - * other template parameter class INTEGRATOR. - * - * In order to use this class, it will be necessary to create an - * INTEGRATOR class following this template: - * - * @code - * template - * class MyLocalOperator : public Subscriptor - * { - * public: - * void cell(typename IntegrationWorker::CellInfo& info) const; - * void bdry(typename IntegrationWorker::FaceInfo& info) const; - * void face(typename IntegrationWorker::FaceInfo& info1, - * typename IntegrationWorker::FaceInfo& info2) const; - * }; - * @endcode - * - * This class will do whatever your problem requires locally on each - * cell and/or face. Once this class is defined, you choose a suitable - * assembler for your problem from the Assembler namespace and set up - * objects: - * - * @code - * MyLocalOperator myop; - * - * AssemblingIntegrator > - * integrator(myop); - * @endcode - * - * You do the necessary initializations of this @p integrator and then - * you have a worker object suitable for integration_loop(). - * - * @ingroup MeshWorker - * @author Guido Kanschat, 2009 - */ - template - class AssemblingIntegrator : public IntegrationWorker, - public ASSEMBLER - { - public: - typedef - typename IntegrationWorker::CellInfo - CellInfo; - - typedef - typename IntegrationWorker::FaceInfo - FaceInfo; - - /** - * Constructor, initializing - * the local data. - * - * The three arguments are - * objects that that will do - * the local computations on - * cells, boundary and interior - * faces. They may be specified - * as pointers to global or - * static functions, pointers - * to objects with an - * operator() that takes the - * required number of - * arguments, or functions with - * bound arguments. - */ - AssemblingIntegrator(const std_cxx1x::function &cell_worker, - const std_cxx1x::function &boundary_worker, - const std_cxx1x::function &face_worker); - - /** - * The cell operator called by - * integration_loop(). - */ - void cell(CellInfo& info); - - /** - * The local boundary operator - * called by - * integration_loop(). - */ - void bdry(FaceInfo& info); - - /** - * The interior face operator - * called by - * integration_loop(). - */ - void face(FaceInfo& info1, - FaceInfo& info2); - - private: - /** - * Pointers to the three - * functions that will do the - * local computations on cells, - * boundary and interior faces. - */ - const std_cxx1x::function cell_worker; - const std_cxx1x::function boundary_worker; - const std_cxx1x::function face_worker; - }; - //----------------------------------------------------------------------// template inline @@ -459,50 +352,9 @@ namespace MeshWorker : interior_fluxes(true), boundary_fluxes(true) {} - -//----------------------------------------------------------------------// - - template - inline - AssemblingIntegrator:: - AssemblingIntegrator(const std_cxx1x::function &cell_worker, - const std_cxx1x::function &boundary_worker, - const std_cxx1x::function &face_worker) - : - cell_worker (cell_worker), - boundary_worker (boundary_worker), - face_worker (face_worker) - {} - - - template - inline void - AssemblingIntegrator::cell(CellInfo& info) - { - cell_worker(info); - ASSEMBLER::assemble(info); - } - - - template - inline void - AssemblingIntegrator::bdry(FaceInfo& info) - { - boundary_worker(info); - ASSEMBLER::assemble(info); - } - - - template - inline void - AssemblingIntegrator::face(FaceInfo& info1, - FaceInfo& info2) - { - face_worker(info1, info2); - ASSEMBLER::assemble(info1, info2); - } } + DEAL_II_NAMESPACE_CLOSE #endif 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 853c7496e2..e6dc2a3c2c 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker_info.h +++ b/deal.II/deal.II/include/numerics/mesh_worker_info.h @@ -47,7 +47,7 @@ namespace MeshWorker void initialize_local(MatrixBlock >& M, const unsigned int row, const unsigned int col); - + public: void initialize_numbers(const unsigned int n); void initialize_vectors(const unsigned int n); @@ -67,7 +67,7 @@ namespace MeshWorker * provide row and column info. */ void initialize_matrices(unsigned int n, bool both); - + /** * Allocate a local matrix * for each of the global @@ -87,7 +87,7 @@ namespace MeshWorker template void initialize_matrices(const std::vector& matrices, bool both); - + /** * Reinitialize matrices for * new cell. Resizes the @@ -95,14 +95,14 @@ namespace MeshWorker * them to zero. */ void reinit(const BlockIndices& local_sizes); - + /** * The local numbers, * computed on a cell or on a * face. */ std::vector J; - + /** * The local vectors. This * field is public, so that @@ -110,7 +110,7 @@ namespace MeshWorker * write to it. */ std::vector > R; - + /** * The local matrices * coupling degrees of @@ -119,7 +119,7 @@ namespace MeshWorker * first cell on a face. */ std::vector > > M1; - + /** * The local matrices * coupling test functions on @@ -133,7 +133,7 @@ namespace MeshWorker std::vector > > M2; }; - + /** * Basic info class only containing information on geometry and * degrees of freedom of the mesh object. @@ -175,13 +175,13 @@ namespace MeshWorker template class DoFInfo : public LocalResults { - public: + public: /// The current cell typename Triangulation::cell_iterator cell; - + /// The current face typename Triangulation::face_iterator face; - + /** * The number of the current * face on the current cell. @@ -191,7 +191,7 @@ namespace MeshWorker * if the info object was * initialized with a cell. */ - + unsigned int face_number; /** * The number of the current @@ -207,7 +207,7 @@ namespace MeshWorker /// The DoF indices of the current cell std::vector indices; - + /** * Constructor setting the * #block_info pointer. @@ -222,7 +222,7 @@ namespace MeshWorker */ template DoFInfo(const DH& dof_handler); - + /** * Set the current cell and * fill #indices. @@ -239,7 +239,7 @@ namespace MeshWorker void reinit(const DHCellIterator& c, const DHFaceIterator& f, const unsigned int n); - + /** * Set the current subface * and fill #indices if the @@ -252,21 +252,21 @@ namespace MeshWorker const unsigned int s); const BlockIndices& local_indices() const; - - + + /// The block structure of the system SmartPointer > block_info; - + private: /// Fill index vector void get_indices(const typename DoFHandler::cell_iterator c); - + /// Fill index vector with level indices void get_indices(const typename MGDoFHandler::cell_iterator c); - + /// Auxiliary vector std::vector indices_org; - + /** * An auxiliary local * BlockIndices object created @@ -277,7 +277,7 @@ namespace MeshWorker */ BlockIndices aux_local_indices; }; - + /** * Class for objects handed to local integration functions. * @@ -340,14 +340,14 @@ namespace MeshWorker * information to DoFInfo. */ IntegrationInfo(const BlockInfo& block_info); - + /** * Constructor forwarding * information to DoFInfo. */ template IntegrationInfo(const DH& dof_handler); - + /** * Build all internal * structures, in particular @@ -384,12 +384,12 @@ namespace MeshWorker * selector. */ void initialize_data(const boost::shared_ptr > data); - + /** * Delete the data created by initialize(). */ void clear(); - + /// This is true if we are assembling for multigrid bool multigrid; /// Access to finite element @@ -403,7 +403,7 @@ namespace MeshWorker * vector of elements. */ const FEVALUESBASE& fe() const; - + /// Access to finite elements /** * This access function must @@ -414,7 +414,7 @@ namespace MeshWorker * @see DGBlockSplitApplication */ const FEVALUESBASE& fe(const unsigned int i) const; - + /** * The vector containing the * values of finite element @@ -430,7 +430,7 @@ namespace MeshWorker * point. */ std::vector > > values; - + /** * The vector containing the * derivatives of finite @@ -446,7 +446,7 @@ namespace MeshWorker * point. */ std::vector > > > gradients; - + /** * The vector containing the * second derivatives of finite @@ -462,7 +462,7 @@ namespace MeshWorker * point. */ std::vector > > > hessians; - + /** * Reinitialize internal data * structures for use on a cell. @@ -478,7 +478,7 @@ namespace MeshWorker void reinit(const DHCellIterator& c, const DHFaceIterator& f, const unsigned int fn); - + /** * Reinitialize internal data * structures for use on a subface. @@ -553,9 +553,10 @@ namespace MeshWorker */ template IntegrationInfoBox(const T&); - - template + + template void initialize(const WORKER&, + ASSEMBLER &assembler, const FiniteElement& el, const Mapping& mapping); @@ -564,20 +565,20 @@ namespace MeshWorker const FiniteElement& el, const Mapping& mapping, const NamedData& data); - + // private: boost::shared_ptr > cell_data; boost::shared_ptr > bdry_data; boost::shared_ptr > face_data; - + CellInfo cell_info; FaceInfo bdry_info; FaceInfo face_info; FaceInfo subface_info; FaceInfo neighbor_info; }; - + //----------------------------------------------------------------------// @@ -588,7 +589,7 @@ namespace MeshWorker J.resize(n); } - + template inline void LocalResults::initialize_vectors(const unsigned int n) @@ -596,7 +597,7 @@ namespace MeshWorker R.resize(n); } - + template template inline void @@ -621,8 +622,8 @@ namespace MeshWorker } } } - - + + template inline void LocalResults::initialize_matrices(const unsigned int n, @@ -642,8 +643,8 @@ namespace MeshWorker } } } - - + + template inline void LocalResults::reinit(const BlockIndices& bi) @@ -658,11 +659,11 @@ namespace MeshWorker for (unsigned int i=0;i template DoFInfo::DoFInfo(const DH& dof_handler) @@ -671,8 +672,8 @@ namespace MeshWorker aux[0] = dof_handler.get_fe().dofs_per_cell; aux_local_indices.reinit(aux); } - - + + template template inline void @@ -685,17 +686,17 @@ namespace MeshWorker if (block_info) LocalResults::reinit(block_info->local()); else - LocalResults::reinit(aux_local_indices); + LocalResults::reinit(aux_local_indices); } - - + + template template inline void DoFInfo::reinit(const DHCellIterator& c, const DHFaceIterator& f, unsigned int n) - { + { if ((cell.state() != IteratorState::valid) || cell != static_cast::cell_iterator> (c)) get_indices(c); @@ -706,10 +707,10 @@ namespace MeshWorker if (block_info) LocalResults::reinit(block_info->local()); else - LocalResults::reinit(aux_local_indices); + LocalResults::reinit(aux_local_indices); } - - + + template template inline void @@ -728,10 +729,10 @@ namespace MeshWorker if (block_info) LocalResults::reinit(block_info->local()); else - LocalResults::reinit(aux_local_indices); + LocalResults::reinit(aux_local_indices); } - - + + template inline const BlockIndices& DoFInfo::local_indices() const @@ -740,10 +741,10 @@ namespace MeshWorker return block_info->local(); return aux_local_indices; } - + //----------------------------------------------------------------------// - + template template IntegrationInfo::IntegrationInfo(const DH& dof_handler) @@ -754,7 +755,7 @@ namespace MeshWorker global_data(boost::shared_ptr >(new VectorDataBase)) {} - + template inline const FVB& IntegrationInfo::fe() const @@ -789,8 +790,8 @@ namespace MeshWorker const bool split_fevalues = this->block_info != 0; fill_local_data(split_fevalues); } - - + + template template inline void @@ -806,12 +807,12 @@ namespace MeshWorker FEFaceValues& fe = dynamic_cast&> (febase); fe.reinit(this->cell, fn); } - + const bool split_fevalues = this->block_info != 0; fill_local_data(split_fevalues); } - - + + template template inline void @@ -828,7 +829,7 @@ namespace MeshWorker FESubfaceValues& fe = dynamic_cast&> (febase); fe.reinit(this->cell, fn, sn); } - + const bool split_fevalues = this->block_info != 0; fill_local_data(split_fevalues); } @@ -845,22 +846,23 @@ namespace MeshWorker subface_info(t), neighbor_info(t) {} - - + + template - template + template void - IntegrationInfoBox::initialize( - const WORKER& integrator, - const FiniteElement& el, - const Mapping& mapping) + IntegrationInfoBox:: + initialize(const WORKER& integrator, + ASSEMBLER &assembler, + const FiniteElement& el, + const Mapping& mapping) { - integrator.initialize_info(cell_info, false); - integrator.initialize_info(bdry_info, false); - integrator.initialize_info(face_info, true); - integrator.initialize_info(subface_info, true); - integrator.initialize_info(neighbor_info, true); - + assembler.initialize_info(cell_info, false); + assembler.initialize_info(bdry_info, false); + assembler.initialize_info(face_info, true); + assembler.initialize_info(subface_info, true); + assembler.initialize_info(neighbor_info, true); + cell_info.initialize >(el, mapping, integrator.cell_quadrature, integrator.cell_flags); bdry_info.initialize >(el, mapping, integrator.bdry_quadrature, @@ -885,17 +887,17 @@ namespace MeshWorker { 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; diff --git a/deal.II/deal.II/include/numerics/mesh_worker_loop.h b/deal.II/deal.II/include/numerics/mesh_worker_loop.h index a4ea53d9d1..fda5bbb7ff 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker_loop.h +++ b/deal.II/deal.II/include/numerics/mesh_worker_loop.h @@ -66,7 +66,7 @@ namespace MeshWorker * @ingroup MeshWorker * @author Guido Kanschat, 2009 */ - template + template void loop(ITERATOR begin, typename identity::type end, CELLINFO& cell_info, FACEINFO& boundary_info, @@ -74,6 +74,7 @@ namespace MeshWorker const std_cxx1x::function &cell_worker, const std_cxx1x::function &boundary_worker, const std_cxx1x::function &face_worker, + ASSEMBLER &assembler, bool cells_first = true) { const bool integrate_cell = (cell_worker != 0); @@ -93,6 +94,7 @@ namespace MeshWorker { cell_info.reinit(cell); cell_worker(cell_info); + assembler.assemble (cell_info); } if (integrate_interior_face || integrate_boundary) @@ -108,6 +110,7 @@ namespace MeshWorker { boundary_info.reinit(cell, face, face_no); boundary_worker(boundary_info); + assembler.assemble (boundary_info); } } else if (integrate_interior_face) @@ -144,6 +147,7 @@ namespace MeshWorker // conform to // old version face_worker(face_info, subface_info); + assembler.assemble (face_info, subface_info); } else { @@ -168,9 +172,9 @@ namespace MeshWorker neighbor_face_no); face_worker(face_info, neighbor_info); + assembler.assemble (face_info, neighbor_info); neighbor->face(neighbor_face_no)->set_user_flag (); } - } } // faces // Execute this, if faces @@ -179,6 +183,7 @@ namespace MeshWorker { cell_info.reinit(cell); cell_worker(cell_info); + assembler.assemble (cell_info); } } } @@ -189,13 +194,14 @@ namespace MeshWorker * @ingroup MeshWorker * @author Guido Kanschat, 2009 */ - template + template void integration_loop(ITERATOR begin, typename identity::type end, IntegrationInfoBox& box, const std_cxx1x::function &cell_worker, const std_cxx1x::function &boundary_worker, const std_cxx1x::function &face_worker, + ASSEMBLER &assembler, bool cells_first = true) { loop @@ -206,6 +212,7 @@ namespace MeshWorker cell_worker, boundary_worker, face_worker, + assembler, cells_first); } } diff --git a/deal.II/examples/step-38/step-38.cc b/deal.II/examples/step-38/step-38.cc index ce44f58800..2a97fa4cda 100644 --- a/deal.II/examples/step-38/step-38.cc +++ b/deal.II/examples/step-38/step-38.cc @@ -327,15 +327,7 @@ void DGMethod::assemble_system () // result of std::bind if the local // integrators were, for example, // non-static member functions. - typedef - MeshWorker::AssemblingIntegrator - , - Vector > > - Integrator; - Integrator integrator(&DGMethod::integrate_cell_term, - &DGMethod::integrate_boundary_term, - &DGMethod::integrate_face_term); + MeshWorker::IntegrationWorker integration_worker; // First, we initialize the // quadrature formulae and the @@ -350,9 +342,9 @@ void DGMethod::assemble_system () // independently, we have to hand // over this value three times. const unsigned int n_gauss_points = dof_handler.get_fe().degree+1; - integrator.initialize_gauss_quadrature(n_gauss_points, - n_gauss_points, - n_gauss_points); + integration_worker.initialize_gauss_quadrature(n_gauss_points, + n_gauss_points, + n_gauss_points); // These are the types of values we // need for integrating our @@ -364,14 +356,16 @@ void DGMethod::assemble_system () UpdateFlags update_flags = update_quadrature_points | update_values | update_gradients; - integrator.add_update_flags(update_flags, true, true, true, true); + integration_worker.add_update_flags(update_flags, true, true, true, true); // Finally, we have to tell the // assembler base class where to // put the local data. These will // be our system matrix and the // right hand side. - integrator.initialize(system_matrix, right_hand_side); + MeshWorker::Assembler::SystemSimple, Vector > + assembler; + assembler.initialize(system_matrix, right_hand_side); // We are now ready to get to the // integration loop. @p info_box is @@ -384,7 +378,7 @@ void DGMethod::assemble_system () // receives all the stuff we // created so far. MeshWorker::IntegrationInfoBox info_box(dof_handler); - info_box.initialize(integrator, fe, mapping); + info_box.initialize(integration_worker, assembler, fe, mapping); // Finally, the integration loop // over all active cells @@ -393,9 +387,10 @@ void DGMethod::assemble_system () MeshWorker::integration_loop (dof_handler.begin_active(), dof_handler.end(), info_box, - std_cxx1x::bind (&Integrator::cell, integrator, _1), - std_cxx1x::bind (&Integrator::bdry, integrator, _1), - std_cxx1x::bind (&Integrator::face, integrator, _1, _2)); + &DGMethod::integrate_cell_term, + &DGMethod::integrate_boundary_term, + &DGMethod::integrate_face_term, + assembler); } -- 2.39.5