Skip to content

GitLab

Commit bd1b3a13 authored by Marcus Mohr's avatar Marcus Mohr
Browse files

Commit starts re-factoring of matrix creation (see issue #125) 

We start introduction member functions toMatrix() to operator classes.
These are intended to replace the createMatrix() free-functions.
parent fa74468d
Pipeline #34415 failed with stages
in 11 minutes and 7 seconds
Hide whitespace changes
Inline Side-by-side
src/hyteg/elementwiseoperators/P2ElementwiseOperator.hpp
View file @ bd1b3a13
......@@ -34,8 +34,8 @@
#include "hyteg/p1functionspace/VertexDoFMacroFace.hpp"
#include "hyteg/p2functionspace/P2Elements.hpp"
#include "hyteg/p2functionspace/P2Function.hpp"
#include "hyteg/sparseassembly/SparseMatrixProxy.hpp"
#include "hyteg/solvers/Smoothables.hpp"
#include "hyteg/sparseassembly/SparseMatrixProxy.hpp"
namespace hyteg {
......@@ -130,6 +130,24 @@ class P2ElementwiseOperator : public Operator< P2Function< real_t >, P2Function<
DoFType flag ) const;
#endif
/// Assemble operator as sparse matrix.
///
/// \param mat a sparse matrix proxy
/// \param src P2Function for determining column indices
/// \param dst P2Function for determining row indices
/// \param level level in mesh hierarchy for which local operator is to be assembled
/// \param flag ignored
///
/// \note src and dst are legal to and often will be the same function object
void toMatrix( const std::shared_ptr< SparseMatrixProxy >& mat,
const P2Function< matIdx_t >& src,
const P2Function< matIdx_t >& dst,
uint_t level,
DoFType flag ) const override
{
this->assembleLocalMatrix( mat, src, dst, level, flag );
}
P2Form getForm() const;
private:
......
src/hyteg/functions/BlockFunction.hpp
View file @ bd1b3a13
......@@ -41,6 +41,9 @@ class BlockFunction
public:
typedef value_t ValueType;
template < typename VType >
using FunctionType = BlockFunction< VType >;
BlockFunction( const std::string name )
: functionName_( name )
{}
......
src/hyteg/operators/BlockOperator.hpp
View file @ bd1b3a13
......@@ -27,6 +27,7 @@
#include "hyteg/p1functionspace/P1ConstantOperator.hpp"
#include "hyteg/p2functionspace/P2ConstantOperator.hpp"
#include "hyteg/solvers/Smoothables.hpp"
#include "hyteg/sparseassembly/SparseMatrixProxy.hpp"
namespace hyteg {
......@@ -96,7 +97,7 @@ class BlockOperator : public Operator< srcBlockFunc_t, dstBlockFunc_t >,
ConstructorArguments... args )
{
const auto op = std::make_shared< OperatorWrapper< OperatorType > >( storage, minLevel, maxLevel, args... );
setSubOperator(i, j, op);
setSubOperator( i, j, op );
}
const std::shared_ptr< GenericOperator< value_t > > getSubOperator( uint_t i, uint_t j )
......@@ -134,16 +135,16 @@ class BlockOperator : public Operator< srcBlockFunc_t, dstBlockFunc_t >,
if ( colIdx == rowIdx )
continue;
subOper_[rowIdx][colIdx]->apply( x[colIdx], (*tmp_)[rowIdx], level, flag, Add );
subOper_[rowIdx][colIdx]->apply( x[colIdx], ( *tmp_ )[rowIdx], level, flag, Add );
}
// calculate new rhs: tmp1_[r] = b_r - sum_{c != r} A_rc x_c
(*tmp_)[rowIdx].assign( { +1, -1 }, { b[rowIdx], (*tmp_)[rowIdx] }, level, flag );
( *tmp_ )[rowIdx].assign( {+1, -1}, {b[rowIdx], ( *tmp_ )[rowIdx]}, level, flag );
// diagonal part:
auto D_rr = subOper_[rowIdx][rowIdx];
if ( auto* op_diag = dynamic_cast< GSSmoothable< GenericFunction< value_t > >* >( D_rr.get() ) )
{
op_diag->smooth_gs( x[rowIdx], (*tmp_)[rowIdx], level, flag );
op_diag->smooth_gs( x[rowIdx], ( *tmp_ )[rowIdx], level, flag );
}
}
}
......@@ -162,6 +163,24 @@ class BlockOperator : public Operator< srcBlockFunc_t, dstBlockFunc_t >,
throw std::runtime_error( "For GaussSeidel src and dst functions must coincide" );
}
void toMatrix( const std::shared_ptr< SparseMatrixProxy >& mat,
const typename srcBlockFunc_t::template FunctionType< matIdx_t >& src,
const typename dstBlockFunc_t::template FunctionType< matIdx_t >& dst,
size_t level,
DoFType flag ) const override
{
for ( uint_t i = 0; i < nRows_; i++ )
{
for ( uint_t j = 0; j < nCols_; j++ )
{
if ( subOper_[i][j] != nullptr )
{
subOper_[i][j]->toMatrix( mat, src[j], dst[i], level, flag );
}
}
}
}
protected:
std::vector< std::vector< std::shared_ptr< GenericOperator< value_t > > > > subOper_;
uint_t nRows_;
......
src/hyteg/operators/GenericOperator.hpp
View file @ bd1b3a13
......@@ -24,6 +24,7 @@
#include "hyteg/boundary/BoundaryConditions.hpp"
#include "hyteg/primitivestorage/PrimitiveStorage.hpp"
#include "hyteg/sparseassembly/SparseMatrixProxy.hpp"
#include "hyteg/types/flags.hpp"
namespace hyteg {
......@@ -64,6 +65,12 @@ class GenericOperator
const GenericFunction< value_t >& dst,
size_t level,
DoFType flag ) const = 0;
virtual void toMatrix( const std::shared_ptr< SparseMatrixProxy >& mat,
const GenericFunction< matIdx_t >& src,
const GenericFunction< matIdx_t >& dst,
size_t level,
DoFType flag ) const = 0;
};
} // namespace hyteg
src/hyteg/operators/Operator.hpp
View file @ bd1b3a13
......@@ -20,89 +20,101 @@
#pragma once
#include <memory>
#include <core/DataTypes.h>
#include <core/timing/TimingTree.h>
#include <memory>
#include "hyteg/functions/FunctionTraits.hpp"
#include "hyteg/primitivestorage/PrimitiveStorage.hpp"
#include "hyteg/sparseassembly/SparseMatrixProxy.hpp"
namespace hyteg {
template< typename SourceFunction, typename DestinationFunction >
template < typename SourceFunction, typename DestinationFunction >
class Operator
{
public:
Operator(const std::shared_ptr<PrimitiveStorage> & storage, uint_t minLevel, uint_t maxLevel)
: storage_(storage), minLevel_(minLevel), maxLevel_(maxLevel)
{
if(storage->getTimingTree()){
enableTiming(storage->getTimingTree());
}
}
typedef SourceFunction srcType;
typedef DestinationFunction dstType;
virtual ~Operator() = default;
void enableTiming( const std::shared_ptr< walberla::WcTimingTree > & timingTree ) { timingTree_ = timingTree; }
const std::shared_ptr< PrimitiveStorage > getStorage() const { return storage_; }
public:
Operator( const std::shared_ptr< PrimitiveStorage >& storage, uint_t minLevel, uint_t maxLevel )
: storage_( storage )
, minLevel_( minLevel )
, maxLevel_( maxLevel )
{
if ( storage->getTimingTree() )
{
enableTiming( storage->getTimingTree() );
}
}
typedef SourceFunction srcType;
typedef DestinationFunction dstType;
virtual ~Operator() = default;
void enableTiming( const std::shared_ptr< walberla::WcTimingTree >& timingTree ) { timingTree_ = timingTree; }
const std::shared_ptr< PrimitiveStorage > getStorage() const { return storage_; }
uint_t getMinLevel() const { return minLevel_; }
uint_t getMaxLevel() const { return maxLevel_; }
// See the free function of the same name. Method might be revived, if we implement a base class for composite
// operators
//
// uint_t getNumberOfGlobalDoFCouplings( uint_t level ) const
// {
// uint_t nCouplings =
// indexing::countLocalDoFCouplings< typename srcType::Tag, typename dstType::Tag >( storage_, level );
// walberla::mpi::allReduceInplace( nCouplings, walberla::mpi::SUM, walberla::mpi::MPIManager::instance()->comm() );
// return nCouplings;
// };
// Needed in VectorToVectorOperator
virtual void apply( const SourceFunction& src,
const DestinationFunction& dst,
size_t level,
DoFType flag,
UpdateType updateType = Replace ) const
{
WALBERLA_ABORT( "Problem with inheritance, this function should have been overridden in the derived class!" );
};
virtual void toMatrix( const std::shared_ptr< SparseMatrixProxy >& mat,
const typename srcType::template FunctionType< matIdx_t >& src,
const typename dstType::template FunctionType< matIdx_t >& dst,
size_t level,
DoFType flag ) const
{
WALBERLA_ABORT( "toMatrix() not implemented in derived class!" );
}
uint_t getMinLevel() const { return minLevel_; }
uint_t getMaxLevel() const { return maxLevel_; }
// See the free function of the same name. Method might be revived, if we implement a base class for composite
// operators
//
// uint_t getNumberOfGlobalDoFCouplings( uint_t level ) const
// {
// uint_t nCouplings =
// indexing::countLocalDoFCouplings< typename srcType::Tag, typename dstType::Tag >( storage_, level );
// walberla::mpi::allReduceInplace( nCouplings, walberla::mpi::SUM, walberla::mpi::MPIManager::instance()->comm() );
// return nCouplings;
// };
// Needed in VectorToVectorOperator
virtual void apply( const SourceFunction& src,
const DestinationFunction& dst,
size_t level,
DoFType flag,
UpdateType updateType = Replace ) const
{
WALBERLA_ABORT( "Problem with inheritance, this function should have been overridden in the derived class!" );
};
protected:
const std::shared_ptr< PrimitiveStorage > storage_;
const uint_t minLevel_;
const uint_t maxLevel_;
protected:
const std::shared_ptr< PrimitiveStorage > storage_;
const uint_t minLevel_;
const uint_t maxLevel_;
std::shared_ptr< walberla::WcTimingTree > timingTree_;
std::shared_ptr< walberla::WcTimingTree > timingTree_;
protected:
void startTiming( const std::string & timerString ) const
{
if ( timingTree_ )
{
timingTree_->start( "Operator " + FunctionTrait< SourceFunction >::getTypeName() + " to " + FunctionTrait< DestinationFunction >::getTypeName() );
timingTree_->start( timerString );
}
}
void stopTiming ( const std::string & timerString ) const
{
if ( timingTree_ )
{
timingTree_->stop( timerString );
timingTree_->stop( "Operator " + FunctionTrait< SourceFunction >::getTypeName() + " to " + FunctionTrait< DestinationFunction >::getTypeName() );
}
}
void startTiming( const std::string& timerString ) const
{
if ( timingTree_ )
{
timingTree_->start( "Operator " + FunctionTrait< SourceFunction >::getTypeName() + " to " +
FunctionTrait< DestinationFunction >::getTypeName() );
timingTree_->start( timerString );
}
}
void stopTiming( const std::string& timerString ) const
{
if ( timingTree_ )
{
timingTree_->stop( timerString );
timingTree_->stop( "Operator " + FunctionTrait< SourceFunction >::getTypeName() + " to " +
FunctionTrait< DestinationFunction >::getTypeName() );
}
}
};
} // namespace hyteg
src/hyteg/operators/OperatorWrapper.hpp
View file @ bd1b3a13
......@@ -25,6 +25,7 @@
#include "hyteg/functions/FunctionTraits.hpp"
#include "hyteg/operators/GenericOperator.hpp"
#include "hyteg/solvers/Smoothables.hpp"
#include "hyteg/sparseassembly/SparseMatrixProxy.hpp"
namespace hyteg {
......@@ -43,7 +44,10 @@ class OperatorWrapper final
/// Constructor that constructs the operator which the class wraps itself around
template < typename... ConstructorArguments >
OperatorWrapper( const std::shared_ptr< PrimitiveStorage >& storage, uint_t minLevel, uint_t maxLevel, ConstructorArguments... args )
OperatorWrapper( const std::shared_ptr< PrimitiveStorage >& storage,
uint_t minLevel,
uint_t maxLevel,
ConstructorArguments... args )
{
wrappedOper_ = std::make_unique< oper_t >( storage, minLevel, maxLevel, args... );
};
......@@ -87,7 +91,7 @@ class OperatorWrapper final
{
const auto& dst_unwrapped = dst.template unwrap< typename oper_t::srcType >();
const auto& rhs_unwrapped = rhs.template unwrap< typename oper_t::dstType >();
op_gs->smooth_gs(dst_unwrapped, rhs_unwrapped, level, flag );
op_gs->smooth_gs( dst_unwrapped, rhs_unwrapped, level, flag );
}
else
{
......@@ -100,6 +104,19 @@ class OperatorWrapper final
}
}
void toMatrix( const std::shared_ptr< SparseMatrixProxy >& mat,
const GenericFunction< matIdx_t >& src,
const GenericFunction< matIdx_t >& dst,
size_t level,
DoFType flag ) const
{
wrappedOper_->toMatrix( mat,
src.template unwrap< typename oper_t::srcType::template FunctionType< matIdx_t > >(),
dst.template unwrap< typename oper_t::dstType::template FunctionType< matIdx_t > >(),
level,
flag );
};
private:
std::unique_ptr< oper_t > wrappedOper_;
};
......
src/hyteg/operators/VectorToVectorOperator.hpp
View file @ bd1b3a13
......@@ -86,6 +86,24 @@ class VectorToVectorOperator : public Operator< SrcVecFuncKind< ValueType >, Dst
return subOper_[i][j];
}
void toMatrix( const std::shared_ptr< SparseMatrixProxy >& mat,
const SrcVecFuncKind< matIdx_t >& src,
const DstVecFuncKind< matIdx_t >& dst,
size_t level,
DoFType flag ) const
{
for ( uint_t i = 0; i < dim_; i++ )
{
for ( uint_t j = 0; j < dim_; j++ )
{
if ( subOper_[i][j] != nullptr )
{
subOper_[i][j]->toMatrix( mat, src[j], dst[i], level, flag );
}
}
}
}
protected:
std::vector< std::vector< std::shared_ptr< scalarOpType > > > subOper_;
uint_t dim_;
......
src/hyteg/p1functionspace/P1ConstantOperator.hpp
View file @ bd1b3a13
......@@ -21,16 +21,22 @@
#include <array>
#include "hyteg/operators/Operator.hpp"
#include "hyteg/forms/P1LinearCombinationForm.hpp"
#include "hyteg/forms/P2LinearCombinationForm.hpp"
#include "hyteg/forms/P2RowSumForm.hpp"
#include "hyteg/forms/form_fenics_base/P1FenicsForm.hpp"
#include "hyteg/memory/LevelWiseMemory.hpp"
#include "hyteg/memory/StencilMemory.hpp"
#include "hyteg/operators/Operator.hpp"
#include "hyteg/p1functionspace/P1Function.hpp"
#include "hyteg/p1functionspace/VertexDoFIndexing.hpp"
#include "hyteg/solvers/Smoothables.hpp"
#include "hyteg/sparseassembly/SparseMatrixProxy.hpp"
// This include can be removed once the implementation of createMatrix() was moved into toMatrix()
#ifdef HYTEG_BUILD_WITH_PETSC
#include "hyteg/p1functionspace/P1Petsc.hpp"
#endif
namespace hyteg {
......@@ -146,6 +152,18 @@ class P1ConstantOperator : public Operator< P1Function< real_t >, P1Function< re
return cellStencilID_;
}
// Remove guard once the implementation of createMatrix() here
#ifdef HYTEG_BUILD_WITH_PETSC
void toMatrix( const std::shared_ptr< SparseMatrixProxy >& mat,
const P1Function< matIdx_t >& src,
const P1Function< matIdx_t >& dst,
size_t level,
DoFType flag ) const override
{
hyteg::petsc::createMatrix( *this, src, dst, mat, level, flag );
}
#endif
private:
void assembleStencils();
......
src/hyteg/p2functionspace/P2ConstantOperator.hpp
View file @ bd1b3a13
......@@ -29,6 +29,11 @@
#include "hyteg/p2functionspace/P2Function.hpp"
#include "hyteg/solvers/Smoothables.hpp"
// This include can be removed once the implementation of createMatrix() was moved into toMatrix()
#ifdef HYTEG_BUILD_WITH_PETSC
#include "hyteg/p2functionspace/P2Petsc.hpp"
#endif
namespace hyteg {
using walberla::real_t;
......@@ -101,6 +106,18 @@ class P2ConstantOperator : public Operator< P2Function< real_t >, P2Function< re
size_t level,
DoFType flag ) const override;
// Remove guard once the implementation of createMatrix() here
#ifdef HYTEG_BUILD_WITH_PETSC
void toMatrix( const std::shared_ptr< SparseMatrixProxy >& mat,
const P2Function< matIdx_t >& src,
const P2Function< matIdx_t >& dst,
size_t level,
DoFType flag ) const override
{
hyteg::petsc::createMatrix( *this, src, dst, mat, level, flag );
}
#endif
private:
void smooth_sor_macro_vertices( const P2Function< real_t >& dst,
const P2Function< real_t >& rhs,
......
src/hyteg/petsc/PETScSparseMatrix.hpp
View file @ bd1b3a13
......@@ -87,6 +87,19 @@ class PETScSparseMatrix
assembled_ = true;
}
inline void createMatrixFromOperator_newAPI( const OperatorType& op,
uint_t level,
const FunctionType< PetscInt >& numerator,
DoFType flag = All )
{
auto proxy = std::make_shared< PETScSparseMatrixProxy >( mat );
op.toMatrix( proxy, numerator, numerator, level, flag );
MatAssemblyBegin( mat, MAT_FINAL_ASSEMBLY );
MatAssemblyEnd( mat, MAT_FINAL_ASSEMBLY );
assembled_ = true;
}
inline bool createMatrixFromOperatorOnce( const OperatorType& op,
uint_t level,
const FunctionType< PetscInt >& numerator,
......
src/hyteg/types/flags.hpp
View file @ bd1b3a13
......@@ -84,5 +84,8 @@ enum class CycleType
WCYCLE
};
typedef int32_t matIdx_t;
} // namespace hyteg
tests/hyteg/Petsc/PetscMatrixAssemblyTest.cpp
View file @ bd1b3a13
......@@ -26,10 +26,10 @@
#include "hyteg/composites/P1EpsilonStokesOperator.hpp"
#include "hyteg/composites/P1P1UzawaDampingFactorEstimationOperator.hpp"
#include "hyteg/dataexport/VTKOutput.hpp"
#include "hyteg/elementwiseoperators/DiagonalNonConstantOperator.hpp"
#include "hyteg/elementwiseoperators/ElementwiseOperatorPetsc.hpp"
#include "hyteg/elementwiseoperators/P1ElementwiseOperator.hpp"
#include "hyteg/elementwiseoperators/P2ElementwiseOperator.hpp"
#include "hyteg/elementwiseoperators/DiagonalNonConstantOperator.hpp"
#include "hyteg/functions/FunctionTraits.hpp"
#include "hyteg/geometry/AnnulusMap.hpp"
#include "hyteg/mesh/MeshInfo.hpp"
......@@ -54,12 +54,14 @@
#include "hyteg/composites/P2P2StabilizedStokesOperator.hpp"
#include "hyteg/composites/P2P2UnstableStokesOperator.hpp"
#include "hyteg/composites/UnsteadyDiffusion.hpp"
#include "hyteg/operators/BlockOperator.hpp"
#include "hyteg/operators/VectorMassOperator.hpp"
#include "hyteg/petsc/PETScExportOperatorMatrix.hpp"
#include "hyteg/petsc/PETScManager.hpp"
#include "hyteg/petsc/PETScSparseMatrix.hpp"
#include "hyteg/primitivestorage/SetupPrimitiveStorage.hpp"
#include "hyteg/primitivestorage/loadbalancing/SimpleBalancer.hpp"
#include "hyteg/composites/P2P1TaylorHoodBlockFunction.hpp"
/// This test checks whether we can assemble a sparse matrix
/// (for/with PETSc) for the given operators. It is mostly a
......@@ -79,6 +81,22 @@ void printTestHdr( std::string msg )
}
template < class operType >
void testAssembly_newAPI( std::shared_ptr< PrimitiveStorage >& storage, uint_t level, std::string tag )
{
WALBERLA_LOG_INFO_ON_ROOT( " * " << tag );
PETScManager petscManager;
PETScSparseMatrix< operType, operType::srcType::template FunctionType > matrix( storage, level, tag.c_str() );
typename operType::srcType::template FunctionType< PetscInt > enumerator( "enumerator", storage, level, level );
enumerator.enumerate( level );
operType oper( storage, level, level );
matrix.createMatrixFromOperator_newAPI( oper, level, enumerator, All );
}
template < class operType >
void testAssembly( std::shared_ptr< PrimitiveStorage >& storage, uint_t level, std::string tag )
{
WALBERLA_LOG_INFO_ON_ROOT( " * " << tag );
......@@ -95,7 +113,10 @@ void testAssembly( std::shared_ptr< PrimitiveStorage >& storage, uint_t level, s
// Version for operators the require a form in their ctors
template < class operType, class formType >
void testAssembly( std::shared_ptr< PrimitiveStorage >& storage, uint_t level, std::shared_ptr< formType >& form, std::string tag )
void testAssembly( std::shared_ptr< PrimitiveStorage >& storage,
uint_t level,
std::shared_ptr< formType >& form,
std::string tag )
{
WALBERLA_LOG_INFO_ON_ROOT( " * " << tag );
......@@ -162,7 +183,9 @@ void testAssembly( uint_t level, std::string tag )
// Specialised one case versions
// -------------------------------
template <>
void testAssembly< P1ConstantUnsteadyDiffusionOperator >( std::shared_ptr< PrimitiveStorage >& storage, uint_t level, std::string tag )
void testAssembly< P1ConstantUnsteadyDiffusionOperator >( std::shared_ptr< PrimitiveStorage >& storage,
uint_t level,
std::string tag )
{
WALBERLA_LOG_INFO_ON_ROOT( " * " << tag );
......@@ -205,7 +228,7 @@ int main( int argc, char* argv[] )
printTestHdr( "Testing P1 Operators" );
testAssembly< P1ConstantMassOperator >( storage, level, "P1ConstantOperator" );
testAssembly_newAPI< P1ConstantMassOperator >( storage, level, "P1ConstantOperator" );
testAssembly< P1ConstantMassOperator_new >( storage, level, "P1ConstantOperator_new" );
testAssembly< P1SurrogateMassOperator >( storage, level, "P1SurrogateOperator" );
testAssembly< P1BlendingMassOperator_new >( storage, level, "P1VariableOperator_new" );
......@@ -220,7 +243,7 @@ int main( int argc, char* argv[] )
printTestHdr( "Testing P2 Operators" );
testAssembly< P2ConstantMassOperator >( storage, level, "P2ConstantOperator" );
testAssembly_newAPI< P2ConstantMassOperator >( storage, level, "P2ConstantOperator" );
WALBERLA_LOG_INFO_ON_ROOT( "Skipping: (since assembly doesn't compile)" );
WALBERLA_LOG_INFO_ON_ROOT( " * P2SurrogateVariableOperator" );
......@@ -234,7 +257,8 @@ int main( int argc, char* argv[] )
// Low-level Operators
// ---------------------
WALBERLA_LOG_INFO_ON_ROOT( "" << rule << "\n LOW-LEVEL OPERATORS\n" << rule );
typedef EdgeDoFOperator< P2FenicsForm< p2_mass_cell_integral_0_otherwise, p2_tet_mass_cell_integral_0_otherwise > > E2EMassOperator;
typedef EdgeDoFOperator< P2FenicsForm< p2_mass_cell_integral_0_otherwise, p2_tet_mass_cell_integral_0_otherwise > >
E2EMassOperator;
testAssembly< E2EMassOperator >( storage, level, "EdgeDoFOperator" );
WALBERLA_LOG_INFO_ON_ROOT( "Skipping: (since assembly doesn't compile)" );
......@@ -255,7 +279,8 @@ int main( int argc, char* argv[] )
auto p2MassFormHyTeG = std::make_shared< P2Form_mass >();
std::shared_ptr< P2RowSumForm > lumpedMassFormP2HyTeG = std::make_shared< P2RowSumForm >( p2MassFormHyTeG );
testAssembly< P2BlendingLumpedDiagonalOperator, P2RowSumForm >( storage, level, lumpedMassFormP2HyTeG, "DiagonalNonConstantOperator" );
testAssembly< P2BlendingLumpedDiagonalOperator, P2RowSumForm >(
storage, level, lumpedMassFormP2HyTeG, "DiagonalNonConstantOperator" );
testAssembly< P2P1ElementwiseBlendingStokesOperator >( storage, level, "P2P1ElementwiseBlendingStokesOperator" );
......@@ -280,6 +305,7 @@ int main( int argc, char* argv[] )
testAssembly< P2P1TaylorHoodStokesOperator >( storage, level, "P2P1TaylorHoodStokesOperator" );
testAssembly< P2P2StabilizedStokesOperator >( storage, level, "P2P2StabilizedStokesOperator" );
testAssembly< P2P1TaylorHoodStokesOperator >( storage, level, "P2P1TaylorHoodStokesOperator" );
WALBERLA_LOG_INFO_ON_ROOT( "Skipping: (since assembly doesn't compile)" );
WALBERLA_LOG_INFO_ON_ROOT( " * P1PolynomialBlendingStokesOperator" );
WALBERLA_LOG_INFO_ON_ROOT( " * P1BlendingStokesOperator" );
......@@ -288,8 +314,8 @@ int main( int argc, char* argv[] )
WALBERLA_LOG_INFO_ON_ROOT( " * P2P2UnstableStokesOperator" );
WALBERLA_LOG_INFO_ON_ROOT( " * P2P1BlendingTaylorHoodStokesOperator" );
WALBERLA_LOG_INFO_ON_ROOT( " * P2P1SurrogateTaylorHoodStokesOperator" );
// testAssembly< P1BlendingStokesOperator >( storage, level, "P1BlendingStokesOperator" );
// testAssembly< P1StokesBlockLaplaceOperator >( storage, level, "P1StokesBlockLaplaceOperator" );
// testAssembly< P1BlendingStokesOperator >( storage, level, "P1BlendingStokesOperator" );
// testAssembly< P1EpsilonStokesOperator >( storage, level, "P1EpsilonStokesOperator" );
// testAssembly< P2P2UnstableStokesOperator >( storage, level, "P2P2UnstableStokesOperator" );
// testAssembly< P2P1BlendingTaylorHoodStokesOperator >( storage, level, "P2P1BlendingTaylorHoodStokesOperator" );
......@@ -321,10 +347,36 @@ int main( int argc, char* argv[] )