Commit cdd1d92d authored by Tomáš Oberhuber's avatar Tomáš Oberhuber
Browse files

Refactoring gradient norm eoc unit test.

parent 279bc71b
Loading
Loading
Loading
Loading
+96 −86
Original line number Diff line number Diff line
@@ -19,6 +19,7 @@
#define	TNLFDMGRADIENTNORM_H

#include<operators/fdm/tnlForwardFiniteDifference.h>
#include<operators/geometric/tnlExactGradientNorm.h>

template< typename Mesh,
          template< typename, int, int, int, typename, typename > class DifferenceOperatorTemplate = tnlForwardFiniteDifference,
@@ -44,6 +45,10 @@ class tnlFDMGradientNorm< tnlGrid< 1,MeshReal, Device, MeshIndex >, DifferenceOp
   typedef Real RealType;
   typedef Device DeviceType;
   typedef Index IndexType;
   typedef tnlExactGradientNorm< 1, RealType > ExactOperatorType;
   
   template< typename MeshEntity = typename MeshType::Cell >
   using XDifferenceOperatorType = DifferenceOperatorTemplate< typename MeshEntity::MeshType, 1, 0, 0, Real, Index >;
   
   tnlFDMGradientNorm()
   : epsSquare( 0.0 ){}
@@ -62,11 +67,8 @@ class tnlFDMGradientNorm< tnlGrid< 1,MeshReal, Device, MeshIndex >, DifferenceOp
                    const MeshEntity& entity,
                    const Real& time = 0.0 ) const
   {
      DifferenceOperatorTemplate< typename MeshEntity::MeshType, 1, 0, 0, Real, Index > XDifferenceOperator;
      //const IndexType& cellIndex = entity.getIndex();
      //const typename MeshEntity::template NeighbourEntities< 1 >& neighbourEntities = entity.getNeighbourEntities();      
      //const typename MeshEntity::MeshType& mesh = entity.getMesh();
      const RealType u_x = XDifferenceOperator( u, entity );
      XDifferenceOperatorType< MeshEntity > XDifference;
      const RealType u_x = XDifference( u, entity );
      return sqrt( this->epsSquare + u_x * u_x );          
   }
                
@@ -97,6 +99,12 @@ class tnlFDMGradientNorm< tnlGrid< 2,MeshReal, Device, MeshIndex >, DifferenceOp
      typedef Real RealType;
      typedef Device DeviceType;
      typedef Index IndexType;
      typedef tnlExactGradientNorm< 2, RealType > ExactOperatorType;
      
      template< typename MeshEntity >
      using XDifferenceOperatorType = DifferenceOperatorTemplate< typename MeshEntity::MeshType, 1, 0, 0, Real, Index >;
      template< typename MeshEntity >
      using YDifferenceOperatorType = DifferenceOperatorTemplate< typename MeshEntity::MeshType, 0, 1, 0, Real, Index >;

      tnlFDMGradientNorm()
      : epsSquare( 0.0 ){}
@@ -117,14 +125,10 @@ class tnlFDMGradientNorm< tnlGrid< 2,MeshReal, Device, MeshIndex >, DifferenceOp
                       const MeshEntity& entity,
                       const Real& time = 0.0 ) const
      {
      DifferenceOperatorTemplate< typename MeshEntity::MeshType, 1, 0, 0, Real, Index > XDifferenceOperator;
      DifferenceOperatorTemplate< typename MeshEntity::MeshType, 0, 1, 0, Real, Index > YDifferenceOperator;
      //const IndexType& cellIndex = entity.getIndex();
      //const typename MeshEntity::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();      
      //const typename MeshEntity::MeshType& mesh = entity.getMesh();
      //const RealType& u_c = u[ cellIndex ];
      const RealType u_x = XDifferenceOperator( u, entity );
      const RealType u_y = YDifferenceOperator( u, entity );
         XDifferenceOperatorType< MeshEntity > XDifference;
         YDifferenceOperatorType< MeshEntity > YDifference;
         const RealType u_x = XDifference( u, entity );
         const RealType u_y = YDifference( u, entity );
         return sqrt( this->epsSquare + u_x * u_x + u_y * u_y );       
      }

@@ -157,6 +161,15 @@ class tnlFDMGradientNorm< tnlGrid< 3, MeshReal, Device, MeshIndex >, DifferenceO
      typedef Real RealType;
      typedef Device DeviceType;
      typedef Index IndexType;
      typedef tnlExactGradientNorm< 3, RealType > ExactOperatorType;
   
      template< typename MeshEntity >
      using XDifferenceOperatorType = DifferenceOperatorTemplate< typename MeshEntity::MeshType, 1, 0, 0, Real, Index >;
      template< typename MeshEntity >
      using YDifferenceOperatorType = DifferenceOperatorTemplate< typename MeshEntity::MeshType, 0, 1, 0, Real, Index >;
      template< typename MeshEntity >
      using ZDifferenceOperatorType = DifferenceOperatorTemplate< typename MeshEntity::MeshType, 0, 0, 1, Real, Index >;

   
      tnlFDMGradientNorm()
      : epsSquare( 0.0 ){}   
