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

Merge branch 'hamilton-jacobi' of... 

Merge branch 'hamilton-jacobi' of geraldine.fjfi.cvut.cz:/local/projects/tnl/tnl into hamilton-jacobi

Conflicts:
	src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/CMakeLists.txt
	src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/MainBuildConfig.h
	src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/tnl-direct-eikonal-solver.h
	src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/tnlDirectEikonalMethodsBase_impl.h
	src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/tnlDirectEikonalProblem_impl.h
	src/Tools/tnl-diff.cpp
parents 0e415b53 223c69a7

src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/CMakeLists.txt

+7 −2
Original line number Diff line number Diff line
@@ -10,8 +10,13 @@ set( tnl_hamilton_jacobi_SOURCES 
#ADD_EXECUTABLE(tnl-hamilton-jacobi main.cpp)

#target_link_libraries (tnl-hamilton-jacobi tnl )



IF( BUILD_CUDA )

    CUDA_ADD_EXECUTABLE(tnl-direct-eikonal-solver tnl-direct-eikonal-solver.cu )

    target_link_libraries (tnl-direct-eikonal-solver tnl )

ELSE(  BUILD_CUDA )               

    ADD_EXECUTABLE(tnl-direct-eikonal-solver tnl-direct-eikonal-solver.cpp )

    target_link_libraries (tnl-direct-eikonal-solver tnl )

ENDIF( BUILD_CUDA )





INSTALL( TARGETS #tnl-hamilton-jacobi

src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/MainBuildConfig.h

+26 −33
Original line number Diff line number Diff line 
/***************************************************************************

                          MainBuildConfig.h  -  description

                          HamiltonJacobiBuildConfigTag.h  -  description

                             -------------------

    begin                : Jul 7, 2014

    copyright            : (C) 2014 by Tomas Oberhuber

    email                : tomas.oberhuber@fjfi.cvut.cz

 ***************************************************************************/



/***************************************************************************

 *                                                                         *

 *   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 2 of the License, or     *

 *   (at your option) any later version.                                   *

 *                                                                         *

 ***************************************************************************/

/* See Copyright Notice in tnl/Copyright */



#pragma once



#include <TNL/Solvers/BuildConfigTags.h>



#ifndef MAINBUILDCONFIG_H_

#define MAINBUILDCONFIG_H_

namespace TNL {



#include <solvers/tnlBuildConfigTags.h>

class HamiltonJacobiBuildConfig {};



class MainBuildConfig

{

   public:

namespace Solvers {   



      static void print() {std::cerr << "MainBuildConfig" <<std::endl; }

};



/****

 * Turn off support for float and long double.

 */

template<> struct tnlConfigTagReal< MainBuildConfig, float > { enum { enabled = false }; };

template<> struct tnlConfigTagReal< MainBuildConfig, long double > { enum { enabled = false }; };

//template<> struct ConfigTagReal< HamiltonJacobiBuildConfig, float > { enum { enabled = false }; };

template<> struct ConfigTagReal< HamiltonJacobiBuildConfig, long double > { enum { enabled = false }; };



/****

 * Turn off support for short int and long int indexing.

 */

template<> struct tnlConfigTagIndex< MainBuildConfig, short int >{ enum { enabled = false }; };

template<> struct tnlConfigTagIndex< MainBuildConfig, long int >{ enum { enabled = false }; };

template<> struct ConfigTagIndex< HamiltonJacobiBuildConfig, short int >{ enum { enabled = false }; };

template<> struct ConfigTagIndex< HamiltonJacobiBuildConfig, long int >{ enum { enabled = false }; };



/****

 * Use of tnlGrid is enabled for allowed dimensions and Real, Device and Index types.

 * Use of Grid is enabled for allowed dimensions and Real, Device and Index types.

 */

template< int Dimensions, typename Real, typename Device, typename Index >

   struct tnlConfigTagMesh< MainBuildConfig, tnlGrid< Dimensions, Real, Device, Index > >

      { enum { enabled = tnlConfigTagDimensions< MainBuildConfig, Dimensions >::enabled  &&

                         tnlConfigTagReal< MainBuildConfig, Real >::enabled &&

                         tnlConfigTagDevice< MainBuildConfig, Device >::enabled &&

                         tnlConfigTagIndex< MainBuildConfig, Index >::enabled }; };

/*template< int Dimension, typename Real, typename Device, typename Index >

   struct ConfigTagMesh< HamiltonJacobiBuildConfig, Meshes::Grid< Dimension, Real, Device, Index > >

      { enum { enabled = ConfigTagDimension< HamiltonJacobiBuildConfig, Dimension >::enabled  &&

                         ConfigTagReal< HamiltonJacobiBuildConfig, Real >::enabled &&

                         ConfigTagDevice< HamiltonJacobiBuildConfig, Device >::enabled &&

                         ConfigTagIndex< HamiltonJacobiBuildConfig, Index >::enabled }; };*/



/****

 * Please, chose your preferred time discretisation  here.

 */

template<> struct tnlConfigTagTimeDiscretisation< MainBuildConfig, tnlExplicitTimeDiscretisationTag >{ enum { enabled = true }; };

template<> struct tnlConfigTagTimeDiscretisation< MainBuildConfig, tnlSemiImplicitTimeDiscretisationTag >{ enum { enabled = false}; };

template<> struct tnlConfigTagTimeDiscretisation< MainBuildConfig, tnlImplicitTimeDiscretisationTag >{ enum { enabled = false }; };

template<> struct ConfigTagTimeDiscretisation< HamiltonJacobiBuildConfig, ExplicitTimeDiscretisationTag >{ enum { enabled = true }; };

template<> struct ConfigTagTimeDiscretisation< HamiltonJacobiBuildConfig, SemiImplicitTimeDiscretisationTag >{ enum { enabled = true }; };

template<> struct ConfigTagTimeDiscretisation< HamiltonJacobiBuildConfig, ImplicitTimeDiscretisationTag >{ enum { enabled = false }; };



/****

 * Only the Runge-Kutta-Merson solver is enabled by default.

 */

template<> struct tnlConfigTagExplicitSolver< MainBuildConfig, tnlExplicitEulerSolverTag >{ enum { enabled = false }; };

//template<> struct ConfigTagExplicitSolver< HamiltonJacobiBuildConfig, ExplicitEulerSolverTag >{ enum { enabled = false }; };



#endif /* MAINBUILDCONFIG_H_ */
 
} // namespace Solvers

} // namespace TNL

src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/tnl-direct-eikonal-solver.h

+6 −4
Original line number Diff line number Diff line
@@ -20,14 +20,16 @@ 
#include <TNL/Functions/MeshFunction.h>

#include <TNL/Meshes/Grid.h>

#include "tnlDirectEikonalProblem.h"

#include "MainBuildConfig.h"



using namespace TNL;



//typedef tnlDefaultBuildMeshConfig BuildConfig;

typedef Solvers::FastBuildConfigTag BuildConfig;

//typedef Solvers::FastBuildConfig BuildConfig;

typedef HamiltonJacobiBuildConfig BuildConfig;



template< typename MeshConfig >

class tnlDirectEikonalSolverConfig

class DirectEikonalSolverConfig

{

   public:

      static void configSetup( Config::ConfigDescription& config )
@@ -44,7 +46,7 @@ template< typename Real, 
          typename MeshConfig,

          typename SolverStarter,

          typename CommunicatorType >

class tnlDirectEikonalSolverSetter

class DirectEikonalSolverSetter

{

   public:
@@ -66,7 +68,7 @@ class tnlDirectEikonalSolverSetter 


int main( int argc, char* argv[] )

{

   if( ! Solvers::Solver< tnlDirectEikonalSolverSetter, tnlDirectEikonalSolverConfig, BuildConfig >::run( argc, argv ) )

   if( ! Solvers::Solver< DirectEikonalSolverSetter, DirectEikonalSolverConfig, BuildConfig >::run( argc, argv ) )

      return EXIT_FAILURE;

   return EXIT_SUCCESS;

}

src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/tnl-run-fsm-eoc-test

+7 −8
Original line number Diff line number Diff line
@@ -2,7 +2,7 @@ 


device="host"

dimensions="2D 3D"

dimensions="2D"

#dimensions="3D"

sizes1D="16 32 64 128 256 512 1024 2048 4096"

#sizes1D="256"

sizes2D="16 32 64 128 256 512 1024"
@@ -12,7 +12,7 @@ testFunctions="paraboloid" 
snapshotPeriod=0.1

finalTime=1.5

solverName="tnl-direct-eikonal-solver"

#solverName="gdb --args tnl-direct-eikonal-solver-dbg"

#solverName="gdb --args tnl-direct-eikonal-solver-dbg --catch-exceptions no"

#



setupTestFunction()
@@ -59,7 +59,7 @@ setupGrid() 
setInitialCondition()

{

   testFunction=$1

   tnl-init --test-function ${testFunction}-sdf \

   tnl-init --test-function ${testFunction} \

            --output-file initial-u.tnl \

            --amplitude ${amplitude} \

            --wave-length ${waveLength} \
@@ -78,7 +78,7 @@ setInitialCondition() 
            --radius ${radius}

   

    tnl-init --test-function ${testFunction}-sdf \

            --output-file final-u.tnl \

            --output-file exact-u.tnl \

            --amplitude ${amplitude} \

            --wave-length ${waveLength} \

            --wave-length-x ${waveLengthX} \
@@ -108,7 +108,6 @@ solve() 
                 --time-step 10 \

                 --time-step-order 2 \

                 --discrete-solver ${discreteSolver} \

                 --merson-adaptivity 1.0e-5 \

                 --min-iterations 20 \

                 --convergence-residue 1.0e-12 \

                 --snapshot-period ${snapshotPeriod} \
@@ -118,7 +117,7 @@ solve() 
computeError()

{

   tnl-diff --mesh mesh.tnl \

            --input-files final-u.tnl u-*.tnl \

            --input-files aux-final.tnl exact-u.tnl \

            --mode sequence \

            --snapshot-period ${snapshotPeriod} \

            --output-file errors.txt \
@@ -172,8 +171,8 @@ runTest() 
               echo "========================================================================="

               echo "===                   STARTING THE SOLVER                             ==="

               echo "========================================================================="

               solve explicit merson

               #solve semi-implicit gmres

               #solve explicit merson

               solve semi-implicit cg

               mv computation-in-progress computation-done

            fi            

            echo "========================================================================="

src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/tnlDirectEikonalMethodsBase.h

+82 −17
Original line number Diff line number Diff line
@@ -9,6 +9,7 @@ 


#include <TNL/Meshes/Grid.h>

#include <TNL/Functions/MeshFunction.h>

#include <TNL/Devices/Cuda.h>



using namespace TNL;
@@ -30,15 +31,21 @@ class tnlDirectEikonalMethodsBase< Meshes::Grid< 1, Real, Device, Index > > 
      typedef Index IndexType;

      typedef Functions::MeshFunction< MeshType > MeshFunctionType;

      typedef Functions::MeshFunction< MeshType, 1, bool > InterfaceMapType;

      using MeshFunctionPointer = SharedPointer< MeshFunctionType >;

      using InterfaceMapPointer = SharedPointer< InterfaceMapType >;

      

      void initInterface( const MeshFunctionType& input,

                          MeshFunctionType& output,

                          InterfaceMapType& interfaceMap );

      void initInterface( const MeshFunctionPointer& input,

                          MeshFunctionPointer& output,

                          InterfaceMapPointer& interfaceMap );

      

      template< typename MeshEntity >

      void updateCell( MeshFunctionType& u,

                       const MeshEntity& cell );

      __cuda_callable__ void updateCell( MeshFunctionType& u,

                                         const MeshEntity& cell,

                                         const RealType velocity = 1.0  );

      

      __cuda_callable__ bool updateCell( volatile Real sArray[18],

                                         int thri, const Real h,

                                         const Real velocity = 1.0 );

};
@@ -48,21 +55,27 @@ template< typename Real, 
class tnlDirectEikonalMethodsBase< Meshes::Grid< 2, Real, Device, Index > >

{

   public:

      

      typedef Meshes::Grid< 2, Real, Device, Index > MeshType;

      typedef Real RealType;

      typedef Device DevcieType;

      typedef Index IndexType;

      typedef Functions::MeshFunction< MeshType > MeshFunctionType;

      typedef Functions::MeshFunction< MeshType, 2, bool > InterfaceMapType;

      using MeshFunctionPointer = SharedPointer< MeshFunctionType >;

      using InterfaceMapPointer = SharedPointer< InterfaceMapType >;      



      void initInterface( const MeshFunctionType& input,

                          MeshFunctionType& output,

                          InterfaceMapType& interfaceMap );

      void initInterface( const MeshFunctionPointer& input,

                          MeshFunctionPointer& output,

                          InterfaceMapPointer& interfaceMap );

      

      template< typename MeshEntity >

      void updateCell( MeshFunctionType& u,

                       const MeshEntity& cell );

      __cuda_callable__ void updateCell( MeshFunctionType& u,

                                         const MeshEntity& cell,

                                         const RealType velocity = 1.0 );

      

      __cuda_callable__ bool updateCell( volatile Real sArray[18][18],

                                         int thri, int thrj, const Real hx, const Real hy,

                                         const Real velocity = 1.0 );

};



template< typename Real,
@@ -71,21 +84,73 @@ template< typename Real, 
class tnlDirectEikonalMethodsBase< Meshes::Grid< 3, Real, Device, Index > >

{

   public:

      

      typedef Meshes::Grid< 3, Real, Device, Index > MeshType;

      typedef Real RealType;

      typedef Device DevcieType;

      typedef Index IndexType;

      typedef Functions::MeshFunction< MeshType > MeshFunctionType;

      typedef Functions::MeshFunction< MeshType, 3, bool > InterfaceMapType;

      using MeshFunctionPointer = SharedPointer< MeshFunctionType >;

      using InterfaceMapPointer = SharedPointer< InterfaceMapType >;      



      void initInterface( const MeshFunctionType& input,

                          MeshFunctionType& output,

                          InterfaceMapType& interfaceMap );

      void initInterface( const MeshFunctionPointer& input,

                          MeshFunctionPointer& output,

                          InterfaceMapPointer& interfaceMap );

      

      template< typename MeshEntity >

      void updateCell( MeshFunctionType& u,

                       const MeshEntity& cell );      

      __cuda_callable__ void updateCell( MeshFunctionType& u,

                                         const MeshEntity& cell,

                                         const RealType velocity = 1.0);

      

      __cuda_callable__ bool updateCell( volatile Real sArray[10][10][10],

                                         int thri, int thrj, int thrk, const Real hx, const Real hy, const Real hz,

                                         const Real velocity = 1.0 );

};



template < typename T1, typename T2 >

T1 meet2DCondition( T1 a, T1 b, const T2 ha, const T2 hb, const T1 value, double v = 1);



template < typename T1 >

__cuda_callable__ void sortMinims( T1 pom[] );





#ifdef HAVE_CUDA

template < typename Real, typename Device, typename Index >

__global__ void CudaInitCaller( const Functions::MeshFunction< Meshes::Grid< 1, Real, Device, Index > >& input, 

                                Functions::MeshFunction< Meshes::Grid< 1, Real, Device, Index > >& output,

                                Functions::MeshFunction< Meshes::Grid< 1, Real, Device, Index >, 1, bool >& interfaceMap  );



template < typename Real, typename Device, typename Index >

__global__ void CudaUpdateCellCaller( tnlDirectEikonalMethodsBase< Meshes::Grid< 1, Real, Device, Index > > ptr,

                                      const Functions::MeshFunction< Meshes::Grid< 1, Real, Device, Index >, 1, bool >& interfaceMap,

                                      Functions::MeshFunction< Meshes::Grid< 1, Real, Device, Index > >& aux,

                                      bool *BlockIterDevice );



template < typename Real, typename Device, typename Index >

__global__ void CudaUpdateCellCaller( tnlDirectEikonalMethodsBase< Meshes::Grid< 2, Real, Device, Index > > ptr,

                                      const Functions::MeshFunction< Meshes::Grid< 2, Real, Device, Index >, 2, bool >& interfaceMap,

                                      Real *aux,

                                      int *BlockIterDevice);

__global__ void CudaParallelReduc( int *BlockIterDevice, int *dBlock, int nBlocks );



template < typename Real, typename Device, typename Index >

__global__ void aux1( Functions::MeshFunction< Meshes::Grid< 2, Real, Device, Index > >& aux, Real *dAux, int a );



template < typename Real, typename Device, typename Index >

__global__ void CudaInitCaller( const Functions::MeshFunction< Meshes::Grid< 2, Real, Device, Index > >& input, 

                                Functions::MeshFunction< Meshes::Grid< 2, Real, Device, Index > >& output,

                                Functions::MeshFunction< Meshes::Grid< 2, Real, Device, Index >, 2, bool >& interfaceMap );



template < typename Real, typename Device, typename Index >

__global__ void CudaInitCaller3d( const Functions::MeshFunction< Meshes::Grid< 3, Real, Device, Index > >& input, 

                                  Functions::MeshFunction< Meshes::Grid< 3, Real, Device, Index > >& output,

                                  Functions::MeshFunction< Meshes::Grid< 3, Real, Device, Index >, 3, bool >& interfaceMap );



template < typename Real, typename Device, typename Index >

__global__ void CudaUpdateCellCaller( tnlDirectEikonalMethodsBase< Meshes::Grid< 3, Real, Device, Index > > ptr,

                                      const Functions::MeshFunction< Meshes::Grid< 3, Real, Device, Index >, 3, bool >& interfaceMap,

                                      Functions::MeshFunction< Meshes::Grid< 3, Real, Device, Index > >& aux,

                                      int *BlockIterDevice );

#endif



#include "tnlDirectEikonalMethodsBase_impl.h"