1#ifndef DUNE_FEM_SOLVER_COMMUNICATION_FEM_HH
2#define DUNE_FEM_SOLVER_COMMUNICATION_FEM_HH
9#include <dune/common/hybridutilities.hh>
11#include <dune/istl/solvercategory.hh>
13#include <dune/fem/space/common/commoperations.hh>
27 template<
class DiscreteFunctionSpace,
class DofVector >
28 class FemCommunicationVector
30 typedef FemCommunicationVector< DiscreteFunctionSpace, DofVector > ThisType;
35 typedef typename DiscreteFunctionSpaceType::LocalBlockIndices BlockIndices;
37 static constexpr std::size_t blockSize = Hybrid::size( BlockIndices() );
39 typedef typename DofVector::field_type DofType;
41 template<
class Operation >
44 typedef typename DiscreteFunctionSpaceType::template CommDataHandle< ThisType, Operation >::Type Type;
48 : dfSpace_( dfSpace ), dofVector_( dofVector )
51 template<
class Operation >
52 typename CommDataHandle< Operation >::Type dataHandle (
const Operation &operation )
54 return space().createDataHandle( *
this, operation );
57 const DofVector &dofVector ()
const {
return dofVector_; }
58 DofVector &dofVector () {
return dofVector_; }
60 const DiscreteFunctionSpaceType &space ()
const {
return dfSpace_; }
64 DofVector &dofVector_;
72 template<
class DiscreteFunctionSpace >
73 class FemCommunication
75 typedef FemCommunication< DiscreteFunctionSpace > ThisType;
81 : dfSpace_( dfSpace ), solverCategory_( solverCategory )
87 void copyOwnerToAll (
const T &x, T &y )
const
91 FemCommunicationVector< DiscreteFunctionSpaceType, T > z( dfSpace_, y );
92 dfSpace_.communicator().exchange( z, DFCommunicationOperation::Add() );
96 void project ( T &x )
const
98 typedef typename T::field_type field_type;
101 for(
int i : dfSpace_.auxiliaryDofs() )
102 x[ i ] = field_type( 0 );
105 template<
class T,
class F >
106 void dot (
const T &x,
const T &y, F &scp )
const
108 const auto &auxiliaryDofs = dfSpace_.auxiliaryDofs();
110 const int numAuxiliarys = auxiliaryDofs.size();
111 for(
int auxiliary = 0, i = 0; auxiliary < numAuxiliarys; ++auxiliary, ++i )
113 const int nextAuxiliary = auxiliaryDofs[ auxiliary ];
114 for( ; i < nextAuxiliary; ++i )
115 scp += x[ i ] * y[ i ];
118 scp = communicator().sum( scp );
122 typename Dune::FieldTraits< typename T::field_type >::real_type norm (
const T &x )
const
125 typename Dune::FieldTraits< typename T::field_type >::real_type norm2( 0 );
127 return sqrt( norm2 );
133 const DiscreteFunctionSpaceType &dfSpace_;
142 template<
class DiscreteFunctionSpace >
145 std::shared_ptr< FemCommunication< DiscreteFunctionSpace > > &communication )
147 communication.reset(
new FemCommunication< DiscreteFunctionSpace >( dfSpace, solverCategory ) );
158 template<
class DiscreteFunctionSpace >
160 :
public std::false_type
Traits::CommunicationType CommunicationType
Collective communication.
Definition: gridpart.hh:97
Type traits to determine the type of reals (when working with complex numbers)
auto dot(const A &a, const B &b) -> typename std::enable_if< IsNumber< A >::value &&!IsVector< A >::value &&!std::is_same< typename FieldTraits< A >::field_type, typename FieldTraits< A >::real_type > ::value, decltype(conj(a) *b)>::type
computes the dot product for fundamental data types according to Petsc's VectDot function: dot(a,...
Definition: dotproduct.hh:42
Dune namespace.
Definition: alignedallocator.hh:13
Category
Definition: solvercategory.hh:23
@ sequential
Category for sequential solvers.
Definition: solvercategory.hh:25