Commit 337fb813 authored by Tomas Sobotik's avatar Tomas Sobotik
Browse files

3D parallel iterative solver almost fully fixed

parent 67e1ee69
Loading
Loading
Loading
Loading
+1 −1
Original line number Diff line number Diff line
@@ -29,7 +29,7 @@ class fastSweepingConfig
         config.addDelimiter( "Parallel Eikonal solver settings:" );
         config.addEntry        < tnlString > ( "problem-name", "This defines particular problem.", "fast-sweeping" );
         config.addRequiredEntry        < tnlString > ( "initial-condition", "Initial condition for solver");
         config.addRequiredEntry        < tnlString > ( "dim", "Dimension of problem.");
         config.addRequiredEntry        < int > ( "dim", "Dimension of problem.");
         config.addEntry       < tnlString > ( "mesh", "Name of mesh.", "mesh.tnl" );
         config.addEntry       < tnlString > ( "exact-input", "Are the function values near the curve equal to the SDF? (yes/no)", "no" );
      }
+3 −3
Original line number Diff line number Diff line
@@ -43,9 +43,9 @@ int main( int argc, char* argv[] )
   if( ! parseCommandLine( argc, argv, configDescription, parameters ) )
      return false;

   const tnlString& dim = parameters.getParameter< tnlString >( "dim" );
   const int& dim = parameters.getParameter< int >( "dim" );

   if(dim == "2")
   if(dim == 2)
   {
		tnlFastSweeping<tnlGrid<2,double,tnlHost, int>, double, int> solver;
		if(!solver.init(parameters))
@@ -58,7 +58,7 @@ int main( int argc, char* argv[] )
	   cout << "Starting solver..." << endl;
	   solver.run();
   }
   else if(dim == "3")
   else if(dim == 3)
   {
		tnlFastSweeping<tnlGrid<3,double,tnlHost, int>, double, int> solver;
		if(!solver.init(parameters))
+87 −38
Original line number Diff line number Diff line
@@ -44,6 +44,54 @@ int main( int argc, char* argv[] )
   tnlDeviceEnum device;
   device = tnlHostDevice;

   const int& dim = parameters.getParameter< int >( "dim" );

  if(dim == 2)
  {

	   typedef parallelGodunovEikonalScheme< tnlGrid<2,double,tnlHost, int>, double, int > SchemeTypeHost;
		/*#ifdef HAVE_CUDA
		   typedef parallelGodunovEikonalScheme< tnlGrid<2,double,tnlCuda, int>, double, int > SchemeTypeDevice;
		#endif
		#ifndef HAVE_CUDA*/
	   typedef parallelGodunovEikonalScheme< tnlGrid<2,double,tnlHost, int>, double, int > SchemeTypeDevice;
		/*#endif*/

	   if(device==tnlHostDevice)
	   {
		   typedef tnlHost Device;


		   tnlParallelEikonalSolver<2,SchemeTypeHost,SchemeTypeDevice, Device> solver;
		   if(!solver.init(parameters))
		   {
			   cerr << "Solver failed to initialize." << endl;
			   return EXIT_FAILURE;
		   }
		   cout << "-------------------------------------------------------------" << endl;
		   cout << "Starting solver loop..." << endl;
		   solver.run();
	   }
	   else if(device==tnlCudaDevice )
	   {
		   typedef tnlCuda Device;
		   //typedef parallelGodunovEikonalScheme< tnlGrid<2,double,Device, int>, double, int > SchemeType;

		   tnlParallelEikonalSolver<2,SchemeTypeHost,SchemeTypeDevice, Device> solver;
		   if(!solver.init(parameters))
		   {
			   cerr << "Solver failed to initialize." << endl;
			   return EXIT_FAILURE;
		   }
		   cout << "-------------------------------------------------------------" << endl;
		   cout << "Starting solver loop..." << endl;
		   solver.run();
	   }
  // }
  }
  else if(dim == 3)
  {

	   typedef parallelGodunovEikonalScheme< tnlGrid<3,double,tnlHost, int>, double, int > SchemeTypeHost;
		/*#ifdef HAVE_CUDA
		   typedef parallelGodunovEikonalScheme< tnlGrid<2,double,tnlCuda, int>, double, int > SchemeTypeDevice;
@@ -83,6 +131,7 @@ int main( int argc, char* argv[] )
		   solver.run();
	   }
 // }
  }

   time(&stop);
   cout << endl;
+1 −0
Original line number Diff line number Diff line
@@ -39,6 +39,7 @@ class parallelEikonalConfig
         config.addRequiredEntry        < double > ( "initial-tau", " initial tau for solver" );
         config.addEntry        < double > ( "cfl-condition", " CFL condition", 0.0 );
         config.addEntry        < int > ( "subgrid-size", "Subgrid size.", 16 );
         config.addRequiredEntry        < int > ( "dim", "Dimension of problem.");
      }
};

+2 −16
Original line number Diff line number Diff line
@@ -21,7 +21,6 @@
#include "tnlParallelEikonalSolver.h"
#include <core/mfilename.h>


template< typename SchemeHost, typename SchemeDevice, typename Device>
tnlParallelEikonalSolver<2,SchemeHost, SchemeDevice, Device, double, int>::tnlParallelEikonalSolver()
{
@@ -1044,7 +1043,7 @@ template< typename SchemeHost, typename SchemeDevice, typename Device>
__device__
void tnlParallelEikonalSolver<2,SchemeHost, SchemeDevice, Device, double, int>::updateSubgridCUDA2D( const int i ,tnlParallelEikonalSolver<2,SchemeHost, SchemeDevice, Device, double, int >* caller, double* a)
{
	int j = threadIdx.x + threadIdx.y * blockDim.x;
//	int j = threadIdx.x + threadIdx.y * blockDim.x;
	int index = (blockIdx.y) * caller->n*caller->n*caller->gridCols
            + (blockIdx.x) * caller->n
            + threadIdx.y * caller->n*caller->gridCols
@@ -1166,7 +1165,6 @@ void tnlParallelEikonalSolver<2,SchemeHost, SchemeDevice, Device, double, int>::
   __syncthreads();
//   if( time + currentTau > finalTime ) currentTau = finalTime - time;


   tnlGridEntity<MeshType, 2, tnlGridEntityNoStencilStorage > Entity(subMesh);
   tnlNeighbourGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighbourEntities(Entity);
   Entity.setCoordinates(tnlStaticVector<2,int>(i,j));
@@ -1312,20 +1310,8 @@ void /*tnlParallelEikonalSolver<2,SchemeHost, SchemeDevice, Device, double, int>
			subgridValue_cmp = cudaSolver->getSubgridValueCUDA2D(blockIdx.y*gridDim.x + blockIdx.x + 1);
			boundary_index = 1;
		}
//		if(threadIdx.y == 0 && (blockIdx.y != 0) && (cudaSolver->currentStep & 1))
//		{
//			u_cmp = cudaSolver->work_u_cuda[gid - blockDim.x*gridDim.x];
//			subgridValue_cmp = cudaSolver->getSubgridValueCUDA2D((blockIdx.y - 1)*gridDim.x + blockIdx.x);
//			boundary_index = 3;
//		}
//		if(threadIdx.y == blockDim.y - 1 && (blockIdx.y != gridDim.y - 1) && (cudaSolver->currentStep & 1))
//		{
//			u_cmp = cudaSolver->work_u_cuda[gid + blockDim.x*gridDim.x];
//			subgridValue_cmp = cudaSolver->getSubgridValueCUDA2D((blockIdx.y + 1)*gridDim.x + blockIdx.x);
//			boundary_index = 0;
//		}

//		__threadfence();
		__threadfence();
		if((subgridValue == INT_MAX || fabs(u_cmp) + cudaSolver->delta < fabs(u) ) && (subgridValue_cmp != INT_MAX && subgridValue_cmp != -INT_MAX))
		{
			cudaSolver->unusedCell_cuda[gid] = 0;
Loading