From: Guido Kanschat Date: Mon, 15 Feb 2010 20:58:04 +0000 (+0000) Subject: apply changes in library interface X-Git-Tag: v8.0.0~6428 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=59a4cf5c78c4d8197c1ba7666b1576f8370efd61;p=dealii.git apply changes in library interface git-svn-id: https://svn.dealii.org/trunk@20640 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/tests/deal.II/mesh_worker_02.cc b/tests/deal.II/mesh_worker_02.cc index 797ae07062..9df0556486 100644 --- a/tests/deal.II/mesh_worker_02.cc +++ b/tests/deal.II/mesh_worker_02.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2000, 2001, 2003, 2004, 2007, 2008, 2009 by the deal.II authors +// Copyright (C) 2000, 2001, 2003, 2004, 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 @@ -40,18 +40,23 @@ template class MatrixIntegrator : public Subscriptor { public: - void cell(typename MeshWorker::IntegrationWorker::CellInfo& info) const; - void bdry(typename MeshWorker::IntegrationWorker::FaceInfo& info) const; - void face(typename MeshWorker::IntegrationWorker::FaceInfo& info1, - typename MeshWorker::IntegrationWorker::FaceInfo& info2) const; + static void cell(MeshWorker::DoFInfo& dinfo, + typename MeshWorker::IntegrationWorker::CellInfo& info); + static void bdry(MeshWorker::DoFInfo& dinfo, + typename MeshWorker::IntegrationWorker::FaceInfo& info); + static void face(MeshWorker::DoFInfo& dinfo1, + MeshWorker::DoFInfo& dinfo2, + typename MeshWorker::IntegrationWorker::FaceInfo& info1, + typename MeshWorker::IntegrationWorker::FaceInfo& info2); }; template -void MatrixIntegrator::cell(typename MeshWorker::IntegrationWorker::CellInfo& info) const +void MatrixIntegrator::cell(MeshWorker::DoFInfo& dinfo, + typename MeshWorker::IntegrationWorker::CellInfo& info) { - const FEValuesBase& fe = info.fe(); - FullMatrix& local_matrix = info.M1[0].matrix; + const FEValuesBase& fe = info.fe_values(); + FullMatrix& local_matrix = dinfo.matrix(0).matrix; for (unsigned int k=0; k::cell(typename MeshWorker::IntegrationWorker::Ce template -void MatrixIntegrator::bdry(typename MeshWorker::IntegrationWorker::FaceInfo& info) const +void MatrixIntegrator::bdry(MeshWorker::DoFInfo& dinfo, + typename MeshWorker::IntegrationWorker::FaceInfo& info) { - const FEFaceValuesBase& fe = info.fe(); - FullMatrix& local_matrix = info.M1[0].matrix; + const FEFaceValuesBase& fe = info.fe_values(); + FullMatrix& local_matrix = dinfo.matrix(0).matrix; const unsigned int deg = fe.get_fe().tensor_degree(); - const double penalty = 2. * deg * (deg+1) * info.face->measure() / info.cell->measure(); + const double penalty = 2. * deg * (deg+1) * dinfo.face->measure() / dinfo.cell->measure(); for (unsigned k=0;k::bdry(typename MeshWorker::IntegrationWorker::Fa template -void MatrixIntegrator::face(typename MeshWorker::IntegrationWorker::FaceInfo& info1, - typename MeshWorker::IntegrationWorker::FaceInfo& info2) const +void MatrixIntegrator::face(MeshWorker::DoFInfo& dinfo1, + MeshWorker::DoFInfo& dinfo2, + typename MeshWorker::IntegrationWorker::FaceInfo& info1, + typename MeshWorker::IntegrationWorker::FaceInfo& info2) { - const FEFaceValuesBase& fe1 = info1.fe(); - const FEFaceValuesBase& fe2 = info2.fe(); - FullMatrix& matrix_v1u1 = info1.M1[0].matrix; - FullMatrix& matrix_v1u2 = info1.M2[0].matrix; - FullMatrix& matrix_v2u1 = info2.M2[0].matrix; - FullMatrix& matrix_v2u2 = info2.M1[0].matrix; + const FEFaceValuesBase& fe1 = info1.fe_values(); + const FEFaceValuesBase& fe2 = info2.fe_values(); + FullMatrix& matrix_v1u1 = dinfo1.matrix(0, false).matrix; + FullMatrix& matrix_v1u2 = dinfo1.matrix(0, true).matrix; + FullMatrix& matrix_v2u1 = dinfo2.matrix(0, true).matrix; + FullMatrix& matrix_v2u2 = dinfo2.matrix(0, false).matrix; const unsigned int deg = fe1.get_fe().tensor_degree(); - const double penalty1 = deg * (deg+1) * info1.face->measure() / info1.cell->measure(); - const double penalty2 = deg * (deg+1) * info2.face->measure() / info2.cell->measure(); + const double penalty1 = deg * (deg+1) * dinfo1.face->measure() / dinfo1.cell->measure(); + const double penalty2 = deg * (deg+1) * dinfo2.face->measure() / dinfo2.cell->measure(); const double penalty = penalty1 + penalty2; for (unsigned k=0;k& dof_handler, SparseMatrix& matrix) const FiniteElement& fe = dof_handler.get_fe(); MappingQ1 mapping; - const MatrixIntegrator local; - MeshWorker::AssemblingIntegrator >, MatrixIntegrator > - integrator(local); + MeshWorker::IntegrationWorker integration_worker; + MeshWorker::Assembler::MatrixSimple > assembler; + const unsigned int n_gauss_points = dof_handler.get_fe().tensor_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); UpdateFlags update_flags = 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); - integrator.initialize(matrix); - MeshWorker::IntegrationInfoBox info_box(dof_handler); - info_box.initialize(integrator, fe, mapping); - MeshWorker::integration_loop(dof_handler.begin_active(), dof_handler.end(), info_box, integrator); + assembler.initialize(matrix); + MeshWorker::IntegrationInfoBox info_box; + info_box.initialize(integration_worker, fe, mapping); + MeshWorker::DoFInfo dof_info(dof_handler); + + info_box.initialize(integration_worker, fe, mapping); + MeshWorker::loop + , MeshWorker::IntegrationInfoBox > + (dof_handler.begin_active(), + dof_handler.end(), + dof_info, + info_box, + &MatrixIntegrator::cell, + &MatrixIntegrator::bdry, + &MatrixIntegrator::face, + assembler); } @@ -152,19 +172,28 @@ assemble(const MGDoFHandler& dof_handler, const FiniteElement& fe = dof_handler.get_fe(); MappingQ1 mapping; - const MatrixIntegrator local; - MeshWorker::AssemblingIntegrator >, MatrixIntegrator > - integrator(local); + MeshWorker::IntegrationWorker integration_worker; + MeshWorker::Assembler::MGMatrixSimple > assembler; + const unsigned int n_gauss_points = dof_handler.get_fe().tensor_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); UpdateFlags update_flags = 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); + + assembler.initialize(matrix); + MeshWorker::IntegrationInfoBox info_box; + info_box.initialize(integration_worker, fe, mapping); + MeshWorker::DoFInfo dof_info(dof_handler); - integrator.initialize(matrix); - integrator.initialize_fluxes(dg_up, dg_down); - MeshWorker::IntegrationInfoBox info_box(dof_handler); - info_box.initialize(integrator, fe, mapping); - MeshWorker::integration_loop(dof_handler.begin_active(), dof_handler.end(), info_box, integrator); + MeshWorker::loop, MeshWorker::IntegrationInfoBox > + (dof_handler.begin(), + dof_handler.end(), + dof_info, + info_box, + &MatrixIntegrator::cell, + &MatrixIntegrator::bdry, + &MatrixIntegrator::face, + assembler); }