@@ -175,16 +188,13 @@ class tnlFDMGradientNorm< tnlGrid< 3, MeshReal, Device, MeshIndex >, DifferenceO
                       const MeshEntity& entity,
                       const Real& time = 0.0 ) const
      {
      DifferenceOperatorTemplate< typename MeshEntity::MeshType, 1, 0, 0, Real, Index > XDifferenceOperator;
      DifferenceOperatorTemplate< typename MeshEntity::MeshType, 0, 1, 0, Real, Index > YDifferenceOperator;
      DifferenceOperatorTemplate< typename MeshEntity::MeshType, 0, 0, 1, Real, Index > ZDifferenceOperator;
      //const IndexType& cellIndex = entity.getIndex();
      //const typename MeshEntity::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();      
      //const typename MeshEntity::MeshType& mesh = entity.getMesh();
      //const RealType& u_c = u[ cellIndex ];
      const RealType u_x = XDifferenceOperator( u, entity );
      const RealType u_y = YDifferenceOperator( u, entity );
      const RealType u_z = ZDifferenceOperator( u, entity );
         XDifferenceOperatorType< MeshEntity > XDifference;
         YDifferenceOperatorType< MeshEntity > YDifference;
         ZDifferenceOperatorType< MeshEntity > ZDifference;

         const RealType u_x = XDifference( u, entity );
         const RealType u_y = YDifference( u, entity );
         const RealType u_z = ZDifference( u, entity );
         return sqrt( this->epsSquare + u_x * u_x + u_y * u_y + u_z * u_z );             
      }

+1 −1
Original line number Diff line number Diff line
ADD_SUBDIRECTORY( diffusion )
ADD_SUBDIRECTORY( fdm )
#ADD_SUBDIRECTORY( geometric )
 No newline at end of file
ADD_SUBDIRECTORY( geometric )
 No newline at end of file
