From 862569ee640e654d48191aa27b953fb09b60bd63 Mon Sep 17 00:00:00 2001 From: kanschat Date: Wed, 29 Sep 2010 00:42:56 +0000 Subject: [PATCH] start making MGLocalBlocksToGlobalBlocks work git-svn-id: https://svn.dealii.org/trunk@22187 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/include/numerics/mesh_worker.h | 47 ++++++++++-- .../include/numerics/mesh_worker_assembler.h | 73 ++++++++++--------- 2 files changed, 81 insertions(+), 39 deletions(-) diff --git a/deal.II/deal.II/include/numerics/mesh_worker.h b/deal.II/deal.II/include/numerics/mesh_worker.h index ee1bfe6213..2321eb8f58 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker.h +++ b/deal.II/deal.II/include/numerics/mesh_worker.h @@ -224,15 +224,23 @@ namespace MeshWorker * column coordinates. The * matrices themselves are * resized by reinit(). - * - * The template parameter @p - * MatrixPtr should point to - * a MatrixBlock - * instantiation in order to - * provide row and column info. */ template void initialize_matrices(const MatrixBlockVector& matrices, + bool both); + + /** + * Allocate a local matrix + * for each of the global + * level objects in @p + * matrices. Additionally, + * set their block row and + * column coordinates. The + * matrices themselves are + * resized by reinit(). + */ + template + void initialize_matrices(const MGMatrixBlockVector& matrices, bool both); /** @@ -426,6 +434,33 @@ namespace MeshWorker } + template + template + inline void + LocalResults::initialize_matrices( + const MGMatrixBlockVector& matrices, + bool both) + { + M1.resize(matrices.size()); + if (both) + M2.resize(matrices.size()); + for (unsigned int i=0;i >& o = matrices.block(i); + const unsigned int row = o[o.min_level()].row; + const unsigned int col = o[o.min_level()].column; + + M1[i].row = row; + M1[i].column = col; + if (both) + { + M2[i].row = row; + M2[i].column = col; + } + } + } + + template inline void LocalResults::initialize_matrices(const unsigned int n, 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 c2c62ee407..922480da07 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker_assembler.h +++ b/deal.II/deal.II/include/numerics/mesh_worker_assembler.h @@ -904,9 +904,10 @@ namespace MeshWorker class MGMatrixLocalBlocksToGlobalBlocks { public: - /// The object that is stored - typedef std_cxx1x::shared_ptr > > MatrixPtr; - + typedef MGMatrixBlockVector MatrixPtrVector; + typedef SmartPointer > + MatrixPtrVectorPtr; + /** * Constructor, initializing * the #threshold, which @@ -923,7 +924,7 @@ namespace MeshWorker * cell matrix vectors. */ void initialize(const BlockInfo* block_info, - std::vector& matrices); + MatrixPtrVector& matrices); /** * Multigrid methods on @@ -936,7 +937,7 @@ namespace MeshWorker * refinement edge, which are * set by this method. */ - void initialize_edge_flux(std::vector& up, std::vector& down); + void initialize_edge_flux(MatrixPtrVector& up, MatrixPtrVector& down); /** * Initialize the local data * in the @@ -995,35 +996,35 @@ namespace MeshWorker * stored as a vector of * pointers. */ - std::vector matrices; + MatrixPtrVectorPtr matrices; /** * The flux matrix between * the fine and the coarse * level at refinement edges. */ - std::vector flux_down; + MatrixPtrVectorPtr flux_down; /** * The flux matrix between * the coarse and the fine * level at refinement edges. */ - std::vector flux_up; + MatrixPtrVectorPtr flux_up; /** * The interface matrix between * the fine and the coarse * level at refinement edges. */ - std::vector interface_out; + MatrixPtrVectorPtr interface_out; /** * The interface matrix between * the coarse and the fine * level at refinement edges. */ - std::vector interface_in; + MatrixPtrVectorPtr interface_in; /** * A pointer to the object containing the block structure. @@ -1939,11 +1940,12 @@ namespace MeshWorker template inline void - MGMatrixLocalBlocksToGlobalBlocks::initialize(const BlockInfo* b, - std::vector& m) + MGMatrixLocalBlocksToGlobalBlocks::initialize( + const BlockInfo* b, + MatrixPtrVector& m) { block_info = b; - matrices = m; + matrices = &m; } @@ -1954,7 +1956,7 @@ namespace MeshWorker DOFINFO& info, bool interior_face) const { - info.initialize_matrices(matrices, interior_face); + info.initialize_matrices(*matrices, interior_face); } @@ -1962,8 +1964,8 @@ namespace MeshWorker template inline void MGMatrixLocalBlocksToGlobalBlocks::initialize_edge_flux( - std::vector& up, - std::vector& down) + MatrixPtrVector& up, + MatrixPtrVector& down) { flux_up = up; flux_down = down; @@ -2020,19 +2022,18 @@ namespace MeshWorker template template inline void - MGMatrixLocalBlocksToGlobalBlocks::assemble( - const DOFINFO& info) + MGMatrixLocalBlocksToGlobalBlocks::assemble(const DOFINFO& info) { const unsigned int level = info.cell->level(); - for (unsigned int i=0;isize();++i) { // Row and column index of // the block we are dealing with - const unsigned int row = matrices[i]->row; - const unsigned int col = matrices[i]->column; + const unsigned int row = matrices->block(i)[level].row; + const unsigned int col = matrices->block(i)[level].column; - assemble(matrices[i]->matrix[level], info.matrix(i,false).matrix, row, col, + assemble(matrices->block(i)[level].matrix, info.matrix(i,false).matrix, row, col, info.indices, info.indices, level, level); } } @@ -2048,33 +2049,39 @@ namespace MeshWorker const unsigned int level1 = info1.cell->level(); const unsigned int level2 = info2.cell->level(); - for (unsigned int i=0;isize();++i) { + MGLevelObject >& o = matrices->block(i); + // Row and column index of // the block we are dealing with - const unsigned int row = matrices[i]->row; - const unsigned int col = matrices[i]->column; + const unsigned int row = o[level1].row; + const unsigned int col = o[level1].column; if (level1 == level2) { - assemble(matrices[i]->matrix[level1], info1.matrix(i,false).matrix, row, col, info1.indices, info1.indices, level1, level1); - assemble(matrices[i]->matrix[level1], info1.matrix(i,true).matrix, row, col, info1.indices, info2.indices, level1, level2); - assemble(matrices[i]->matrix[level1], info2.matrix(i,false).matrix, row, col, info2.indices, info2.indices, level2, level2); - assemble(matrices[i]->matrix[level1], info2.matrix(i,true).matrix, row, col, info2.indices, info1.indices, level2, level1); + assemble(o[level1].matrix, info1.matrix(i,false).matrix, row, col, + info1.indices, info1.indices, level1, level1); + assemble(o[level1].matrix, info1.matrix(i,true).matrix, row, col, + info1.indices, info2.indices, level1, level2); + assemble(o[level1].matrix, info2.matrix(i,false).matrix, row, col, + info2.indices, info2.indices, level2, level2); + assemble(o[level1].matrix, info2.matrix(i,true).matrix, row, col, + info2.indices, info1.indices, level2, level1); } else { Assert(level1 > level2, ExcNotImplemented()); - if (flux_up.size() != 0) + if (flux_up->size() != 0) { // Do not add M22, // which is done by // the coarser cell - assemble(matrices[i]->matrix[level1], info1.matrix(i,false).matrix, row, col, + assemble(o[level1].matrix, info1.matrix(i,false).matrix, row, col, info1.indices, info1.indices, level1, level1); - assemble(flux_up[i]->matrix[level1], info1.matrix(i,true).matrix, row, col, + assemble(flux_up->block(i)[level1].matrix, info1.matrix(i,true).matrix, row, col, info1.indices, info2.indices, level1, level2, true); - assemble(flux_down[i]->matrix[level1], info2.matrix(i,true).matrix, row, col, + assemble(flux_down->block(i)[level1].matrix, info2.matrix(i,true).matrix, row, col, info2.indices, info1.indices, level2, level1); } } -- 2.39.5