From 4682e1f420b573308ace8fe5bc4533fbe97a1709 Mon Sep 17 00:00:00 2001 From: kanschat Date: Wed, 16 Jun 2010 22:15:20 +0000 Subject: [PATCH] remove MeshWorker::IntegrationWorker git-svn-id: https://svn.dealii.org/trunk@21218 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/include/numerics/mesh_worker.h | 168 ------------- .../include/numerics/mesh_worker_info.h | 220 +++++++++++++----- .../numerics/mesh_worker_info.templates.h | 81 +++++++ .../deal.II/source/numerics/mesh_worker.cc | 68 ------ .../source/numerics/mesh_worker_info.cc | 1 + deal.II/examples/step-12/step-12.cc | 55 ++--- deal.II/examples/step-39/step-39.cc | 132 +++++------ 7 files changed, 345 insertions(+), 380 deletions(-) diff --git a/deal.II/deal.II/include/numerics/mesh_worker.h b/deal.II/deal.II/include/numerics/mesh_worker.h index 135617028e..aca24c21bc 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker.h +++ b/deal.II/deal.II/include/numerics/mesh_worker.h @@ -174,174 +174,6 @@ template class MGDoFHandler; */ 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, #boundary_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, #boundary_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 - */ - template - class IntegrationWorker - { - public: - /** - * The info type expected by a - * cell integrator. - */ - typedef IntegrationInfo CellInfo; - - /** - * Initialize default values. - */ - IntegrationWorker(); - /** - * Initialize the - * VectorSelector objects - * #cell_selector, - * #boundary_selector and - * #face_selector in order to - * save computational - * eeffort. If no selectors - * are used, then values for - * all named vectors in - * DoFInfo::global_data will be - * computed in all quadrature - * points. - * - * This function will also - * add UpdateFlags to the - * flags stored in this class. - */ - void initialize_update_flags(); - - /** - * Add additional values for update. - */ - void add_update_flags(const UpdateFlags flags, bool cell = true, - bool boundary = true, bool face = true, - bool neighbor = true); - - /** Assign n-point Gauss - * quadratures to each of the - * quadrature rules. Here, a - * size of zero points means - * that no loop over these grid - * entities should be - * performed. - */ - void initialize_gauss_quadrature(unsigned int n_cell_points, - unsigned int n_boundary_points, - unsigned int n_face_points); - - /** - * Select the vectors from - * DoFInfo::global_data - * that should be computed in - * the quadrature points on cells. - */ - MeshWorker::VectorSelector cell_selector; - - /** - * Select the vectors from - * DoFInfo::global_data - * that should be computed in - * the quadrature points on - * boundary faces. - */ - MeshWorker::VectorSelector boundary_selector; - - /** - * Select the vectors from - * DoFInfo::global_data - * that should be computed in - * the quadrature points on - * interior faces. - */ - MeshWorker::VectorSelector face_selector; - - /** - * The set of update flags - * for boundary cell integration. - * - * Defaults to - * #update_JxW_values. - */ - UpdateFlags cell_flags; - /** - * The set of update flags - * for boundary face integration. - * - * Defaults to - * #update_JxW_values and - * #update_normal_vectors. - */ - UpdateFlags boundary_flags; - - /** - * The set of update flags - * for interior face integration. - * - * Defaults to - * #update_JxW_values and - * #update_normal_vectors. - */ - UpdateFlags face_flags; - - /** - * The set of update flags - * for interior face integration. - * - * Defaults to - * #update_default, since - * quadrature weights are - * taken from the other cell. - */ - UpdateFlags neighbor_flags; - - /** - * The quadrature rule used - * on cells. - */ - Quadrature cell_quadrature; - - /** - * The quadrature rule used - * on boundary faces. - */ - Quadrature boundary_quadrature; - - /** - * The quadrature rule used - * on interior faces. - */ - Quadrature face_quadrature; - }; } 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 dd854287bd..cfa71204d6 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker_info.h +++ b/deal.II/deal.II/include/numerics/mesh_worker_info.h @@ -726,6 +726,31 @@ namespace MeshWorker * functions at quadrature points. * * + * In order to allow for sufficient generality, a few steps have to be + * undertaken to use this class. + * + * First, you should consider if you need values from any vectors in a + * NamedData object. If so, fill the VectorSelector objects + * #cell_selector, #boundary_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, #boundary_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 */ @@ -737,29 +762,102 @@ namespace MeshWorker /// The type of the info object for cells typedef IntegrationInfo CellInfo; - /** - * @name FEValues setup - */ - /* @{ */ - template - void initialize(const WORKER&, - const FiniteElement& el, + void initialize(const FiniteElement& el, const Mapping& mapping, const BlockInfo* block_info = 0); - template - void initialize(const WORKER&, - const FiniteElement& el, + template + void initialize(const FiniteElement& el, const Mapping& mapping, const NamedData& data, const BlockInfo* block_info = 0); - template - void initialize(const WORKER&, - const FiniteElement& el, + template + void initialize(const FiniteElement& el, const Mapping& mapping, const NamedData*>& data, const BlockInfo* block_info = 0); + /** + * @name FEValues setup + */ + /* @{ */ + void initialize_update_flags(); + + /** + * Add additional values for update. + */ + void add_update_flags(const UpdateFlags flags, bool cell = true, + bool boundary = true, bool face = true, + bool neighbor = true); + + /** Assign n-point Gauss + * quadratures to each of the + * quadrature rules. Here, a + * size of zero points means + * that no loop over these grid + * entities should be + * performed. + */ + void initialize_gauss_quadrature(unsigned int n_cell_points, + unsigned int n_boundary_points, + unsigned int n_face_points); + + /** + * The set of update flags + * for boundary cell integration. + * + * Defaults to + * #update_JxW_values. + */ + UpdateFlags cell_flags; + /** + * The set of update flags + * for boundary face integration. + * + * Defaults to + * #update_JxW_values and + * #update_normal_vectors. + */ + UpdateFlags boundary_flags; + + /** + * The set of update flags + * for interior face integration. + * + * Defaults to + * #update_JxW_values and + * #update_normal_vectors. + */ + UpdateFlags face_flags; + + /** + * The set of update flags + * for interior face integration. + * + * Defaults to + * #update_default, since + * quadrature weights are + * taken from the other cell. + */ + UpdateFlags neighbor_flags; + + /** + * The quadrature rule used + * on cells. + */ + Quadrature cell_quadrature; + + /** + * The quadrature rule used + * on boundary faces. + */ + Quadrature boundary_quadrature; + + /** + * The quadrature rule used + * on interior faces. + */ + Quadrature face_quadrature; /* @} */ /** @@ -767,6 +865,50 @@ namespace MeshWorker */ /* @{ */ + /** + * Initialize the + * VectorSelector objects + * #cell_selector, + * #boundary_selector and + * #face_selector in order to + * save computational + * eeffort. If no selectors + * are used, then values for + * all named vectors in + * DoFInfo::global_data will be + * computed in all quadrature + * points. + * + * This function will also + * add UpdateFlags to the + * flags stored in this class. + */ + /** + * Select the vectors from + * DoFInfo::global_data + * that should be computed in + * the quadrature points on cells. + */ + MeshWorker::VectorSelector cell_selector; + + /** + * Select the vectors from + * DoFInfo::global_data + * that should be computed in + * the quadrature points on + * boundary faces. + */ + MeshWorker::VectorSelector boundary_selector; + + /** + * Select the vectors from + * DoFInfo::global_data + * that should be computed in + * the quadrature points on + * interior faces. + */ + MeshWorker::VectorSelector face_selector; + boost::shared_ptr > cell_data; boost::shared_ptr > boundary_data; boost::shared_ptr > face_data; @@ -1324,51 +1466,28 @@ namespace MeshWorker //----------------------------------------------------------------------// template - template - void - IntegrationInfoBox:: - initialize(const WORKER& integrator, - const FiniteElement& el, - const Mapping& mapping, - const BlockInfo* block_info) - { - cell.initialize >(el, mapping, integrator.cell_quadrature, - integrator.cell_flags, block_info); - boundary.initialize >(el, mapping, integrator.boundary_quadrature, - integrator.boundary_flags, block_info); - face.initialize >(el, mapping, integrator.face_quadrature, - integrator.face_flags, block_info); - subface.initialize >(el, mapping, integrator.face_quadrature, - integrator.face_flags, block_info); - neighbor.initialize >(el, mapping, integrator.face_quadrature, - integrator.neighbor_flags, block_info); - } - - - template - template + template void IntegrationInfoBox::initialize( - const WORKER& integrator, const FiniteElement& el, const Mapping& mapping, const NamedData& data, const BlockInfo* block_info) { - initialize(integrator, el, mapping, block_info); + initialize(el, mapping, block_info); boost::shared_ptr > p; - p = boost::shared_ptr >(new VectorData (integrator.cell_selector)); + p = boost::shared_ptr >(new VectorData (cell_selector)); p->initialize(data); cell_data = p; cell.initialize_data(p); - p = boost::shared_ptr >(new VectorData (integrator.boundary_selector)); + p = boost::shared_ptr >(new VectorData (boundary_selector)); p->initialize(data); boundary_data = p; boundary.initialize_data(p); - p = boost::shared_ptr >(new VectorData (integrator.face_selector)); + p = boost::shared_ptr >(new VectorData (face_selector)); p->initialize(data); face_data = p; face.initialize_data(p); @@ -1378,51 +1497,48 @@ namespace MeshWorker template - template + template void IntegrationInfoBox::initialize( - const WORKER& integrator, const FiniteElement& el, const Mapping& mapping, const NamedData*>& data, const BlockInfo* block_info) { - initialize(integrator, el, mapping, block_info); + initialize(el, mapping, block_info); boost::shared_ptr > p; - p = boost::shared_ptr >(new MGVectorData (integrator.cell_selector)); + p = boost::shared_ptr >(new MGVectorData (cell_selector)); p->initialize(data); cell_data = p; cell.initialize_data(p); - p = boost::shared_ptr >(new MGVectorData (integrator.boundary_selector)); + p = boost::shared_ptr >(new MGVectorData (boundary_selector)); p->initialize(data); boundary_data = p; boundary.initialize_data(p); - p = boost::shared_ptr >(new MGVectorData (integrator.face_selector)); + p = boost::shared_ptr >(new MGVectorData (face_selector)); p->initialize(data); face_data = p; face.initialize_data(p); subface.initialize_data(p); neighbor.initialize_data(p); } - - + + template template void IntegrationInfoBox::post_cell(const DoFInfoBox&) {} - + template template void IntegrationInfoBox::post_faces(const DoFInfoBox&) {} - - } DEAL_II_NAMESPACE_CLOSE 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 409643d926..493e361660 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 @@ -13,6 +13,7 @@ #include +#include DEAL_II_NAMESPACE_OPEN @@ -220,6 +221,86 @@ namespace MeshWorker 0, n_comp, 0, info.indices.size()); } } + + +//----------------------------------------------------------------------// + + + template + void + IntegrationInfoBox::initialize_update_flags () + { + cell_flags = update_JxW_values; + boundary_flags = UpdateFlags(update_JxW_values | update_normal_vectors); + face_flags = boundary_flags; + neighbor_flags = update_default; + + 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 (boundary_selector.has_values() != 0) boundary_flags |= update_values; + if (boundary_selector.has_gradients() != 0) boundary_flags |= update_gradients; + if (boundary_selector.has_hessians() != 0) boundary_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) neighbor_flags |= update_values; + if (face_selector.has_gradients() != 0) neighbor_flags |= update_gradients; + if (face_selector.has_hessians() != 0) neighbor_flags |= update_hessians; + } + + + template + void + IntegrationInfoBox::add_update_flags( + const UpdateFlags flags, + bool cell, + bool boundary, + bool face, + bool neighbor) + { + if (cell) cell_flags |= flags; + if (boundary) boundary_flags |= flags; + if (face) face_flags |= flags; + if (neighbor) neighbor_flags |= flags; + } + + + template + void + IntegrationInfoBox::initialize_gauss_quadrature( + unsigned int cp, + unsigned int bp, + unsigned int fp) + { + cell_quadrature = QGauss(cp); + boundary_quadrature = QGauss(bp); + face_quadrature = QGauss(fp); + + } + + + template + void + IntegrationInfoBox:: + initialize(const FiniteElement& el, + const Mapping& mapping, + const BlockInfo* block_info) + { + cell.initialize >(el, mapping, cell_quadrature, + cell_flags, block_info); + boundary.initialize >(el, mapping, boundary_quadrature, + boundary_flags, block_info); + face.initialize >(el, mapping, face_quadrature, + face_flags, block_info); + subface.initialize >(el, mapping, face_quadrature, + face_flags, block_info); + neighbor.initialize >(el, mapping, face_quadrature, + neighbor_flags, block_info); + } } diff --git a/deal.II/deal.II/source/numerics/mesh_worker.cc b/deal.II/deal.II/source/numerics/mesh_worker.cc index b430dd3230..803c7e5a83 100644 --- a/deal.II/deal.II/source/numerics/mesh_worker.cc +++ b/deal.II/deal.II/source/numerics/mesh_worker.cc @@ -18,72 +18,4 @@ DEAL_II_NAMESPACE_OPEN -namespace MeshWorker -{ - - template - IntegrationWorker::IntegrationWorker () - { - cell_flags = update_JxW_values; - boundary_flags = UpdateFlags(update_JxW_values | update_normal_vectors); - face_flags = boundary_flags; - neighbor_flags = update_default; - } - - - template - void - IntegrationWorker::initialize_update_flags () - { - 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 (boundary_selector.has_values() != 0) boundary_flags |= update_values; - if (boundary_selector.has_gradients() != 0) boundary_flags |= update_gradients; - if (boundary_selector.has_hessians() != 0) boundary_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) neighbor_flags |= update_values; - if (face_selector.has_gradients() != 0) neighbor_flags |= update_gradients; - if (face_selector.has_hessians() != 0) neighbor_flags |= update_hessians; - } - - - template - void - IntegrationWorker::add_update_flags( - const UpdateFlags flags, bool cell, bool boundary, bool face, bool neighbor) - { - if (cell) cell_flags |= flags; - if (boundary) boundary_flags |= flags; - if (face) face_flags |= flags; - if (neighbor) neighbor_flags |= flags; - } - - - template - void - IntegrationWorker::initialize_gauss_quadrature( - unsigned int cp, - unsigned int bp, - unsigned int fp) - { - cell_quadrature = QGauss(cp); - boundary_quadrature = QGauss(bp); - face_quadrature = QGauss(fp); - - } -} - -#if deal_II_dimension > 1 -namespace MeshWorker -{ - template class IntegrationWorker; -} -#endif - DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/deal.II/source/numerics/mesh_worker_info.cc b/deal.II/deal.II/source/numerics/mesh_worker_info.cc index c70e915fe4..9943e4c195 100644 --- a/deal.II/deal.II/source/numerics/mesh_worker_info.cc +++ b/deal.II/deal.II/source/numerics/mesh_worker_info.cc @@ -27,6 +27,7 @@ DEAL_II_NAMESPACE_OPEN namespace MeshWorker { template class IntegrationInfo; + template class IntegrationInfoBox; template class LocalResults; template class DoFInfo; diff --git a/deal.II/examples/step-12/step-12.cc b/deal.II/examples/step-12/step-12.cc index 42c8692ba1..bd07fd93a6 100644 --- a/deal.II/examples/step-12/step-12.cc +++ b/deal.II/examples/step-12/step-12.cc @@ -199,7 +199,7 @@ class Step12 // order to make our life easier // below. typedef MeshWorker::DoFInfo DoFInfo; - typedef typename MeshWorker::IntegrationWorker::CellInfo CellInfo; + typedef MeshWorker::IntegrationInfo CellInfo; // The following three functions // are then the ones that get called @@ -309,8 +309,8 @@ void Step12::assemble_system () // object distributes these into // the global sparse matrix and the // right hand side vector. - MeshWorker::IntegrationWorker integration_worker; - + MeshWorker::IntegrationInfoBox info_box; + // First, we initialize the // quadrature formulae and the // update flags in the worker base @@ -324,45 +324,46 @@ void Step12::assemble_system () // independently, we have to hand // over this value three times. const unsigned int n_gauss_points = dof_handler.get_fe().degree+1; - integration_worker.initialize_gauss_quadrature(n_gauss_points, + info_box.initialize_gauss_quadrature(n_gauss_points, n_gauss_points, n_gauss_points); - + // These are the types of values we // need for integrating our // system. They are added to the // flags used on cells, boundary // and interior faces, as well as // interior neighbor faces, which is - // forced by the four @p true values. + // forced by the four @p true + // values. + info_box.initialize_update_flags(); UpdateFlags update_flags = update_quadrature_points | update_values | update_gradients; - integration_worker.add_update_flags(update_flags, true, true, true, true); + info_box.add_update_flags(update_flags, true, true, true, true); + + // After preparing all data in + // info_box, we initialize + // the FEValus objects in there. + info_box.initialize(fe, mapping); + + // The object created so far helps + // us do the local integration on + // each cell and face. Now, we need + // an object which receives the + // integrated (local) data and + // forwards them to the assembler. + MeshWorker::DoFInfo dof_info(dof_handler); - // 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. + // Now, we have to create the + // assembler object and tell it, + // where to put the local + // data. These will be our system + // matrix and the 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 - // an object that generates the - // extended iterators for cells and - // faces of type - // MeshWorker::IntegrationInfo. Since - // we need five different of them, - // this is a handy shortcut. It - // receives all the stuff we - // created so far. - MeshWorker::IntegrationInfoBox info_box; - MeshWorker::DoFInfo dof_info(dof_handler); - info_box.initialize(integration_worker, fe, mapping); - + // Finally, the integration loop // over all active cells // (determined by the first diff --git a/deal.II/examples/step-39/step-39.cc b/deal.II/examples/step-39/step-39.cc index a686a13596..8563ea4101 100644 --- a/deal.II/examples/step-39/step-39.cc +++ b/deal.II/examples/step-39/step-39.cc @@ -87,7 +87,7 @@ Functions::SlitSingularityFunction<2> exact_solution; // interior faces, respectively. All // the information needed for the // local integration is provided by - // MeshWorker::IntegrationWorker::CellInfo. Note + // MeshWorker::IntegrationInfo. Note // that this public interface cannot // be changed, because it is expected // by MeshWorker::integration_loop(). @@ -96,13 +96,13 @@ class MatrixIntegrator : public Subscriptor { public: static void cell(MeshWorker::DoFInfo& dinfo, - typename MeshWorker::IntegrationWorker::CellInfo& info); + typename MeshWorker::IntegrationInfo& info); static void bdry(MeshWorker::DoFInfo& dinfo, - typename MeshWorker::IntegrationWorker::CellInfo& info); + typename MeshWorker::IntegrationInfo& info); static void face(MeshWorker::DoFInfo& dinfo1, MeshWorker::DoFInfo& dinfo2, - typename MeshWorker::IntegrationWorker::CellInfo& info1, - typename MeshWorker::IntegrationWorker::CellInfo& info2); + typename MeshWorker::IntegrationInfo& info1, + typename MeshWorker::IntegrationInfo& info2); }; @@ -116,7 +116,7 @@ class MatrixIntegrator : public Subscriptor template void MatrixIntegrator::cell( MeshWorker::DoFInfo& dinfo, - typename MeshWorker::IntegrationWorker::CellInfo& info) + typename MeshWorker::IntegrationInfo& info) { const FEValuesBase& fe = info.fe_values(); FullMatrix& local_matrix = dinfo.matrix(0,false).matrix; @@ -133,7 +133,7 @@ void MatrixIntegrator::cell( template void MatrixIntegrator::bdry( MeshWorker::DoFInfo& dinfo, - typename MeshWorker::IntegrationWorker::CellInfo& info) + typename MeshWorker::IntegrationInfo& info) { const FEValuesBase& fe = info.fe_values(); FullMatrix& local_matrix = dinfo.matrix(0,false).matrix; @@ -155,8 +155,8 @@ template void MatrixIntegrator::face( MeshWorker::DoFInfo& dinfo1, MeshWorker::DoFInfo& dinfo2, - typename MeshWorker::IntegrationWorker::CellInfo& info1, - typename MeshWorker::IntegrationWorker::CellInfo& info2) + typename MeshWorker::IntegrationInfo& info1, + typename MeshWorker::IntegrationInfo& info2) { const FEValuesBase& fe1 = info1.fe_values(); const FEValuesBase& fe2 = info2.fe_values(); @@ -205,22 +205,22 @@ template class RHSIntegrator : public Subscriptor { public: - static void cell(MeshWorker::DoFInfo& dinfo, typename MeshWorker::IntegrationWorker::CellInfo& info); - static void bdry(MeshWorker::DoFInfo& dinfo, typename MeshWorker::IntegrationWorker::CellInfo& info); + static void cell(MeshWorker::DoFInfo& dinfo, typename MeshWorker::IntegrationInfo& info); + static void bdry(MeshWorker::DoFInfo& dinfo, typename MeshWorker::IntegrationInfo& info); static void face(MeshWorker::DoFInfo& dinfo1, MeshWorker::DoFInfo& dinfo2, - typename MeshWorker::IntegrationWorker::CellInfo& info1, - typename MeshWorker::IntegrationWorker::CellInfo& info2); + typename MeshWorker::IntegrationInfo& info1, + typename MeshWorker::IntegrationInfo& info2); }; template -void RHSIntegrator::cell(MeshWorker::DoFInfo&, typename MeshWorker::IntegrationWorker::CellInfo&) +void RHSIntegrator::cell(MeshWorker::DoFInfo&, typename MeshWorker::IntegrationInfo&) {} template -void RHSIntegrator::bdry(MeshWorker::DoFInfo& dinfo, typename MeshWorker::IntegrationWorker::CellInfo& info) +void RHSIntegrator::bdry(MeshWorker::DoFInfo& dinfo, typename MeshWorker::IntegrationInfo& info) { const FEValuesBase& fe = info.fe_values(); Vector& local_vector = dinfo.vector(0).block(0); @@ -242,8 +242,8 @@ void RHSIntegrator::bdry(MeshWorker::DoFInfo& dinfo, typename MeshWork template void RHSIntegrator::face(MeshWorker::DoFInfo&, MeshWorker::DoFInfo&, - typename MeshWorker::IntegrationWorker::CellInfo&, - typename MeshWorker::IntegrationWorker::CellInfo&) + typename MeshWorker::IntegrationInfo&, + typename MeshWorker::IntegrationInfo&) {} @@ -251,17 +251,17 @@ template class Estimator : public Subscriptor { public: - static void cell(MeshWorker::DoFInfo& dinfo, typename MeshWorker::IntegrationWorker::CellInfo& info); - static void bdry(MeshWorker::DoFInfo& dinfo, typename MeshWorker::IntegrationWorker::CellInfo& info); + static void cell(MeshWorker::DoFInfo& dinfo, typename MeshWorker::IntegrationInfo& info); + static void bdry(MeshWorker::DoFInfo& dinfo, typename MeshWorker::IntegrationInfo& info); static void face(MeshWorker::DoFInfo& dinfo1, MeshWorker::DoFInfo& dinfo2, - typename MeshWorker::IntegrationWorker::CellInfo& info1, - typename MeshWorker::IntegrationWorker::CellInfo& info2); + typename MeshWorker::IntegrationInfo& info1, + typename MeshWorker::IntegrationInfo& info2); }; template -void Estimator::cell(MeshWorker::DoFInfo& dinfo, typename MeshWorker::IntegrationWorker::CellInfo& info) +void Estimator::cell(MeshWorker::DoFInfo& dinfo, typename MeshWorker::IntegrationInfo& info) { const FEValuesBase& fe = info.fe_values(); @@ -276,7 +276,7 @@ void Estimator::cell(MeshWorker::DoFInfo& dinfo, typename MeshWorker:: template -void Estimator::bdry(MeshWorker::DoFInfo& dinfo, typename MeshWorker::IntegrationWorker::CellInfo& info) +void Estimator::bdry(MeshWorker::DoFInfo& dinfo, typename MeshWorker::IntegrationInfo& info) { const FEValuesBase& fe = info.fe_values(); @@ -298,8 +298,8 @@ void Estimator::bdry(MeshWorker::DoFInfo& dinfo, typename MeshWorker:: template void Estimator::face(MeshWorker::DoFInfo& dinfo1, MeshWorker::DoFInfo& dinfo2, - typename MeshWorker::IntegrationWorker::CellInfo& info1, - typename MeshWorker::IntegrationWorker::CellInfo& info2) + typename MeshWorker::IntegrationInfo& info1, + typename MeshWorker::IntegrationInfo& info2) { const FEValuesBase& fe = info1.fe_values(); const std::vector& uh1 = info1.values[0][0]; @@ -330,7 +330,7 @@ template class Step39 { public: - typedef typename MeshWorker::IntegrationWorker::CellInfo CellInfo; + typedef typename MeshWorker::IntegrationInfo CellInfo; Step39(const FiniteElement& fe); @@ -434,18 +434,19 @@ template void Step39::assemble_matrix() { - MeshWorker::IntegrationWorker integration_worker; - MeshWorker::Assembler::MatrixSimple > assembler; - + MeshWorker::IntegrationInfoBox info_box; const unsigned int n_gauss_points = dof_handler.get_fe().tensor_degree()+1; - integration_worker.initialize_gauss_quadrature(n_gauss_points, n_gauss_points, n_gauss_points); + info_box.initialize_gauss_quadrature(n_gauss_points, n_gauss_points, n_gauss_points); + info_box.initialize_update_flags(); UpdateFlags update_flags = update_values | update_gradients; - integration_worker.add_update_flags(update_flags, true, true, true, true); - - assembler.initialize(matrix); - MeshWorker::IntegrationInfoBox info_box; + info_box.add_update_flags(update_flags, true, true, true, true); + info_box.initialize(fe, mapping); + MeshWorker::DoFInfo dof_info(dof_handler); - info_box.initialize(integration_worker, fe, mapping); + + MeshWorker::Assembler::MatrixSimple > assembler; + assembler.initialize(matrix); + MeshWorker::integration_loop( dof_handler.begin_active(), dof_handler.end(), dof_info, info_box, @@ -460,19 +461,20 @@ template void Step39::assemble_mg_matrix() { - MeshWorker::IntegrationWorker integration_worker; - MeshWorker::Assembler::MGMatrixSimple > assembler; - + MeshWorker::IntegrationInfoBox info_box; const unsigned int n_gauss_points = mg_dof_handler.get_fe().tensor_degree()+1; - integration_worker.initialize_gauss_quadrature(n_gauss_points, n_gauss_points, n_gauss_points); + info_box.initialize_gauss_quadrature(n_gauss_points, n_gauss_points, n_gauss_points); + info_box.initialize_update_flags(); UpdateFlags update_flags = update_values | update_gradients; - integration_worker.add_update_flags(update_flags, true, true, true, true); + info_box.add_update_flags(update_flags, true, true, true, true); + info_box.initialize(fe, mapping); + MeshWorker::DoFInfo dof_info(mg_dof_handler); + + MeshWorker::Assembler::MGMatrixSimple > assembler; assembler.initialize(mg_matrix); assembler.initialize_fluxes(mg_matrix_dg_up, mg_matrix_dg_down); - MeshWorker::IntegrationInfoBox info_box; - MeshWorker::DoFInfo dof_info(mg_dof_handler); - info_box.initialize(integration_worker, fe, mapping); + MeshWorker::integration_loop ( mg_dof_handler.begin(), mg_dof_handler.end(), dof_info, info_box, @@ -487,21 +489,21 @@ template void Step39::assemble_right_hand_side() { - MeshWorker::IntegrationWorker integration_worker; - MeshWorker::Assembler::ResidualSimple > assembler; - + MeshWorker::IntegrationInfoBox info_box; const unsigned int n_gauss_points = dof_handler.get_fe().tensor_degree()+1; - integration_worker.initialize_gauss_quadrature(n_gauss_points, n_gauss_points+1, n_gauss_points); + info_box.initialize_gauss_quadrature(n_gauss_points, n_gauss_points+1, n_gauss_points); + info_box.initialize_update_flags(); UpdateFlags update_flags = update_quadrature_points | update_values | update_gradients; - integration_worker.add_update_flags(update_flags, true, true, true, true); - + info_box.add_update_flags(update_flags, true, true, true, true); + info_box.initialize(fe, mapping); + + MeshWorker::DoFInfo dof_info(dof_handler); + + MeshWorker::Assembler::ResidualSimple > assembler; NamedData* > data; Vector* rhs = &right_hand_side; data.add(rhs, "RHS"); assembler.initialize(data); - MeshWorker::IntegrationInfoBox info_box; - MeshWorker::DoFInfo dof_info(dof_handler); - info_box.initialize(integration_worker, fe, mapping); MeshWorker::integration_loop( dof_handler.begin_active(), dof_handler.end(), @@ -597,11 +599,9 @@ Step39::estimate() cell != triangulation.end();++cell,++i) cell->set_user_index(i); - MeshWorker::IntegrationWorker integration_worker; - MeshWorker::Assembler::CellsAndFaces assembler; - + MeshWorker::IntegrationInfoBox info_box; const unsigned int n_gauss_points = dof_handler.get_fe().tensor_degree()+1; - integration_worker.initialize_gauss_quadrature(n_gauss_points, n_gauss_points+1, n_gauss_points); + info_box.initialize_gauss_quadrature(n_gauss_points, n_gauss_points+1, n_gauss_points); NamedData* > solution_data; solution_data.add(&solution, "solution"); @@ -611,19 +611,21 @@ Step39::estimate() cs.add("solution", true, true, true); fs.add("solution", true, true, false); - integration_worker.cell_selector = cs; - integration_worker.boundary_selector = fs; - integration_worker.face_selector = fs; - integration_worker.initialize_update_flags(); - integration_worker.add_update_flags(update_quadrature_points, false, true, false, false); - + info_box.cell_selector = cs; + info_box.boundary_selector = fs; + info_box.face_selector = fs; + info_box.initialize_update_flags(); + info_box.add_update_flags(update_quadrature_points, false, true, false, false); + info_box.initialize(fe, mapping, solution_data); + + MeshWorker::DoFInfo dof_info(dof_handler); + + MeshWorker::Assembler::CellsAndFaces assembler; NamedData* > out_data; BlockVector* est = &estimates; out_data.add(est, "cells"); assembler.initialize(out_data, false); - MeshWorker::IntegrationInfoBox info_box; - MeshWorker::DoFInfo dof_info(dof_handler); - info_box.initialize(integration_worker, fe, mapping, solution_data); + MeshWorker::integration_loop ( dof_handler.begin_active(), dof_handler.end(), dof_info, info_box, -- 2.39.5