+6 −6
Original line number Diff line number Diff line
@@ -2,13 +2,13 @@ ADD_EXECUTABLE( tnlFDMGradientNormTest${mpiExt}${debugExt} tnlFDMGradientNormTes
TARGET_LINK_LIBRARIES( tnlFDMGradientNormTest${mpiExt}${debugExt} ${CPPUNIT_LIBRARIES}
                                                                 tnl${mpiExt}${debugExt}-0.1 )

ADD_EXECUTABLE( tnlTwoSidedGradientNormTest${mpiExt}${debugExt} tnlTwoSidedGradientNormTest.cpp )
TARGET_LINK_LIBRARIES( tnlTwoSidedGradientNormTest${mpiExt}${debugExt} ${CPPUNIT_LIBRARIES}
                                                                 tnl${mpiExt}${debugExt}-0.1 )
#ADD_EXECUTABLE( tnlTwoSidedGradientNormTest${mpiExt}${debugExt} tnlTwoSidedGradientNormTest.cpp )
#TARGET_LINK_LIBRARIES( tnlTwoSidedGradientNormTest${mpiExt}${debugExt} ${CPPUNIT_LIBRARIES}
#                                                                 tnl${mpiExt}${debugExt}-0.1 )

ADD_EXECUTABLE( tnlCoFVMGradientNormTest${mpiExt}${debugExt} tnlCoFVMGradientNormTest.cpp )
TARGET_LINK_LIBRARIES( tnlCoFVMGradientNormTest${mpiExt}${debugExt} ${CPPUNIT_LIBRARIES}
                                                                 tnl${mpiExt}${debugExt}-0.1 )
#ADD_EXECUTABLE( tnlCoFVMGradientNormTest${mpiExt}${debugExt} tnlCoFVMGradientNormTest.cpp )
#TARGET_LINK_LIBRARIES( tnlCoFVMGradientNormTest${mpiExt}${debugExt} ${CPPUNIT_LIBRARIES}
#                                                                 tnl${mpiExt}${debugExt}-0.1 )

if( BUILD_CUDA )                                                           
   CUDA_ADD_EXECUTABLE( tnlFDMGradientNormTest-cuda${mpiExt}${debugExt} ${headers} tnlFDMGradientNormTest.cu
+107 −113
Original line number Diff line number Diff line
@@ -24,156 +24,150 @@
#include <operators/fdm/tnlForwardFiniteDifference.h>
#include <operators/fdm/tnlCentralFiniteDifference.h>
#include "../../tnlUnitTestStarter.h"
#include "../tnlPDEOperatorEocTest.h"
#include "../tnlPDEOperatorEocUnitTest.h"

template< int Dimensions,
          typename Real,
          typename Device,
          typename Index,
          typename TestFunction >
class tnlPDEOperatorEocTestResult< 
   tnlFDMGradientNorm< tnlGrid< Dimensions, Real, Device, Index >,
                       tnlForwardFiniteDifference,
                       Real,
                       Index >,
   TestFunction >

template< typename ApproximateOperator >
class tnlFDMGradientNormEocTestResults
{
   public:       
      static Real getL1Eoc() { return ( Real ) 1.0; };
      static Real getL1Tolerance() { return ( Real ) 0.05; };
 
      static Real getL2Eoc() { return ( Real ) 1.0; };
      static Real getL2Tolerance() { return ( Real ) 0.05; };
      typedef typename ApproximateOperator::RealType RealType;
           
      static Real getMaxEoc() { return ( Real ) 1.0; };
      static Real getMaxTolerance() { return ( Real ) 0.05; };
      const RealType eoc[ 3 ] =       { 1.0,  1.0,  1.0 };
      const RealType tolerance[ 3 ] = { 0.05, 0.05, 0.05 };      
};

template< typename MeshType,
          typename RealType,
          typename IndexType >
class tnlFDMGradientNormEocTestResults< tnlCentralFiniteDifference< MeshType, 1, 0, 0, RealType, IndexType > >
{
   public:
      
      const RealType eoc[ 3 ] =       { 2.0,  2.0,  2.0 };
      const RealType tolerance[ 3 ] = { 0.05, 0.05, 0.05 };
};

template< int Dimensions,
          typename Real,
          typename Device,
          typename Index,
          typename TestFunction >
class tnlPDEOperatorEocTestResult< 
   tnlFDMGradientNorm< tnlGrid< Dimensions, Real, Device, Index >,
                       tnlBackwardFiniteDifference,
                       Real,
                       Index >,
   TestFunction >
template< typename ApproximateOperator,
          typename TestFunction,
          bool write = false,
          bool verbose = false >
class tnlFDMGradientNormTest
   : public tnlPDEOperatorEocTest< ApproximateOperator, TestFunction >,
     public tnlFDMGradientNormEocTestResults< typename ApproximateOperator::template XDifferenceOperatorType< typename ApproximateOperator::MeshType::Cell > >
{
   public:
      static Real getL1Eoc() { return ( Real ) 1.0; };
      static Real getL1Tolerance() { return ( Real ) 0.05; };
      
      static Real getL2Eoc() { return ( Real ) 1.0; };
      static Real getL2Tolerance() { return ( Real ) 0.05; };
      typedef ApproximateOperator ApproximateOperatorType;
      typedef typename ApproximateOperatorType::ExactOperatorType ExactOperatorType;
      typedef typename ApproximateOperator::MeshType MeshType;
      typedef typename ApproximateOperator::RealType RealType;
      typedef typename ApproximateOperator::IndexType IndexType;
      
      static Real getMaxEoc() { return ( Real ) 1.0; };
      static Real getMaxTolerance() { return ( Real ) 0.05; };
      const IndexType coarseMeshSize[ 3 ] = { 1024, 256, 64 };
      
};
      static tnlString getType()
      { 
         return tnlString( "tnlLinearDiffusionTest< " ) + 
                ApproximateOperator::getType() + ", " +
                TestFunction::getType() + " >";
      }
      
template< int Dimensions,
          typename Real,
          typename Device,
          typename Index,
          typename TestFunction >
class tnlPDEOperatorEocTestResult< 
   tnlFDMGradientNorm< tnlGrid< Dimensions, Real, Device, Index >,
                       tnlCentralFiniteDifference,
                       Real,
                       Index >,
   TestFunction >
      void setupTest()
      {
   public:
      static Real getL1Eoc() { return ( Real ) 2.0; };
      static Real getL1Tolerance() { return ( Real ) 0.05; };
         this->setupFunction();
      }
            
      static Real getL2Eoc() { return ( Real ) 2.0; };
      static Real getL2Tolerance() { return ( Real ) 0.05; };
      void getApproximationError( const IndexType meshSize,
                                  RealType errors[ 3 ] )
      {
         this->setupMesh( meshSize );
         this->performTest( this->approximateOperator,
                            this->exactOperator,
                            errors,
                            write,
                            verbose );

      static Real getMaxEoc() { return ( Real ) 2.0; };
      static Real getMaxTolerance() { return ( Real ) 0.1; };
      }
      
};
      void runUnitTest()
      {  
         RealType coarseErrors[ 3 ], fineErrors[ 3 ];
         this->getApproximationError( coarseMeshSize[ MeshType::getDimensionsCount() - 1 ], coarseErrors );
         this->getApproximationError( 2 * coarseMeshSize[ MeshType::getDimensionsCount() - 1 ], fineErrors );
         this->checkEoc( coarseErrors, fineErrors, this->eoc, this->tolerance, verbose );                            
      }
      
   protected:

      ApproximateOperator approximateOperator;
      
template< typename Mesh,
      ExactOperatorType exactOperator;

};

