Commit 50cab7c3 authored by Tomas Sobotik's avatar Tomas Sobotik
Browse files

Godunov 3D parallel iterative solver working properly

parent 031804b1

examples/hamilton-jacobi-parallel/tnlParallelEikonalSolver3D_impl.h

+6 −6
Original line number Diff line number Diff line
@@ -949,37 +949,37 @@ void tnlParallelEikonalSolver<3,SchemeHost, SchemeDevice, Device, double, int>:: 
			{

				Ent.setCoordinates(tnlStaticVector<3,int>(0,j,k));

			   	Ent.refresh();

				u[l] = u[Ent.getIndex()] + Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 1, 0, 0 >()*(threadIdx.x) ;//+  2*Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 1, 0, 0 >()*(threadIdx.x+this->n);

				u[l] = u[Ent.getIndex()];// + Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 1, 0, 0 >()*(threadIdx.x) ;//+  2*Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 1, 0, 0 >()*(threadIdx.x+this->n);

			}

			else if(boundaryCondition == 2)

			{

				Ent.setCoordinates(tnlStaticVector<3,int>(blockDim.x - 1,j,k));

			   	Ent.refresh();

				u[l] = u[Ent.getIndex()] + Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 1, 0, 0 >()*(this->n - 1 - threadIdx.x);//+ 2*Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 1, 0, 0 >()*(blockDim.x - threadIdx.x - 1+this->n);

				u[l] = u[Ent.getIndex()];// + Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 1, 0, 0 >()*(this->n - 1 - threadIdx.x);//+ 2*Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 1, 0, 0 >()*(blockDim.x - threadIdx.x - 1+this->n);

			}

			else if(boundaryCondition == 8)

			{

				Ent.setCoordinates(tnlStaticVector<3,int>(i,0,k));

			   	Ent.refresh();

				u[l] = u[Ent.getIndex()] + Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 0, 1, 0 >()*(threadIdx.y) ;//+ 2*Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 1, 0, 0 >()*(threadIdx.y+this->n);

				u[l] = u[Ent.getIndex()];// + Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 0, 1, 0 >()*(threadIdx.y) ;//+ 2*Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 1, 0, 0 >()*(threadIdx.y+this->n);

			}

			else if(boundaryCondition == 1)

			{

				Ent.setCoordinates(tnlStaticVector<3,int>(i,blockDim.y - 1,k));

			   	Ent.refresh();

				u[l] = u[Ent.getIndex()] + Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 0, 1, 0 >()*(this->n - 1 - threadIdx.y) ;//+ Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 1, 0, 0 >()*(blockDim.y - threadIdx.y  - 1 +this->n);

				u[l] = u[Ent.getIndex()];// + Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 0, 1, 0 >()*(this->n - 1 - threadIdx.y) ;//+ Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 1, 0, 0 >()*(blockDim.y - threadIdx.y  - 1 +this->n);

			}

			else if(boundaryCondition == 32)

			{

				Ent.setCoordinates(tnlStaticVector<3,int>(i,j,0));

			   	Ent.refresh();

				u[l] = u[Ent.getIndex()] + Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 0, 0, 1 >()*(threadIdx.z);

				u[l] = u[Ent.getIndex()];// + Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 0, 0, 1 >()*(threadIdx.z);

			}

			else if(boundaryCondition == 16)

			{

				Ent.setCoordinates(tnlStaticVector<3,int>(i,j,blockDim.z - 1));

			   	Ent.refresh();

				u[l] = u[Ent.getIndex()] + Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 0, 0, 1 >()*(this->n - 1 - threadIdx.z) ;

				u[l] = u[Ent.getIndex()];// + Sign(u[0])*this->subMesh.template getSpaceStepsProducts< 0, 0, 1 >()*(this->n - 1 - threadIdx.z) ;

			}

		}

	}

src/operators/godunov-eikonal/parallelGodunovEikonal3D_impl.h

+22 −22
Original line number Diff line number Diff line
@@ -39,7 +39,7 @@ template< typename MeshReal, 
Real  parallelGodunovEikonalScheme< tnlGrid< 3,MeshReal, Device, MeshIndex >, Real, Index > :: negativePart(const Real arg) const

{

	if(arg < 0.0)

		return arg;

		return -arg;

	return 0.0;

}
@@ -272,10 +272,10 @@ Real parallelGodunovEikonalScheme< tnlGrid< 3, MeshReal, Device, MeshIndex >, Re 
          	  	  	  	  	  	  	  	  	  	  	  	  	  	  	  	  	  	  	  	 ) const

{

	RealType signui;

	if(boundaryCondition == 0)

//	if(boundaryCondition == 0)

		signui = sign(u[cellIndex], this->epsilon);

	else

		signui = Sign(u[cellIndex]);

//	else

//		signui = Sign(u[cellIndex]);





	RealType xb = u[cellIndex];
@@ -358,24 +358,24 @@ Real parallelGodunovEikonalScheme< tnlGrid< 3, MeshReal, Device, MeshIndex >, Re 
	   }





//	   if(xb - xf > 0.0)

//		   a = xb;

//	   else

//		   a = xf;

//

//	   if(yb - yf > 0.0)

//		   b = yb;

//	   else

//		   b = yf;

//

//	   if(zb - zf > 0.0)

//		   c = zb;

//	   else

//		   c = zf;

//

//	   d = ( 1.0 - sqrt(a*a + b*b + c*c)*ihx );

	   if(xb - xf > 0.0)

		   a = xb;

	   else

		   a = xf;



	   if(yb - yf > 0.0)

		   b = yb;

	   else

		   b = yf;



	   if(zb - zf > 0.0)

		   c = zb;

	   else

		   c = zf;



	   d = ( 1.0 - sqrt(a*a + b*b + c*c)*ihx );



	   d = 1.0 - sqrt(xf*xf + xb*xb + yf*yf + yb*yb + zf*zf + zb*zb)*ihx; /*upwind*/

//	   d = 1.0 - sqrt(xf*xf + xb*xb + yf*yf + yb*yb + zf*zf + zb*zb)*ihx; /*upwind*/



	   if(Sign(d) > 0.0 )

		   return Sign(u[cellIndex])*d;