From: kanschat Date: Mon, 12 Mar 2012 13:43:15 +0000 (+0000) Subject: implement new block structure in DoFInfo X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=17ddd8c4dc62f846889fa43bcb4b42af1ed2987b;p=dealii-svn.git implement new block structure in DoFInfo git-svn-id: https://svn.dealii.org/trunk@25263 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/include/deal.II/meshworker/dof_info.h b/deal.II/include/deal.II/meshworker/dof_info.h index 93055eddda..40addebb6c 100644 --- a/deal.II/include/deal.II/meshworker/dof_info.h +++ b/deal.II/include/deal.II/meshworker/dof_info.h @@ -1,7 +1,7 @@ //--------------------------------------------------------------------------- // $Id$ // -// Copyright (C) 2006, 2007, 2008, 2009, 2010, 2011 by the deal.II authors +// Copyright (C) 2006, 2007, 2008, 2009, 2010, 2011, 2012 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -93,9 +93,19 @@ namespace MeshWorker */ unsigned int sub_number; - /// The DoF indices of the current cell + /* + * The DoF indices of the + * current cell + */ std::vector indices; + /** + * The DoF indices on the + * current cell, organized by + * local blocks + */ + std::vector > indices_by_block; + /** * Constructor setting the * #block_info pointer. diff --git a/deal.II/include/deal.II/meshworker/dof_info.templates.h b/deal.II/include/deal.II/meshworker/dof_info.templates.h index faa56529ce..d43af1b649 100644 --- a/deal.II/include/deal.II/meshworker/dof_info.templates.h +++ b/deal.II/include/deal.II/meshworker/dof_info.templates.h @@ -1,7 +1,7 @@ //--------------------------------------------------------------------------- // $Id$ // -// Copyright (C) 2009, 2010, 2011 by the deal.II authors +// Copyright (C) 2009, 2010, 2011, 2012 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -22,7 +22,11 @@ namespace MeshWorker template DoFInfo::DoFInfo(const BlockInfo& info) : block_info(&info, typeid(*this).name()) - {} + { + indices_by_block.resize(info.local().size()); + for (unsigned int i=0;i @@ -44,6 +48,14 @@ namespace MeshWorker { indices_org.resize(c->get_fe().dofs_per_cell); c->get_dof_indices(indices_org); + for (unsigned int i=0;i + bi = block_info->local().global_to_local(i); + indices_by_block[bi.first][bi.second] = indices_org[i]; + } + // Remove this after + // changing block codes for (unsigned int i=0;iblock_info->renumber(i)] = indices_org[i]; } diff --git a/deal.II/include/deal.II/meshworker/simple.h b/deal.II/include/deal.II/meshworker/simple.h index 56f90d1010..6e40780534 100644 --- a/deal.II/include/deal.II/meshworker/simple.h +++ b/deal.II/include/deal.II/meshworker/simple.h @@ -115,6 +115,16 @@ namespace MeshWorker * this class can be used in a MeshWorker::loop() to assemble the cell * and face matrices into the global matrix. * + * The assembler can handle two different types of local data. First, + * by default, the obvious choice of taking a single local matrix with + * dimensions equal to the number of degrees of freedom of the + * cell. Alternatively, a local block structure can be initialized + * with initialize_local_blocks(). After this, the local data will be + * arranged as an array of n by n FullMatrix blocks, which are + * ordered lexicographically in DoFInfo. Note that + * initialize_local_blocks() has to be called before initialize_info() + * to take the desired effect. + * * @todo On locally refined meshes, a ConstraintMatrix should be used * to automatically eliminate hanging nodes. * @@ -140,9 +150,45 @@ namespace MeshWorker */ void initialize(MATRIX& m); /** - * Initialize the constraints. + * Initialize the + * constraints for laer use + * by the assemble functions. */ void initialize(const ConstraintMatrix& constraints); + + /** + * Store information on the + * local block structure. If + * the assembler is + * inititialized with this + * function, + * initialize_info() will + * generate one local matrix + * for each block row and + * column, whch will be + * numbered + * lexicographically, row by + * row. + * + * In spite of using local + * block structure, all + * blocks will be enteres + * into the same global + * matrix, disregarding any + * global block structure. + * + * @note The argument of this + * function will be copied + * into the member object + * #local_indices. Thus, + * every subsequent change in + * the block structure must + * be initialzied or will not + * be used by the assembler. + */ + + void initialize_local_blocks(const BlockIndices& local_indices); + /** * Initialize the local data * in the @@ -191,10 +237,21 @@ namespace MeshWorker */ SmartPointer > matrix; /** - * A pointer to the object containing constraints. + * A pointer to the object + * containing constraints. */ SmartPointer > constraints; + /** + * The object containing the + * local block structure. Set + * by + * initialize_local_blocks() + * and used by assembling + * functions. + */ + BlockIndices local_indices; + /** * The smallest positive * number that will be @@ -579,7 +636,15 @@ namespace MeshWorker constraints = &c; } + + template + inline void + MatrixSimple::initialize_local_blocks(const BlockIndices& b) + { + local_indices = b; + } + template template inline void @@ -598,7 +663,7 @@ namespace MeshWorker { AssertDimension(M.m(), i1.size()); AssertDimension(M.n(), i2.size()); - + if(constraints == 0) { for (unsigned int j=0; j::assemble(const DOFINFO& info) { - assemble(info.matrix(0,false).matrix, info.indices, info.indices); + if (local_indices.size() == 0) + assemble(info.matrix(0,false).matrix, info.indices, info.indices); + else + { + } } @@ -627,12 +696,18 @@ namespace MeshWorker MatrixSimple::assemble(const DOFINFO& info1, const DOFINFO& info2) { - assemble(info1.matrix(0,false).matrix, info1.indices, info1.indices); - assemble(info1.matrix(0,true).matrix, info1.indices, info2.indices); - assemble(info2.matrix(0,false).matrix, info2.indices, info2.indices); - assemble(info2.matrix(0,true).matrix, info2.indices, info1.indices); + if (local_indices.size() == 0) + { + assemble(info1.matrix(0,false).matrix, info1.indices, info1.indices); + assemble(info1.matrix(0,true).matrix, info1.indices, info2.indices); + assemble(info2.matrix(0,false).matrix, info2.indices, info2.indices); + assemble(info2.matrix(0,true).matrix, info2.indices, info1.indices); + } + else + { + } } - + //----------------------------------------------------------------------//