2 #ifndef DUNE_FEM_PETSCDOFMAPPING_HH 3 #define DUNE_FEM_PETSCDOFMAPPING_HH 39 template<
class SlaveDofs >
40 class PetscDofMappings
43 typedef PetscDofMappings< SlaveDofs > ThisType;
47 typedef SlaveDofs SlaveDofsType;
49 typedef std::vector< int > DofMappingType;
50 typedef DofMappingType::size_type IndexType;
51 typedef typename DofMappingType::value_type DofType;
53 typedef PetscInt GlobalDofType ;
54 typedef DynamicVector< GlobalDofType > GlobalDofMappingType ;
56 PetscDofMappings ( SlaveDofs *slaveDofs )
58 : slaveDofs_( *slaveDofs ),
59 numOwnedDofBlocks_( 0 ),
61 processStartIndex_( 0 ),
73 const int sequence = slaveDofs_.space().sequence();
74 if( sequence_ != sequence )
76 initialize( slaveDofs_ );
77 sequence_ = sequence ;
84 GlobalDofType numOwnedDofBlocks ()
const {
return numOwnedDofBlocks_; }
86 GlobalDofType numSlaveBlocks ()
const {
return numSlaveBlocks_; }
88 GlobalDofType processStartIndex ()
const {
return processStartIndex_; }
90 size_t size ()
const {
return localSlaveMapping_.size(); }
92 const DofType& localSlaveMapping ( IndexType index )
const 94 return localSlaveMapping_[ index ];
97 const GlobalDofType& globalMapping (
const IndexType index )
const 99 return globalDofMapping_[ index ];
103 bool isSlave ( IndexType i )
const 105 return static_cast< int >( localSlaveMapping( i ) ) >= numOwnedDofBlocks();
112 template <
class Stream >
113 void write(
const Stream& )
const {}
115 template <
class Stream >
116 void read(
const Stream& ) {}
118 void resize () { update (); }
119 bool compress() {
return update (); }
126 PetscDofMappings (
const ThisType& );
127 PetscDofMappings& operator= (
const ThisType& );
132 template<
typename SlaveDofProv
ider >
133 void initializeMappings ( SlaveDofProvider& slaveDofs )
142 typedef typename SlaveDofProvider :: DiscreteFunctionSpaceType SpaceType;
143 const SpaceType& space = slaveDofs.space();
146 int ownedDofBlocks = 0;
149 localSlaveMapping_.resize( space.blockMapper().size(), -1 );
150 globalDofMapping_.resize( space.blockMapper().size(), 0 );
153 GlobalDofType index = processStartIndex_ ;
154 for(
int slave = 0, i = 0; slave < slaveDofs.size(); ++slave )
156 const int nextSlave = slaveDofs[ slave ];
157 for(; i < nextSlave; ++i )
159 localSlaveMapping_[ i ] = i - slave;
160 globalDofMapping_[ i ] = index++;
168 if( static_cast< size_t >( i ) == localSlaveMapping_.size() )
171 assert( static_cast< size_t >( i ) < localSlaveMapping_.size() );
173 localSlaveMapping_[ i ] = numOwnedDofBlocks_ + slave;
174 globalDofMapping_[ i ] = 0;
179 checkDofMappingConsistency();
181 assert( numOwnedDofBlocks_ == ownedDofBlocks );
183 typedef typename SpaceType :: template ToNewDimRange< 1 > :: Type DofSpaceType ;
185 if( space.continuous() )
187 typedef DofSpaceType GlobalDofSpaceType ;
188 GlobalDofSpaceType dofSpace( space.gridPart() );
194 VectorDiscreteFunction< GlobalDofSpaceType, GlobalDofMappingType >
195 dofMappingFunc(
"globalDofs", dofSpace, globalDofMapping_ );
198 dofMappingFunc.communicate();
203 typedef FiniteVolumeSpace<
typename DofSpaceType :: FunctionSpaceType,
204 typename DofSpaceType :: GridPartType > GlobalDofSpaceType ;
205 GlobalDofSpaceType dofSpace( space.gridPart() );
211 VectorDiscreteFunction< GlobalDofSpaceType, GlobalDofMappingType >
212 dofMappingFunc(
"globalDofs", dofSpace, globalDofMapping_ );
215 dofMappingFunc.communicate();
219 void initialize (
const SlaveDofsType& slaveDofs )
221 numSlaveBlocks_ = slaveDofs.size() - 1;
222 numOwnedDofBlocks_ = slaveDofs.space().blockMapper().size() - numSlaveBlocks_;
225 unsigned long processStartIndex = 0;
226 unsigned long numOwnedDofBlocks = numOwnedDofBlocks_;
229 MPI_Exscan( &numOwnedDofBlocks, &processStartIndex, 1, MPI_UNSIGNED_LONG,
MPI_SUM, PETSC_COMM_WORLD );
232 processStartIndex_ = processStartIndex ;
234 initializeMappings( slaveDofs );
237 checkSlaveConsistency( slaveDofs );
244 template<
typename SlavesType >
245 void checkSlaveConsistency (
const SlavesType& slaves )
const 247 for(
size_t i = 0; i < size(); ++i )
249 assert( isSlave( i ) == slaves.isSlave( i ) );
253 void checkDofMappingConsistency ()
const 257 std::map< DofType, bool > buf;
258 for( std::vector< int >::const_iterator it = localSlaveMapping_.begin(); it != localSlaveMapping_.end(); ++it )
262 std::cerr <<
"localSlaveMapping_ has not been initialized completely on rank " <<
MPIManager::rank() << std::endl;
266 if( buf.find( *it ) != buf.end() )
268 std::cerr << *it <<
" was found more than once in localSlaveMapping_ (on rank " <<
MPIManager::rank() <<
")\n";
281 SlaveDofsType& slaveDofs_;
282 GlobalDofType numOwnedDofBlocks_;
283 GlobalDofType numSlaveBlocks_;
284 GlobalDofType processStartIndex_;
285 DofMappingType localSlaveMapping_;
286 GlobalDofMappingType globalDofMapping_;
294 #endif // #if HAVE_PETSC 296 #endif // DUNE_FEM_PETSCDOFMAPPING_HH
static int rank()
Definition: mpimanager.hh:116
static const int MPI_SUM
Definition: odeobjectfiles.cc:18
Definition: coordinate.hh:4