template< typename Operator,
          typename Function, 
          typename Operator,
          int MeshSize,
          bool WriteFunctions,
          bool Verbose >
bool testDifferenceOperator()
          bool write,
          bool verbose >
bool runTest()
{
   typedef tnlExactGradientNorm< Mesh::meshDimensions > ExactOperator;
   return tnlUnitTestStarter::run<
            tnlPDEOperatorEocTester< 
                Operator,
                ExactOperator,
                Function,
                typename Mesh::Cell,
                MeshSize,
                WriteFunctions,
                Verbose > >();
   
   typedef tnlFDMGradientNormTest< Operator, Function, write, verbose > OperatorTest;
#ifdef HAVE_CPPUNIT   
   if( ! tnlUnitTestStarter::run< tnlPDEOperatorEocUnitTest< OperatorTest > >() )
      return false;
   return true;
#endif      
}

template< typename Mesh,
          typename Function,
          int MeshSize,
          bool WriteFunctions,
          bool Verbose >
          bool write,
          bool verbose >
bool setDifferenceOperator()
{
   typedef tnlFDMGradientNorm< Mesh, tnlForwardFiniteDifference > ForwardGradientNorm;
   typedef tnlFDMGradientNorm< Mesh, tnlBackwardFiniteDifference > BackwardGradientNorm;
   typedef tnlFDMGradientNorm< Mesh, tnlCentralFiniteDifference > CentralGradientNorm;
   return ( testDifferenceOperator< Mesh, Function, ForwardGradientNorm,  MeshSize, WriteFunctions, Verbose >() &&
            testDifferenceOperator< Mesh, Function, BackwardGradientNorm, MeshSize, WriteFunctions, Verbose >() &&
            testDifferenceOperator< Mesh, Function, CentralGradientNorm,  MeshSize, WriteFunctions, Verbose >() );
   return ( runTest< ForwardGradientNorm, Function, write, verbose >() &&
            runTest< BackwardGradientNorm, Function, write, verbose >() &&
            runTest< CentralGradientNorm, Function, write, verbose >() );
}

template< typename Mesh,
          int MeshSize,
          bool WriteFunctions,
          bool Verbose >
bool setFunction()
          bool write,
          bool verbose >
bool setTestFunction()
{
   const int Dimensions = Mesh::meshDimensions;
   typedef tnlExpBumpFunction< Dimensions, double >  Function;
   return setDifferenceOperator< Mesh, Function, MeshSize, WriteFunctions, Verbose >();
   return setDifferenceOperator< Mesh, tnlExpBumpFunction< Mesh::getDimensionsCount(), double >, write, verbose >();
}

template< typename Device,
          int MeshSize,
          bool WriteFunctions,
          bool Verbose >
bool setGrid()
          bool write,
          bool verbose >
bool setMesh()
{
   typedef double MeshReal;
   typedef int MeshIndex;
   typedef tnlGrid< 1, MeshReal, Device, MeshIndex > Grid1D;
   typedef tnlGrid< 2, MeshReal, Device, MeshIndex > Grid2D;
   typedef tnlGrid< 3, MeshReal, Device, MeshIndex > Grid3D;
   return ( setFunction< Grid1D, MeshSize, WriteFunctions, Verbose >() &&
            setFunction< Grid2D, MeshSize, WriteFunctions, Verbose >() &&
            setFunction< Grid3D, MeshSize, WriteFunctions, Verbose >() );
   return ( setTestFunction< tnlGrid< 1, double, Device, int >, write, verbose >() &&
            setTestFunction< tnlGrid< 2, double, Device, int >, write, verbose >() &&
            setTestFunction< tnlGrid< 3, double, Device, int >, write, verbose >() );
}

int main( int argc, char* argv[] )
{
   const int meshSize( 32 );
   const bool writeFunctions( false );
   const bool verbose( true );
#ifdef HAVE_CPPUNIT
    return setGrid< tnlHost, meshSize, writeFunctions, verbose >();
#else
   const bool verbose( false );
   const bool write( false );
   
   if( ! setMesh< tnlHost, write, verbose  >() )
      return EXIT_FAILURE;
#ifdef HAVE_CUDA
   if( ! setMesh< tnlCuda, write, verbose >() )
      return EXIT_FAILURE;
#endif   
   return EXIT_SUCCESS;
}


#endif	/* TNLFDMGRADIENTNORMTEST_H */