Commit 322859ec authored by Nils Kohl's avatar Nils Kohl 🌝
Browse files

Implemented the P1DGE constant vector Laplace operator.

Includes some fixes in the P0toP1 and P1toP0 mixed operators.
Also switched to P0Function for DGE.
parent 1ff584e1
Pipeline #39117 failed with stages
in 30 minutes and 27 seconds
......@@ -6,6 +6,7 @@ target_sources( hyteg
DGDiffusionForm_Example.cpp
DGForm.hpp
DGForm2D.hpp
DGFormAbort.hpp
DGFormVolume.hpp
DGFunction.cpp
DGFunction.hpp
......@@ -18,6 +19,7 @@ target_sources( hyteg
P0_to_P1_divt_form.cpp
P0_to_P1_divt_form.hpp
DGVectorLaplaceForm.hpp
DGVectorMassForm.hpp
DGVectorFunction.hpp
DGVectorOperators.hpp
)
......
......@@ -54,10 +54,10 @@ class P0ScalarToP1VectorOperator : public Operator< P0Function< real_t >, P1Vect
}
void toMatrix( const std::shared_ptr< SparseMatrixProxy >& mat,
const P0Function< idx_t >& src,
const P1VectorFunction< idx_t >& dst,
size_t level,
DoFType flag ) const
const P0Function< idx_t >& src,
const P1VectorFunction< idx_t >& dst,
size_t level,
DoFType flag ) const
{
operX.toMatrix( mat, src, dst[0], level, flag );
operY.toMatrix( mat, src, dst[1], level, flag );
......@@ -75,4 +75,9 @@ class P0ScalarToP1VectorOperator : public Operator< P0Function< real_t >, P1Vect
typedef P0ScalarToP1VectorOperator< P0ToP1ConstantDivTxOperator, P0ToP1ConstantDivTyOperator, P0ToP1ConstantDivTzOperator >
P0ToP1ConstantDivTOperator;
typedef P0ScalarToP1VectorOperator< P0ToP1ConstantP1EDGVectorLaplaceXCouplingOperator,
P0ToP1ConstantP1EDGVectorLaplaceYCouplingOperator,
P0ToP1ConstantP1EDGVectorLaplaceZCouplingOperator >
P0ToP1ConstantP1EDGVectorLaplaceCouplingOperator;
} // namespace hyteg
......@@ -23,6 +23,7 @@
#include <hyteg/communication/Syncing.hpp>
#include "hyteg/celldofspace/CellDoFIndexing.hpp"
#include "hyteg/dgfunctionspace/DGFormAbort.hpp"
#include "hyteg/dgfunctionspace/DGFunction.hpp"
#include "hyteg/dgfunctionspace/P0_to_P1_divt_form.hpp"
#include "hyteg/functions/Function.hpp"
......@@ -125,15 +126,15 @@ class P0ToP1Operator : public Operator< P0Function< real_t >, P1Function< real_t
std::vector< PrimitiveID > pids;
if ( dim == 2 )
{
src.template communicate< Face, Face >( level );
pids = storage->getFaceIDs();
}
else
{
src.template communicate< Cell, Cell >( level );
pids = storage->getCellIDs();
}
src.communicate( level );
for ( const auto& pid : pids )
{
const auto srcPolyDegree = 0;
......@@ -560,14 +561,14 @@ class P0ToP1Operator : public Operator< P0Function< real_t >, P1Function< real_t
if ( dim == 2 )
{
auto nVertexDoFIndicesArray = facedof::macroface::getMicroVerticesFromMicroFace(
neighborInfo.neighborElementIndices( n ), faceType );
neighborInfo.neighborElementIndices( n ), neighborInfo.neighborFaceType( n ) );
nVertexDoFIndices.insert(
nVertexDoFIndices.begin(), nVertexDoFIndicesArray.begin(), nVertexDoFIndicesArray.end() );
}
else
{
auto nVertexDoFIndicesArray = celldof::macrocell::getMicroVerticesFromMicroCell(
neighborInfo.neighborElementIndices( n ), cellType );
neighborInfo.neighborElementIndices( n ), neighborInfo.neighborCellType( n ) );
nVertexDoFIndices.insert(
nVertexDoFIndices.begin(), nVertexDoFIndicesArray.begin(), nVertexDoFIndicesArray.end() );
}
......@@ -629,7 +630,7 @@ class P0ToP1Operator : public Operator< P0Function< real_t >, P1Function< real_t
const auto globalColIdx = srcDofMemory[volumedofspace::indexing::index(
neighborInfo.neighborElementIndices( n ).x(),
neighborInfo.neighborElementIndices( n ).y(),
faceType,
neighborInfo.neighborFaceType( n ),
0,
1,
level,
......@@ -673,120 +674,6 @@ class P0ToP1Operator : public Operator< P0Function< real_t >, P1Function< real_t
WALBERLA_UNUSED( flag );
}
/// \brief This is similar to the implementation in the dg::DGOperator class.
template < typename VType >
inline void assembleAndOrApplyOld( const P0Function< VType >& src,
const P1Function< VType >& dst,
size_t level,
DoFType flag,
const std::shared_ptr< SparseMatrixProxy >& mat,
UpdateType updateType = Replace ) const
{
using indexing::Index;
WALBERLA_CHECK( updateType == Replace );
auto srcDGF = src.getDGFunction();
communication::syncFunctionBetweenPrimitives( dst, level );
if ( this->getStorage()->hasGlobalCells() )
{
WALBERLA_ABORT( "P0 to P1 apply not implemented in 3D." );
}
else
{
const int dim = 2;
for ( const auto& faceIt : this->getStorage()->getFaces() )
{
const auto faceId = faceIt.first;
const auto& face = *faceIt.second;
const auto srcPolyDegree = 0;
const auto dstPolyDegree = 1;
const auto numSrcDofs = 1;
const auto numDstDofs = 3;
const auto srcDofMemory = srcDGF->volumeDoFFunction()->dofMemory( faceId, level );
auto dstDofMemory = face.getData( dst.getFaceDataID() )->getPointer( level );
const auto srcMemLayout = srcDGF->volumeDoFFunction()->memoryLayout();
for ( auto faceType : facedof::allFaceTypes )
{
for ( auto elementIdx : facedof::macroface::Iterator( level, faceType ) )
{
// TODO: all these coord computations can be executed _once_ and then the coordinates can be incremented by h
// TODO: blending
// This object does the heavy lifting of computing all required coordinates and normals.
volumedofspace::indexing::ElementNeighborInfo neighborInfo(
elementIdx, faceType, level, src.getBoundaryCondition(), faceId, this->getStorage() );
// We only write to the DoFs in the current volume, let's prepare a temporary vector for that.
Eigen::Matrix< real_t, Eigen::Dynamic, 1 > dstDofs;
dstDofs.resize( numDstDofs, Eigen::NoChange_t::NoChange );
dstDofs.setZero();
/////////////////////////
// Volume contribution //
/////////////////////////
std::array< Point3D, 3 > coords;
for ( uint_t i = 0; i < neighborInfo.elementVertexCoords().size(); i++ )
{
coords[i][0] = neighborInfo.elementVertexCoords().at( i )( 0 );
coords[i][1] = neighborInfo.elementVertexCoords().at( i )( 1 );
coords[i][2] = neighborInfo.elementVertexCoords().at( i )( 2 );
}
Matrixr< 3, 1 > elMat;
form_->integrateAll( coords, elMat );
// P0 has only one DoF
const auto srcDoF = srcDofMemory[volumedofspace::indexing::index(
elementIdx.x(), elementIdx.y(), faceType, 0, 1, level, srcMemLayout )];
// Micro-vertex indices
auto vertexDoFIndices = facedof::macroface::getMicroVerticesFromMicroFace( elementIdx, faceType );
if ( mat == nullptr )
{
// Matrix-vector multiplication.
// dstDofs += localMat * srcDofs;
WALBERLA_ABORT( "Not implemented." );
}
else
{
// Sparse assembly.
for ( uint_t dstDofIdx = 0; dstDofIdx < numDstDofs; dstDofIdx++ )
{
const auto globalRowIdx = dstDofMemory[vertexdof::macroface::index(
level, vertexDoFIndices.at( dstDofIdx ).x(), vertexDoFIndices.at( dstDofIdx ).y() )];
const auto globalColIdx = srcDofMemory[volumedofspace::indexing::index(
elementIdx.x(), elementIdx.y(), faceType, 0, 1, level, srcMemLayout )];
mat->addValue( globalRowIdx, globalColIdx, elMat( dstDofIdx, 0 ) );
}
}
if ( mat == nullptr )
{
// Write DoFs.
for ( uint_t dstDofIdx = 0; dstDofIdx < numDstDofs; dstDofIdx++ )
{
WALBERLA_ABORT( "Not implemented." );
}
}
}
}
}
}
WALBERLA_UNUSED( flag );
}
std::shared_ptr< Form > form_;
};
......@@ -794,4 +681,8 @@ typedef P0ToP1Operator< dg::p0_to_p1_divt_0_affine_q0 > P0ToP1ConstantDivTxOpera
typedef P0ToP1Operator< dg::p0_to_p1_divt_1_affine_q0 > P0ToP1ConstantDivTyOperator;
typedef P0ToP1Operator< dg::p0_to_p1_divt_2_affine_q0 > P0ToP1ConstantDivTzOperator;
typedef P0ToP1Operator< dg::DGVectorLaplaceFormP1EDG_0 > P0ToP1ConstantP1EDGVectorLaplaceXCouplingOperator;
typedef P0ToP1Operator< dg::DGVectorLaplaceFormP1EDG_1 > P0ToP1ConstantP1EDGVectorLaplaceYCouplingOperator;
typedef P0ToP1Operator< dg::DGFormAbort > P0ToP1ConstantP1EDGVectorLaplaceZCouplingOperator;
} // namespace hyteg
\ No newline at end of file
......@@ -23,6 +23,7 @@
#include <hyteg/communication/Syncing.hpp>
#include "hyteg/celldofspace/CellDoFIndexing.hpp"
#include "hyteg/dgfunctionspace/DGFormAbort.hpp"
#include "hyteg/dgfunctionspace/DGFunction.hpp"
#include "hyteg/dgfunctionspace/P1_to_P0_div_form.hpp"
#include "hyteg/functions/Function.hpp"
......@@ -550,14 +551,14 @@ class P1ToP0Operator : public Operator< P1Function< real_t >, P0Function< real_t
if ( dim == 2 )
{
auto nVertexDoFIndicesArray = facedof::macroface::getMicroVerticesFromMicroFace(
neighborInfo.neighborElementIndices( n ), faceType );
neighborInfo.neighborElementIndices( n ), neighborInfo.neighborFaceType( n ) );
nVertexDoFIndices.insert(
nVertexDoFIndices.begin(), nVertexDoFIndicesArray.begin(), nVertexDoFIndicesArray.end() );
}
else
{
auto nVertexDoFIndicesArray = celldof::macrocell::getMicroVerticesFromMicroCell(
neighborInfo.neighborElementIndices( n ), cellType );
neighborInfo.neighborElementIndices( n ), neighborInfo.neighborCellType( n ) );
nVertexDoFIndices.insert(
nVertexDoFIndices.begin(), nVertexDoFIndicesArray.begin(), nVertexDoFIndicesArray.end() );
}
......@@ -608,7 +609,6 @@ class P1ToP0Operator : public Operator< P1Function< real_t >, P0Function< real_t
{
neighborMicroVolType = invCellTypeMap[neighborInfo.neighborCellType( n )];
}
// Sparse assembly.
for ( uint_t srcDofIdx = 0; srcDofIdx < numSrcDofs; srcDofIdx++ )
{
......@@ -664,120 +664,6 @@ class P1ToP0Operator : public Operator< P1Function< real_t >, P0Function< real_t
WALBERLA_UNUSED( flag );
}
/// \brief This is similar to the implementation in the dg::DGOperator class.
template < typename VType >
inline void assembleAndOrApplyOld( const P1Function< VType >& src,
const P0Function< VType >& dst,
size_t level,
DoFType flag,
const std::shared_ptr< SparseMatrixProxy >& mat,
UpdateType updateType = Replace ) const
{
using indexing::Index;
WALBERLA_CHECK( updateType == Replace );
auto dstDGF = dst.getDGFunction();
communication::syncFunctionBetweenPrimitives( src, level );
if ( this->getStorage()->hasGlobalCells() )
{
WALBERLA_ABORT( "P1 to P0 apply not implemented in 3D." );
}
else
{
const int dim = 2;
for ( const auto& faceIt : this->getStorage()->getFaces() )
{
const auto faceId = faceIt.first;
const auto& face = *faceIt.second;
const auto srcPolyDegree = 1;
const auto dstPolyDegree = 0;
const auto numSrcDofs = 3;
const auto numDstDofs = 1;
auto dstDofMemory = dstDGF->volumeDoFFunction()->dofMemory( faceId, level );
auto srcDofMemory = face.getData( src.getFaceDataID() )->getPointer( level );
const auto dstMemLayout = dstDGF->volumeDoFFunction()->memoryLayout();
for ( auto faceType : facedof::allFaceTypes )
{
for ( auto elementIdx : facedof::macroface::Iterator( level, faceType ) )
{
// TODO: all these coord computations can be executed _once_ and then the coordinates can be incremented by h
// TODO: blending
// This object does the heavy lifting of computing all required coordinates and normals.
volumedofspace::indexing::ElementNeighborInfo neighborInfo(
elementIdx, faceType, level, dst.getBoundaryCondition(), faceId, this->getStorage() );
// We only write to the DoFs in the current volume, let's prepare a temporary vector for that.
Eigen::Matrix< real_t, Eigen::Dynamic, 1 > dstDofs;
dstDofs.resize( numDstDofs, Eigen::NoChange_t::NoChange );
dstDofs.setZero();
/////////////////////////
// Volume contribution //
/////////////////////////
std::array< Point3D, 3 > coords;
for ( uint_t i = 0; i < neighborInfo.elementVertexCoords().size(); i++ )
{
coords[i][0] = neighborInfo.elementVertexCoords().at( i )( 0 );
coords[i][1] = neighborInfo.elementVertexCoords().at( i )( 1 );
coords[i][2] = neighborInfo.elementVertexCoords().at( i )( 2 );
}
Matrixr< 1, 3 > elMat;
form_->integrateAll( coords, elMat );
// P0 has only one DoF
// const auto dstDoF = dstDofMemory[volumedofspace::indexing::index(
// elementIdx.x(), elementIdx.y(), faceType, 0, 1, level, dstMemLayout )];
// Micro-vertex indices
auto vertexDoFIndices = facedof::macroface::getMicroVerticesFromMicroFace( elementIdx, faceType );
if ( mat == nullptr )
{
// Matrix-vector multiplication.
// dstDofs += localMat * srcDofs;
WALBERLA_ABORT( "Not implemented." );
}
else
{
// Sparse assembly.
for ( uint_t srcDofIdx = 0; srcDofIdx < numSrcDofs; srcDofIdx++ )
{
const auto globalRowIdx = dstDofMemory[volumedofspace::indexing::index(
elementIdx.x(), elementIdx.y(), faceType, 0, 1, level, dstMemLayout )];
const auto globalColIdx = srcDofMemory[vertexdof::macroface::index(
level, vertexDoFIndices.at( srcDofIdx ).x(), vertexDoFIndices.at( srcDofIdx ).y() )];
mat->addValue( globalRowIdx, globalColIdx, elMat( 0, srcDofIdx ) );
}
}
if ( mat == nullptr )
{
// Write DoFs.
for ( uint_t dstDofIdx = 0; dstDofIdx < numDstDofs; dstDofIdx++ )
{
WALBERLA_ABORT( "Not implemented." );
}
}
}
}
}
}
WALBERLA_UNUSED( flag );
}
std::shared_ptr< Form > form_;
};
......@@ -785,4 +671,8 @@ typedef P1ToP0Operator< dg::p1_to_p0_div_0_affine_q0 > P1ToP0ConstantDivxOperato
typedef P1ToP0Operator< dg::p1_to_p0_div_1_affine_q0 > P1ToP0ConstantDivyOperator;
typedef P1ToP0Operator< dg::p1_to_p0_div_2_affine_q0 > P1ToP0ConstantDivzOperator;
typedef P1ToP0Operator< dg::DGVectorLaplaceFormEDGP1_0 > P1ToP0ConstantP1EDGVectorLaplaceXCouplingOperator;
typedef P1ToP0Operator< dg::DGVectorLaplaceFormEDGP1_1 > P1ToP0ConstantP1EDGVectorLaplaceYCouplingOperator;
typedef P1ToP0Operator< dg::DGFormAbort > P1ToP0ConstantP1EDGVectorLaplaceZCouplingOperator;
} // namespace hyteg
\ No newline at end of file
......@@ -76,4 +76,9 @@ class P1VectorToP0ScalarOperator : public Operator< P1VectorFunction< real_t >,
typedef P1VectorToP0ScalarOperator< P1ToP0ConstantDivxOperator, P1ToP0ConstantDivyOperator, P1ToP0ConstantDivzOperator >
P1ToP0ConstantDivOperator;
typedef P1VectorToP0ScalarOperator< P1ToP0ConstantP1EDGVectorLaplaceXCouplingOperator,
P1ToP0ConstantP1EDGVectorLaplaceYCouplingOperator,
P1ToP0ConstantP1EDGVectorLaplaceZCouplingOperator >
P1ToP0ConstantP1EDGVectorLaplaceCouplingOperator;
} // namespace hyteg
target_sources( hyteg
PRIVATE
P0Function.hpp
P0Operator.hpp
)
......@@ -21,7 +21,7 @@
#pragma once
#include "hyteg/dgfunctionspace/DGBasisLinearLagrange_Example.hpp"
#include "hyteg/dgfunctionspace/DGFunction.hpp"
#include "hyteg/p0functionspace/P0Function.hpp"
#include "hyteg/p1dgefunctionspace/DGBasisEnriched.hpp"
#include "hyteg/p1functionspace/P1VectorFunction.hpp"
#include "hyteg/primitivestorage/PrimitiveStorage.hpp"
......@@ -45,7 +45,7 @@ class P1DGEFunction final : public Function< P1DGEFunction< ValueType > >
std::shared_ptr< P1VectorFunction< ValueType > > getConformingPart() const { return u_conforming_; }
std::shared_ptr< dg::DGFunction< ValueType > > getDiscontinuousPart() const { return u_discontinuous_; }
std::shared_ptr< P0Function< ValueType > > getDiscontinuousPart() const { return u_discontinuous_; }
void setBoundaryCondition( BoundaryCondition bc )
{
......@@ -58,8 +58,7 @@ class P1DGEFunction final : public Function< P1DGEFunction< ValueType > >
template < typename SenderType, typename ReceiverType >
void communicate( const uint_t& level ) const
{
u_conforming_->communicate( level );
u_discontinuous_->communicate( level );
WALBERLA_ABORT( "Not implemented." );
}
void add( const ValueType scalar, uint_t level, DoFType flag = All ) const
......@@ -135,7 +134,7 @@ class P1DGEFunction final : public Function< P1DGEFunction< ValueType > >
uint_t level,
DoFType flag = All ) const
{
std::vector< std::reference_wrapper< const dg::DGFunction< ValueType > > > dg_list;
std::vector< std::reference_wrapper< const P0Function< ValueType > > > dg_list;
std::vector< std::reference_wrapper< const P1VectorFunction< ValueType > > > p1_list;
for ( auto f : functions )
......@@ -163,7 +162,12 @@ class P1DGEFunction final : public Function< P1DGEFunction< ValueType > >
return result;
}
void enumerate( uint_t level ) const { enumerate( level, 0 ); }
void enumerate( uint_t level ) const
{
ValueType start = 0;
enumerate( level, start );
WALBERLA_UNUSED( start );
}
void enumerate( uint_t level, ValueType& offset ) const
{
......@@ -267,8 +271,8 @@ class P1DGEFunction final : public Function< P1DGEFunction< ValueType > >
std::vector< uint_t > vertexDoFIndices( numDofs );
std::vector< real_t > dofValues( numDofs );
auto dofs = u_discontinuous_->volumeDoFFunction()->dofMemory( faceID, level );
const auto memLayout = u_discontinuous_->volumeDoFFunction()->memoryLayout();
auto dofs = u_discontinuous_->getDGFunction()->volumeDoFFunction()->dofMemory( faceID, level );
const auto memLayout = u_discontinuous_->getDGFunction()->volumeDoFFunction()->memoryLayout();
for ( auto faceType : facedof::allFaceTypes )
{
......@@ -300,7 +304,7 @@ class P1DGEFunction final : public Function< P1DGEFunction< ValueType > >
std::shared_ptr< dg::DGBasisLinearLagrange_Example > basis_;
std::shared_ptr< P1VectorFunction< ValueType > > u_conforming_;
std::shared_ptr< dg::DGFunction< ValueType > > u_discontinuous_;
std::shared_ptr< P0Function< ValueType > > u_discontinuous_;
dg::DGBasisLinearLagrange_Example basis_conforming_;
DGBasisEnriched basis_discontinuous_;
......@@ -315,7 +319,7 @@ P1DGEFunction< ValueType >::P1DGEFunction( const std::string&
: Function< P1DGEFunction< ValueType > >( name, storage, minLevel, maxLevel )
, basis_{ std::make_shared< dg::DGBasisLinearLagrange_Example >() }
, u_conforming_{ std::make_shared< P1VectorFunction< ValueType > >( name, storage, minLevel, maxLevel, bc ) }
, u_discontinuous_{ std::make_shared< dg::DGFunction< ValueType > >( name, storage, minLevel, maxLevel, basis_, 0, bc ) }
, u_discontinuous_{ std::make_shared< P0Function< ValueType > >( name, storage, minLevel, maxLevel, bc ) }
, basis_conforming_()
, basis_discontinuous_()
{}
......
......@@ -20,40 +20,85 @@
#pragma once
#include "hyteg/dgfunctionspace/DGVectorLaplaceForm.hpp"
#include "hyteg/dgfunctionspace/DGVectorMassForm.hpp"
#include "hyteg/mixedoperators/P0ScalarToP1VectorOperator.hpp"
#include "hyteg/mixedoperators/P1VectorToP0ScalarOperator.hpp"
#include "hyteg/operators/Operator.hpp"
#include "hyteg/operators/VectorLaplaceOperator.hpp"
#include "hyteg/p0functionspace/P0Operator.hpp"
#include "hyteg/p1dgefunctionspace/P1DGEFunction.hpp"
#include "hyteg/dgfunctionspace/DGVectorMassForm.hpp"
#include "hyteg/p1functionspace/P1ConstantOperator.hpp"
namespace hyteg {
template < typename ValueType >
class P1DGEMassOperator final : public Operator< P1DGEFunction< real_t >, P1DGEFunction< real_t > > {
class P1DGEMassOperator final : public Operator< P1DGEFunction< real_t >, P1DGEFunction< real_t > >
{
public:
P1DGEMassOperator( const std::shared_ptr< PrimitiveStorage >& storage,
uint_t minLevel,
uint_t maxLevel);
P1DGEMassOperator( const std::shared_ptr< PrimitiveStorage >& storage, uint_t minLevel, uint_t maxLevel ) {}
void apply( const P1DGEFunction< real_t >& src,
const P1DGEFunction< real_t >& dst,
size_t level,
DoFType flag,
UpdateType updateType ) const override
size_t level,
DoFType flag,
UpdateType updateType ) const override
{
WALBERLA_ABORT( "Not implemented." );
}
void toMatrix( const std::shared_ptr< SparseMatrixProxy >& mat,
const P1DGEFunction< idx_t >& src,
const P1DGEFunction< idx_t >& dst,
const P1DGEFunction< idx_t >& src,
const P1DGEFunction< idx_t >& dst,
size_t level,
DoFType flag ) const override
{
WALBERLA_ABORT( "Not implemented." );
}
};
template < typename ValueType >
class P1DGELaplaceOperator final : public Operator< P1DGEFunction< real_t >, P1DGEFunction< real_t > >
{
public:
P1DGELaplaceOperator( const std::shared_ptr< PrimitiveStorage >& storage, uint_t minLevel, uint_t maxLevel )
: Operator< P1DGEFunction< real_t >, P1DGEFunction< real_t > >( storage, minLevel, maxLevel )
, cg_to_cg_coupling_( storage, minLevel, maxLevel )
, eg_to_cg_coupling_( storage, minLevel, maxLevel )
, cg_to_eg_coupling_( storage, minLevel, maxLevel )
, eg_to_eg_coupling_( storage, minLevel, maxLevel, std::make_shared< dg::DGVectorLaplaceFormEDGEDG >() )
{}
void apply( const P1DGEFunction< real_t >& src,
const P1DGEFunction< real_t >& dst,
size_t level,
DoFType flag,
UpdateType updateType ) const override
{
WALBERLA_ABORT( "Not implemented." );
}
private:
void toMatrix( const std::shared_ptr< SparseMatrixProxy >& mat,
const P1DGEFunction< idx_t >& src,
const P1DGEFunction< idx_t >& dst,
size_t level,
DoFType flag ) const override
{
communication::syncVectorFunctionBetweenPrimitives( *src.getConformingPart(), level );
communication::syncVectorFunctionBetweenPrimitives( *dst.getConformingPart(), level );
src.getDiscontinuousPart()->communicate( level );
dst.getDiscontinuousPart()->communicate( level );
cg_to_cg_coupling_.toMatrix( mat, *src.getConformingPart(), *dst.getConformingPart(), level, flag );
eg_to_cg_coupling_.toMatrix( mat, *src.getDiscontinuousPart(), *dst.getConformingPart(), level, flag );
cg_to_eg_coupling_.toMatrix( mat, *src.getConformingPart(), *dst.getDiscontinuousPart(), level, flag );
eg_to_eg_coupling_.toMatrix( mat, *src.getDiscontinuousPart(), *dst.getDiscontinuousPart(), level, flag );
}
private:
P1ConstantVectorLaplaceOperator cg_to_cg_coupling_;
P1ToP0ConstantP1EDGVectorLaplaceCouplingOperator cg_to_eg_coupling_;
P0ToP1ConstantP1EDGVectorLaplaceCouplingOperator eg_to_cg_coupling_;
P0Operator< dg::DGVectorLaplaceFormEDGEDG > eg_to_eg_coupling_;
};
}
} // namespace hyteg
......@@ -804,6 +804,9 @@ if (HYTEG_BUILD_WITH_PETSC)
waLBerla_compile_test(FILES p1dgefunctionspace/P1DGEFunctionspaceTest.cpp DEPENDS hyteg core)
waLBerla_execute_test(NAME P1DGEFunctionspaceTest)
waLBerla_compile_test(FILES p1dgefunctionspace/P1DGESymmetryTest.cpp DEPENDS hyteg core)
waLBerla_execute_test(NAME P1DGESymmetryTest)
endif()