Commit 068368aa authored by Marcus Mohr's avatar Marcus Mohr
Browse files

Implements PETSc vector <-> function conversion for GenericFunction

The commit introduces implementations of
- createVectorFromFunction()
- createFunctionFromVector()
for functions of type GenericFunction< PetscReal >.

The two free functions reside in the new file GenericFunctionPetsc.hpp.
Because of the peculiarities of GenericFunctions these free functions
are, however, only wrappers that call the corresponding one of the two
new member functions
- GenericFunction::toVector()
- GenericFunction::fromVector()
As usual these are pure virtual and the actual work happens inside
the member functions of the FunctionWrapper class.

We extend the FunctionWrapperTest to check compilation and execution
for various kinds of wrapped functions.

Note that due to instantiation requirements the commit also implements
- createVectorFromFunction()
- createFunctionFromVector()
for functions of type DGFunction. However, these implementations are
pseudo only and will abort.
parent f7da1fc5
/*
* Copyright (c) 2021 Marcus Mohr.
*
* This file is part of HyTeG
* (see https://i10git.cs.fau.de/hyteg/hyteg).
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*/
#pragma once
#include "hyteg/petsc/PETScWrapper.hpp"
#include "hyteg/sparseassembly/SparseMatrixProxy.hpp"
#include "hyteg/sparseassembly/VectorProxy.hpp"
#ifdef HYTEG_BUILD_WITH_PETSC
namespace hyteg {
namespace petsc {
inline void createVectorFromFunction( const DGFunction< PetscReal >& function,
const DGFunction< PetscInt >& numerator,
const std::shared_ptr< VectorProxy >& vec,
uint_t level,
DoFType flag )
{
WALBERLA_ABORT( "Congrats :( You have detected another unimplemented feature of DGFunction" );
}
inline void createFunctionFromVector( const DGFunction< PetscReal >& function,
const DGFunction< PetscInt >& numerator,
const std::shared_ptr< VectorProxy >& vec,
uint_t level,
DoFType flag )
{
WALBERLA_ABORT( "Congrats :( You have detected another unimplemented feature of DGFunction" );
}
} // namespace petsc
} // namespace hyteg
#endif
......@@ -24,6 +24,17 @@
#include "hyteg/functions/FunctionTraits.hpp"
#include "hyteg/functions/GenericFunction.hpp"
#include "hyteg/sparseassembly/VectorProxy.hpp"
// A whole lot of includes, so that createVectorFromFunction below has
// a valid prototype for all possible cases
#include "hyteg/dgfunctionspace/DGPetsc.hpp"
#include "hyteg/edgedofspace/EdgeDoFPetsc.hpp"
#include "hyteg/p1functionspace/P1Petsc.hpp"
#include "hyteg/p2functionspace/P2Petsc.hpp"
// only needed for using PetscInt in to/fromVector() below!
#include "hyteg/petsc/PETScWrapper.hpp"
namespace hyteg {
......@@ -172,6 +183,43 @@ class FunctionWrapper final : public GenericFunction< typename FunctionTrait< fu
return numberOfLocalDoFs< typename FunctionTrait< WrappedFuncType >::Tag >( *storage, level );
}
#ifdef HYTEG_BUILD_WITH_PETSC
/// conversion to/from linear algebra representation
/// @{
void toVector( const GenericFunction< PetscInt >& numerator,
const std::shared_ptr< VectorProxy >& vec,
uint_t level,
DoFType flag ) const
{
if constexpr ( std::is_same< value_t, PetscReal >::value )
{
using numer_t = typename func_t::template FunctionType< PetscInt >;
petsc::createVectorFromFunction( *wrappedFunc_, numerator.template unwrap< numer_t >(), vec, level, flag );
}
else
{
WALBERLA_ABORT( "FunctionWrapper::toVector() only works for ValueType being identical to PetscReal" );
}
};
void fromVector( const GenericFunction< PetscInt >& numerator,
const std::shared_ptr< VectorProxy >& vec,
uint_t level,
DoFType flag ) const
{
if constexpr ( std::is_same< value_t, PetscReal >::value )
{
using numer_t = typename func_t::template FunctionType< PetscInt >;
petsc::createFunctionFromVector( *wrappedFunc_, numerator.template unwrap< numer_t >(), vec, level, flag );
}
else
{
WALBERLA_ABORT( "FunctionWrapper::fromVector() only works for ValueType being identical to PetscReal" );
}
};
/// @}
#endif
private:
std::unique_ptr< func_t > wrappedFunc_;
};
......
......@@ -23,7 +23,9 @@
#include "core/DataTypes.h"
#include "hyteg/boundary/BoundaryConditions.hpp"
#include "hyteg/petsc/PETScWrapper.hpp"
#include "hyteg/primitivestorage/PrimitiveStorage.hpp"
#include "hyteg/sparseassembly/VectorProxy.hpp"
#include "hyteg/types/flags.hpp"
namespace hyteg {
......@@ -125,6 +127,21 @@ class GenericFunction
virtual void enumerate( uint_t level ) const = 0;
virtual void enumerate( uint_t level, value_t& offset ) const = 0;
#ifdef HYTEG_BUILD_WITH_PETSC
/// conversion to/from linear algebra representation
/// @{
virtual void toVector( const GenericFunction< PetscInt >& numerator,
const std::shared_ptr< VectorProxy >& vec,
uint_t level,
DoFType flag ) const = 0;
virtual void fromVector( const GenericFunction< PetscInt >& numerator,
const std::shared_ptr< VectorProxy >& vec,
uint_t level,
DoFType flag ) const = 0;
/// @}
#endif
};
// Special version of numberOfLocalDoFs for GenericFunctions
......
/*
* Copyright (c) 2021 Marcus Mohr.
*
* This file is part of HyTeG
* (see https://i10git.cs.fau.de/hyteg/hyteg).
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*/
#pragma once
#include "core/DataTypes.h"
#include "hyteg/functions/GenericFunction.hpp"
#include "hyteg/petsc/PETScWrapper.hpp"
#include "hyteg/sparseassembly/VectorProxy.hpp"
#include "hyteg/types/flags.hpp"
namespace hyteg {
namespace petsc {
inline void createVectorFromFunction( const GenericFunction< PetscReal >& function,
const GenericFunction< PetscInt >& numerator,
const std::shared_ptr< VectorProxy >& vec,
uint_t level,
DoFType flag )
{
function.toVector( numerator, vec, level, flag );
}
inline void createFunctionFromVector( const GenericFunction< PetscReal >& function,
const GenericFunction< PetscInt >& numerator,
const std::shared_ptr< VectorProxy >& vec,
uint_t level,
DoFType flag )
{
function.fromVector( numerator, vec, level, flag );
}
} // namespace petsc
} // namespace hyteg
......@@ -114,6 +114,7 @@ class PETScVector
protected:
MPI_Comm petscCommunicator_;
Vec vec;
};
} // namespace hyteg
......
......@@ -18,6 +18,8 @@
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*/
#include "hyteg/functions/FunctionWrapper.hpp"
#include "core/DataTypes.h"
#include "core/Environment.h"
#include "core/debug/CheckFunctions.h"
......@@ -26,32 +28,59 @@
#include "hyteg/facedofspace/FaceDoFFunction.hpp"
#include "hyteg/functions/FunctionProperties.hpp"
#include "hyteg/functions/GenericFunctionPetsc.hpp"
#include "hyteg/mesh/MeshInfo.hpp"
#include "hyteg/p1functionspace/P1Function.hpp"
#include "hyteg/p1functionspace/P1VectorFunction.hpp"
#include "hyteg/p2functionspace/P2Function.hpp"
#include "hyteg/p2functionspace/P2VectorFunction.hpp"
#include "hyteg/petsc/PETScManager.hpp"
#include "hyteg/petsc/PETScVector.hpp"
#include "hyteg/primitivestorage/SetupPrimitiveStorage.hpp"
// Perform some basic test on the FunctionWrapper class
#include "hyteg/functions/FunctionWrapper.hpp"
using namespace hyteg;
template < typename func_t >
void wrapFunction( const std::shared_ptr< PrimitiveStorage >& storage,
size_t minLevel,
size_t maxLevel,
uint_t dim )
void wrapFunction( const std::shared_ptr< PrimitiveStorage >& storage, size_t minLevel, size_t maxLevel, uint_t dim )
{
std::string kind = FunctionTrait< func_t >::getTypeName();
FunctionWrapper< func_t > wrappedFunc( kind, storage, minLevel, maxLevel );
WALBERLA_LOG_INFO_ON_ROOT( "-> wrapped a '" << kind << "'" );
WALBERLA_CHECK_EQUAL( dim, wrappedFunc.getDimension() );
std::string kind = FunctionTrait< func_t >::getTypeName();
FunctionWrapper< func_t > wrappedFunc( kind, storage, minLevel, maxLevel );
WALBERLA_LOG_INFO_ON_ROOT( "-> wrapped a '" << kind << "'" );
WALBERLA_CHECK_EQUAL( dim, wrappedFunc.getDimension() );
}
#ifdef HYTEG_BUILD_WITH_PETSC
template < template < typename > class FunctionKind >
void testPETScConversion( const std::shared_ptr< PrimitiveStorage >& storage )
{
WALBERLA_LOG_INFO_ON_ROOT( "-> Running test for " << FunctionTrait< FunctionKind< real_t > >::getTypeName() );
uint_t level = 4;
FunctionWrapper< FunctionKind< real_t > > src( "testing", storage, level, level );
FunctionWrapper< FunctionKind< real_t > > dst( "testing", storage, level, level );
FunctionWrapper< FunctionKind< PetscInt > > numerator( "numerator", storage, level, level );
std::function< real_t( const hyteg::Point3D& ) > expression = []( const hyteg::Point3D& x ) {
real_t value;
value = std::sin( real_c( 3 ) * x[0] ) + real_c( 0.5 ) * x[1] * x[1];
return value;
};
numerator.enumerate( level );
src.interpolate( expression, level );
PETScVector< real_t, GenericFunction > vector( src, numerator, level, All );
vector.createFunctionFromVector( dst, numerator, level, All );
dst.assign( {real_c( 1 ), real_c( -1 )}, {dst, src}, level, All );
// real_t diff = dst.getMaxComponentMagnitude( level, All );
// WALBERLA_CHECK_FLOAT_EQUAL( diff, real_c( 0 ) );
// WALBERLA_LOG_INFO_ON_ROOT( " max. difference = " << diff );
}
#endif
int main( int argc, char* argv[] )
{
......@@ -171,5 +200,17 @@ int main( int argc, char* argv[] )
WALBERLA_CHECK_FLOAT_EQUAL( p2vec[0].getMaxMagnitude( maxLevel ), real_c( 4 ) );
WALBERLA_LOG_INFO_ON_ROOT( "P2VecFunc.interpolate() -> check" );
#ifdef HYTEG_BUILD_WITH_PETSC
WALBERLA_LOG_INFO_ON_ROOT( "--------------------------------" );
WALBERLA_LOG_INFO_ON_ROOT( " Function <-> Vector Conversion" );
WALBERLA_LOG_INFO_ON_ROOT( "--------------------------------" );
PETScManager petscManager( &argc, &argv );
testPETScConversion< P1Function >( storage );
testPETScConversion< P2Function >( storage );
testPETScConversion< EdgeDoFFunction >( storage );
testPETScConversion< P1VectorFunction >( storage );
testPETScConversion< P2VectorFunction >( storage );
#endif
return EXIT_SUCCESS;
}
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment