Commit b291ac30 authored by Marcel Koch's avatar Marcel Koch
Browse files

allow only one cuda exec

parent 9cf90afa
......@@ -22,23 +22,18 @@
#include "core/math/Random.h"
#include "core/timing/Timer.h"
#include "hyteg/dataexport/VTKOutput.hpp"
#include "hyteg/mesh/MeshInfo.hpp"
#include "hyteg/misc/ExactStencilWeights.hpp"
#include "hyteg/p1functionspace/P1ConstantOperator.hpp"
#include "hyteg/p1functionspace/P1Function.hpp"
#include "hyteg/composites/P2P1TaylorHoodFunction.hpp"
#include "hyteg/composites/P2P1TaylorHoodStokesOperator.hpp"
#include "hyteg/ginkgo/GinkgoBlockSolver.hpp"
#include "hyteg/ginkgo/GinkgoSparseMatrixProxy.hpp"
#include "hyteg/ginkgo/GinkgoUtilities.hpp"
#include "hyteg/ginkgo/GinkgoVectorProxy.hpp"
#include "hyteg/mesh/MeshInfo.hpp"
#include "hyteg/p1functionspace/P1Function.hpp"
#include "hyteg/primitivestorage/SetupPrimitiveStorage.hpp"
#include "hyteg/primitivestorage/Visualization.hpp"
#include "hyteg/primitivestorage/loadbalancing/SimpleBalancer.hpp"
#include "hyteg/composites/P2P1TaylorHoodFunction.hpp"
#include "hyteg/composites/P2P1TaylorHoodStokesOperator.hpp"
#ifndef HYTEG_BUILD_WITH_GINKGO
WALBERLA_ABORT( "This test only works with Ginkgo enabled. Please enable it via -DHYTEG_BUILD_WITH_GINKGO=ON" )
WALBERLA_ABORT( "This test only works with Ginkgo enabled. Please enable it via -DHYTEG_BUILD_WITH_GINKGO=ON" )
#endif
using walberla::real_t;
......@@ -63,77 +58,81 @@ void petscSolveTest( std::shared_ptr< const gko::Executor > exec,
std::shared_ptr< PrimitiveStorage > storage = std::make_shared< PrimitiveStorage >( setupStorage );
writeDomainPartitioningVTK( storage, "../../output", "P2P1Stokes3DPetscSolve_Domain" );
hyteg::P2P1TaylorHoodFunction< real_t > x( "x", storage, level, level );
hyteg::P2P1TaylorHoodFunction< real_t > x_exact( "x_exact", storage, level, level );
hyteg::P2P1TaylorHoodFunction< real_t > b( "b", storage, level, level );
hyteg::P2P1TaylorHoodFunction< real_t > btmp( "btmp", storage, level, level );
hyteg::P2P1TaylorHoodFunction< real_t > err( "err", storage, level, level );
hyteg::P2P1TaylorHoodFunction< real_t > residuum( "res", storage, level, level );
hyteg::P2P1TaylorHoodFunction< real_t > nullspace( "nullspace", storage, level, level );
hyteg::P2P1TaylorHoodFunction< real_t > x( "x", storage, level, level );
hyteg::P2P1TaylorHoodFunction< real_t > x_exact( "x_exact", storage, level, level );
hyteg::P2P1TaylorHoodFunction< real_t > b( "b", storage, level, level );
hyteg::P2P1TaylorHoodFunction< real_t > btmp( "btmp", storage, level, level );
hyteg::P2P1TaylorHoodFunction< real_t > err( "err", storage, level, level );
hyteg::P2P1TaylorHoodFunction< real_t > residuum( "res", storage, level, level );
hyteg::P2P1TaylorHoodFunction< real_t > nullspace( "nullspace", storage, level, level );
hyteg::P2P1TaylorHoodFunction< int32_t > numerator( "numerator", storage, level, level );
numerator.enumerate( level );
hyteg::P2P1TaylorHoodStokesOperator A( storage, level, level );
std::function< real_t( const hyteg::Point3D& ) > exactU = []( const hyteg::Point3D& xx ) { return real_c(20) * xx[0] * xx[1] * xx[1] * xx[1]; };
std::function< real_t( const hyteg::Point3D& ) > exactV = []( const hyteg::Point3D& xx ) { return real_c(5) * xx[0] * xx[0] * xx[0] * xx[0] - real_c(5) * xx[1] * xx[1] * xx[1] * xx[1]; };
std::function< real_t( const hyteg::Point3D& ) > exactW = []( const hyteg::Point3D& ) { return real_c(0); };
std::function< real_t( const hyteg::Point3D& ) > exactP = []( const hyteg::Point3D& xx ) { return real_c(60) * std::pow( xx[0], 2.0 ) * xx[1] - real_c(20) * std::pow( xx[1], 3.0 ); };
std::function< real_t( const hyteg::Point3D& ) > zero = []( const hyteg::Point3D& ) { return real_c(0); };
std::function< real_t( const hyteg::Point3D& ) > ones = []( const hyteg::Point3D& ) { return real_c(1); };
hyteg::P2P1TaylorHoodStokesOperator A( storage, level, level );
b.uvw.interpolate( { exactU, exactV, exactW }, level, DirichletBoundary );
b.p.interpolate( zero, level, All );
std::function< real_t( const hyteg::Point3D& ) > exactU = []( const hyteg::Point3D& xx ) {
return real_c( 20 ) * xx[0] * xx[1] * xx[1] * xx[1];
};
std::function< real_t( const hyteg::Point3D& ) > exactV = []( const hyteg::Point3D& xx ) {
return real_c( 5 ) * xx[0] * xx[0] * xx[0] * xx[0] - real_c( 5 ) * xx[1] * xx[1] * xx[1] * xx[1];
};
std::function< real_t( const hyteg::Point3D& ) > exactW = []( const hyteg::Point3D& ) { return real_c( 0 ); };
std::function< real_t( const hyteg::Point3D& ) > exactP = []( const hyteg::Point3D& xx ) {
return real_c( 60 ) * std::pow( xx[0], 2.0 ) * xx[1] - real_c( 20 ) * std::pow( xx[1], 3.0 );
};
std::function< real_t( const hyteg::Point3D& ) > zero = []( const hyteg::Point3D& ) { return real_c( 0 ); };
std::function< real_t( const hyteg::Point3D& ) > ones = []( const hyteg::Point3D& ) { return real_c( 1 ); };
x.uvw.interpolate( { exactU, exactV, exactW }, level, DirichletBoundary );
b.uvw.interpolate( { exactU, exactV, exactW }, level, DirichletBoundary );
b.p.interpolate( zero, level, All );
x_exact.uvw.interpolate( { exactU, exactV, exactW }, level );
x_exact.p.interpolate( exactP, level );
x.uvw.interpolate( { exactU, exactV, exactW }, level, DirichletBoundary );
hyteg::vertexdof::projectMean( x_exact.p, level );
x_exact.uvw.interpolate( { exactU, exactV, exactW }, level );
x_exact.p.interpolate( exactP, level );
nullspace.p.interpolate( ones, level, All );
hyteg::vertexdof::projectMean( x_exact.p, level );
uint_t localDoFs1 = hyteg::numberOfLocalDoFs< P2P1TaylorHoodFunctionTag >( *storage, level );
uint_t globalDoFs1 = hyteg::numberOfGlobalDoFs< P2P1TaylorHoodFunctionTag >( *storage, level );
uint_t globalDoFsvelocity = 3 * hyteg::numberOfGlobalDoFs< P2FunctionTag >( *storage, level );
nullspace.p.interpolate( ones, level, All );
WALBERLA_LOG_INFO( "localDoFs: " << localDoFs1 << " globalDoFs: " << globalDoFs1
<< ", global velocity dofs: " << globalDoFsvelocity );
uint_t localDoFs1 = hyteg::numberOfLocalDoFs< P2P1TaylorHoodFunctionTag >( *storage, level );
uint_t globalDoFs1 = hyteg::numberOfGlobalDoFs< P2P1TaylorHoodFunctionTag >( *storage, level );
uint_t globalDoFsvelocity = 3 * hyteg::numberOfGlobalDoFs< P2FunctionTag >( *storage, level );
GinkgoBlockSolver< P2P1TaylorHoodStokesOperator > solver( storage, level, exec );
WALBERLA_LOG_INFO_ON_ROOT( "localDoFs: " << localDoFs1 << " globalDoFs: " << globalDoFs1
<< ", global velocity dofs: " << globalDoFsvelocity );
walberla::WcTimer timer;
solver.solve(A, x, b, level);
timer.end();
GinkgoBlockSolver< P2P1TaylorHoodStokesOperator > solver( storage, level, exec );
hyteg::vertexdof::projectMean( x.p, level );
walberla::WcTimer timer;
solver.solve( A, x, b, level );
timer.end();
WALBERLA_LOG_INFO_ON_ROOT( "time was: " << timer.last() );
A.apply( x, residuum, level, hyteg::Inner );
hyteg::vertexdof::projectMean( x.p, level );
err.assign( {1.0, -1.0}, {x, x_exact}, level );
WALBERLA_LOG_INFO_ON_ROOT( "time was: " << timer.last() );
A.apply( x, residuum, level, hyteg::Inner );
real_t discr_l2_err = std::sqrt( err.dotGlobal( err, level ) / (real_t) globalDoFs1 );
real_t discr_l2_err_1_u = std::sqrt( err.uvw[0].dotGlobal( err.uvw[0], level ) / (real_t) globalDoFs1 );
real_t discr_l2_err_1_v = std::sqrt( err.uvw[1].dotGlobal( err.uvw[1], level ) / (real_t) globalDoFs1 );
real_t discr_l2_err_1_w = std::sqrt( err.uvw[2].dotGlobal( err.uvw[2], level ) / (real_t) globalDoFs1 );
real_t discr_l2_err_1_p = std::sqrt( err.p.dotGlobal( err.p, level ) / (real_t) globalDoFs1 );
real_t residuum_l2_1 = std::sqrt( residuum.dotGlobal( residuum, level ) / (real_t) globalDoFs1 );
err.assign( { 1.0, -1.0 }, { x, x_exact }, level );
WALBERLA_LOG_INFO_ON_ROOT( "discrete L2 error = " << discr_l2_err );
WALBERLA_LOG_INFO_ON_ROOT( "discrete L2 error u = " << discr_l2_err_1_u );
WALBERLA_LOG_INFO_ON_ROOT( "discrete L2 error v = " << discr_l2_err_1_v );
WALBERLA_LOG_INFO_ON_ROOT( "discrete L2 error w = " << discr_l2_err_1_w );
WALBERLA_LOG_INFO_ON_ROOT( "discrete L2 error p = " << discr_l2_err_1_p );
WALBERLA_LOG_INFO_ON_ROOT( "residuum 1 = " << residuum_l2_1 );
real_t discr_l2_err = std::sqrt( err.dotGlobal( err, level ) / (real_t) globalDoFs1 );
real_t discr_l2_err_1_u = std::sqrt( err.uvw[0].dotGlobal( err.uvw[0], level ) / (real_t) globalDoFs1 );
real_t discr_l2_err_1_v = std::sqrt( err.uvw[1].dotGlobal( err.uvw[1], level ) / (real_t) globalDoFs1 );
real_t discr_l2_err_1_w = std::sqrt( err.uvw[2].dotGlobal( err.uvw[2], level ) / (real_t) globalDoFs1 );
real_t discr_l2_err_1_p = std::sqrt( err.p.dotGlobal( err.p, level ) / (real_t) globalDoFs1 );
real_t residuum_l2_1 = std::sqrt( residuum.dotGlobal( residuum, level ) / (real_t) globalDoFs1 );
//WALBERLA_CHECK_LESS( residuum_l2_1, resEps );
//WALBERLA_CHECK_LESS( discr_l2_err_1_u + discr_l2_err_1_v + discr_l2_err_1_w, errEpsUSum );
//WALBERLA_CHECK_LESS( discr_l2_err_1_p, errEpsP);
WALBERLA_LOG_INFO_ON_ROOT( "discrete L2 error = " << discr_l2_err );
WALBERLA_LOG_INFO_ON_ROOT( "discrete L2 error u = " << discr_l2_err_1_u );
WALBERLA_LOG_INFO_ON_ROOT( "discrete L2 error v = " << discr_l2_err_1_v );
WALBERLA_LOG_INFO_ON_ROOT( "discrete L2 error w = " << discr_l2_err_1_w );
WALBERLA_LOG_INFO_ON_ROOT( "discrete L2 error p = " << discr_l2_err_1_p );
WALBERLA_LOG_INFO_ON_ROOT( "residuum 1 = " << residuum_l2_1 );
auto tt = storage->getTimingTree()->getReduced().getCopyWithRemainder();
//WALBERLA_CHECK_LESS( residuum_l2_1, resEps );
//WALBERLA_CHECK_LESS( discr_l2_err_1_u + discr_l2_err_1_v + discr_l2_err_1_w, errEpsUSum );
//WALBERLA_CHECK_LESS( discr_l2_err_1_p, errEpsP);
}
}
......@@ -145,20 +144,19 @@ int main( int argc, char* argv[] )
walberla::Environment walberlaEnv( argc, argv );
walberla::MPIManager::instance()->useWorldComm();
// PETScManager petscManager( &argc, &argv );
//
// printPETScVersionNumberString();
auto level = argc > 2 ? std::stoi( argv[2] ) : 0;
for ( auto tag : tag_list )
{
if(walberla::MPIManager::instance()->rank() != 0 && (tag == exec_tag::cuda || tag == exec_tag::hip || tag == exec_tag::dpcpp)){
tag = exec_tag::reference;
}
auto exec = get_executor( tag );
if ( exec )
if ( walberla::mpi::allReduce(bool(exec), walberla::mpi::Operation::LOGICAL_AND) )
{
try
{
WALBERLA_LOG_INFO_ON_ROOT("Running test for " << get_executor_name(exec) << " executor");
WALBERLA_LOG_INFO_ON_ROOT( "Running test for " << get_executor_name( exec ) << " executor" );
petscSolveTest( exec,
level,
hyteg::MeshInfo::fromGmshFile( "../../data/meshes/3D/cube_center_at_origin_24el.msh" ),
......@@ -167,7 +165,7 @@ int main( int argc, char* argv[] )
0.33 );
} catch ( const gko::NotImplemented& e )
{
std::cout << e.what() << " for executor " << get_executor_name( exec ) << std::endl;
WALBERLA_LOG_INFO_ON_ROOT( e.what() << " for executor " << get_executor_name( exec ) );
}
}
}
......
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