From f2abb6bc7c7ddb550545288a132969fb02e718fc Mon Sep 17 00:00:00 2001 From: kanschat Date: Wed, 16 Jun 2010 22:19:40 +0000 Subject: [PATCH] remove MeshWorker::IntegrationWorker git-svn-id: https://svn.dealii.org/trunk@21219 0785d39b-7218-0410-832d-ea1e28bc413d --- tests/deal.II/mesh_worker_01.cc | 32 ++++++++-------- tests/deal.II/mesh_worker_02.cc | 57 ++++++++++++++--------------- tests/multigrid/mg_renumbered_01.cc | 12 +++++- tests/multigrid/renumbering_01.cc | 3 +- tests/multigrid/renumbering_02.cc | 3 +- 5 files changed, 56 insertions(+), 51 deletions(-) diff --git a/tests/deal.II/mesh_worker_01.cc b/tests/deal.II/mesh_worker_01.cc index b4cbe552c5..dd241600a8 100644 --- a/tests/deal.II/mesh_worker_01.cc +++ b/tests/deal.II/mesh_worker_01.cc @@ -41,7 +41,7 @@ template class Local : public Subscriptor { public: - typedef typename MeshWorker::IntegrationWorker::CellInfo CellInfo; + typedef MeshWorker::IntegrationInfo CellInfo; void cell(MeshWorker::DoFInfo& dinfo, CellInfo& info) const; void bdry(MeshWorker::DoFInfo& dinfo, CellInfo& info) const; @@ -144,29 +144,27 @@ test_simple(MGDoFHandler& mgdofs) local.cells = true; local.faces = false; - MeshWorker::IntegrationWorker integrator; - integrator.initialize_gauss_quadrature(1, 1, 1); + MappingQ1 mapping; + + MeshWorker::IntegrationInfoBox info_box; + info_box.initialize_gauss_quadrature(1, 1, 1); + info_box.initialize_update_flags(); + info_box.initialize(fe, mapping); + + MeshWorker::DoFInfo dof_info(dofs); MeshWorker::Assembler::SystemSimple, Vector > assembler; assembler.initialize(matrix, v); - - { - MappingQ1 mapping; - - MeshWorker::IntegrationInfoBox info_box; - info_box.initialize(integrator, fe, mapping); - MeshWorker::DoFInfo dof_info(dofs); - - MeshWorker::integration_loop - (dofs.begin_active(), dofs.end(), - dof_info, info_box, - std_cxx1x::bind (&Local::cell, local, _1, _2), + + MeshWorker::integration_loop + (dofs.begin_active(), dofs.end(), + dof_info, info_box, + std_cxx1x::bind (&Local::cell, local, _1, _2), std_cxx1x::bind (&Local::bdry, local, _1, _2), std_cxx1x::bind (&Local::face, local, _1, _2, _3, _4), assembler, true); - } - + for (unsigned int i=0;i& dinfo, - typename MeshWorker::IntegrationWorker::CellInfo& info); + MeshWorker::IntegrationInfo& info); static void bdry(MeshWorker::DoFInfo& dinfo, - typename MeshWorker::IntegrationWorker::CellInfo& info); + MeshWorker::IntegrationInfo& info); static void face(MeshWorker::DoFInfo& dinfo1, MeshWorker::DoFInfo& dinfo2, - typename MeshWorker::IntegrationWorker::CellInfo& info1, - typename MeshWorker::IntegrationWorker::CellInfo& info2); + MeshWorker::IntegrationInfo& info1, + MeshWorker::IntegrationInfo& info2); }; template void MatrixIntegrator::cell(MeshWorker::DoFInfo& dinfo, - typename MeshWorker::IntegrationWorker::CellInfo& info) + MeshWorker::IntegrationInfo& info) { const FEValuesBase& fe = info.fe_values(); FullMatrix& local_matrix = dinfo.matrix(0).matrix; @@ -68,7 +68,7 @@ void MatrixIntegrator::cell(MeshWorker::DoFInfo& dinfo, template void MatrixIntegrator::bdry(MeshWorker::DoFInfo& dinfo, - typename MeshWorker::IntegrationWorker::CellInfo& info) + MeshWorker::IntegrationInfo& info) { const FEValuesBase& fe = info.fe_values(); FullMatrix& local_matrix = dinfo.matrix(0).matrix; @@ -89,8 +89,8 @@ void MatrixIntegrator::bdry(MeshWorker::DoFInfo& dinfo, template void MatrixIntegrator::face(MeshWorker::DoFInfo& dinfo1, MeshWorker::DoFInfo& dinfo2, - typename MeshWorker::IntegrationWorker::CellInfo& info1, - typename MeshWorker::IntegrationWorker::CellInfo& info2) + MeshWorker::IntegrationInfo& info1, + MeshWorker::IntegrationInfo& info2) { const FEValuesBase& fe1 = info1.fe_values(); const FEValuesBase& fe2 = info2.fe_values(); @@ -133,22 +133,21 @@ void assemble(const DoFHandler& dof_handler, SparseMatrix& matrix) { const FiniteElement& fe = dof_handler.get_fe(); - MappingQ1 mapping; - - MeshWorker::IntegrationWorker integration_worker; - MeshWorker::Assembler::MatrixSimple > assembler; + MappingQ1 mapping; + 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.initialize(integration_worker, fe, mapping); + 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(), @@ -171,19 +170,19 @@ assemble(const MGDoFHandler& dof_handler, const FiniteElement& fe = dof_handler.get_fe(); MappingQ1 mapping; - MeshWorker::IntegrationWorker integration_worker; - MeshWorker::Assembler::MGMatrixSimple > 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.initialize(integration_worker, fe, mapping); + info_box.add_update_flags(update_flags, true, true, true, true); + info_box.initialize(fe, mapping); + MeshWorker::DoFInfo dof_info(dof_handler); - + + MeshWorker::Assembler::MGMatrixSimple > assembler; + assembler.initialize(matrix); + MeshWorker::loop, MeshWorker::IntegrationInfoBox > (dof_handler.begin(), dof_handler.end(), diff --git a/tests/multigrid/mg_renumbered_01.cc b/tests/multigrid/mg_renumbered_01.cc index 5aa8fdc4f4..f6c6015790 100644 --- a/tests/multigrid/mg_renumbered_01.cc +++ b/tests/multigrid/mg_renumbered_01.cc @@ -1,4 +1,14 @@ -#define baerbel_mg_test 1 +//---------------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 2010 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------------------------------------------------------- #include #include diff --git a/tests/multigrid/renumbering_01.cc b/tests/multigrid/renumbering_01.cc index 50489f46f7..1622896689 100644 --- a/tests/multigrid/renumbering_01.cc +++ b/tests/multigrid/renumbering_01.cc @@ -1,8 +1,7 @@ //---------------------------------------------------------------------------- // $Id$ -// Version: $Name$ // -// Copyright (C) 2000 - 2006 by the deal.II authors +// Copyright (C) 2000 - 2006, 2010 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer diff --git a/tests/multigrid/renumbering_02.cc b/tests/multigrid/renumbering_02.cc index 0fd7541f13..70a829cc8e 100644 --- a/tests/multigrid/renumbering_02.cc +++ b/tests/multigrid/renumbering_02.cc @@ -1,8 +1,7 @@ //---------------------------------------------------------------------------- // $Id$ -// Version: $Name$ // -// Copyright (C) 2000 - 2006 by the deal.II authors +// Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 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 -- 2.39.5