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

Merge branch 'hamilton-jacobi' into 'develop' 

Hamilton jacobi

See merge request !27
parents 3f05fea7 7ce90375

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

+0 −3
Original line number Diff line number Diff line
@@ -127,9 +127,6 @@ class tnlDirectEikonalMethodsBase< Meshes::Grid< 3, Real, Device, Index > > 
            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[] );

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

+53 −73
Original line number Diff line number Diff line
@@ -102,7 +102,7 @@ updateBlocks( InterfaceMapType interfaceMap, 
        ArrayContainer BlockIterHost, int numThreadsPerBlock/*, Real **sArray*/ )

{

#pragma omp parallel for schedule( dynamic )

  for( int i = 0; i < BlockIterHost.getSize(); i++ )

  for( IndexType i = 0; i < BlockIterHost.getSize(); i++ )

  {

    if( BlockIterHost[ i ] )

    {
@@ -142,14 +142,14 @@ updateBlocks( InterfaceMapType interfaceMap, 
      if( sArray == nullptr )

        std::cout << "Error while allocating memory for sArray." << std::endl;

      

      for( int thri = 0; thri < sizeSArray; thri++ ){

        for( int thrj = 0; thrj < sizeSArray; thrj++ )

      for( IndexType thri = 0; thri < sizeSArray; thri++ ){

        for( IndexType thrj = 0; thrj < sizeSArray; thrj++ )

          sArray[ thri * sizeSArray + thrj ] = std::numeric_limits< Real >::max();

      }

      

      

      //printf("numThreadsPerBlock = %d\n", numThreadsPerBlock);

      for( int thrj = 0; thrj < numThreadsPerBlock + 1; thrj++ )

      for( IndexType thrj = 0; thrj < numThreadsPerBlock + 1; thrj++ )

      {        

        if( dimX > (blIdx+1) * numThreadsPerBlock  && thrj+1 < ykolik )

          sArray[ ( thrj+1 )* sizeSArray +xkolik] = aux[ blIdy*numThreadsPerBlock*dimX - dimX + blIdx*numThreadsPerBlock - 1 + (thrj+1)*dimX + xkolik ];
@@ -165,14 +165,14 @@ updateBlocks( InterfaceMapType interfaceMap, 
          sArray[thrj+1] = aux[ blIdy*numThreadsPerBlock*dimX - dimX + blIdx*numThreadsPerBlock - 1 + thrj+1 ];

      }

      

      for( int k = 0; k < numThreadsPerBlock; k++ ){

        for( int l = 0; l < numThreadsPerBlock; l++ )

      for( IndexType k = 0; k < numThreadsPerBlock; k++ ){

        for( IndexType l = 0; l < numThreadsPerBlock; l++ )

          if( blIdy * numThreadsPerBlock + k < dimY && blIdx * numThreadsPerBlock + l < dimX )

            sArray[(k+1) * sizeSArray + l+1] = aux[ blIdy * numThreadsPerBlock * dimX + numThreadsPerBlock * blIdx  + k*dimX + l ];

      }

      

      for( int k = 0; k < numThreadsPerBlock; k++ ){ 

        for( int l = 0; l < numThreadsPerBlock; l++ ){

      for( IndexType k = 0; k < numThreadsPerBlock; k++ ){ 

        for( IndexType l = 0; l < numThreadsPerBlock; l++ ){

          if( blIdy * numThreadsPerBlock + k < dimY && blIdx * numThreadsPerBlock + l < dimX ){

            //std::cout << "proslo i = " << k * numThreadsPerBlock + l << std::endl;

            if( ! interfaceMap[ blIdy * numThreadsPerBlock * dimX + numThreadsPerBlock * blIdx  + k*dimX + l ] )
@@ -191,8 +191,8 @@ updateBlocks( InterfaceMapType interfaceMap, 
       std::cout << std::endl;

       }*/

      

      for( int k = 0; k < numThreadsPerBlock; k++ ) 

        for( int l = numThreadsPerBlock-1; l >-1; l-- ) { 

      for( IndexType k = 0; k < numThreadsPerBlock; k++ ) 

        for( IndexType l = numThreadsPerBlock-1; l >-1; l-- ) { 

          if( blIdy * numThreadsPerBlock + k < dimY && blIdx * numThreadsPerBlock + l < dimX )

          {

            if( ! interfaceMap[ blIdy * numThreadsPerBlock * dimX + numThreadsPerBlock * blIdx  + k*dimX + l ] )
@@ -209,8 +209,8 @@ updateBlocks( InterfaceMapType interfaceMap, 
       std::cout << std::endl;

       }*/

      

      for( int k = numThreadsPerBlock-1; k > -1; k-- ) 

        for( int l = 0; l < numThreadsPerBlock; l++ ) {

      for( IndexType k = numThreadsPerBlock-1; k > -1; k-- ) 

        for( IndexType l = 0; l < numThreadsPerBlock; l++ ) {

          if( blIdy * numThreadsPerBlock + k < dimY && blIdx * numThreadsPerBlock + l < dimX )

          {

            if( ! interfaceMap[ blIdy * numThreadsPerBlock * dimX + numThreadsPerBlock * blIdx  + k*dimX + l ] )
@@ -227,8 +227,8 @@ updateBlocks( InterfaceMapType interfaceMap, 
       std::cout << std::endl;

       }*/

      

      for( int k = numThreadsPerBlock-1; k > -1; k-- ){

        for( int l = numThreadsPerBlock-1; l >-1; l-- ) { 

      for( IndexType k = numThreadsPerBlock-1; k > -1; k-- ){

        for( IndexType l = numThreadsPerBlock-1; l >-1; l-- ) { 

          if( blIdy * numThreadsPerBlock + k < dimY && blIdx * numThreadsPerBlock + l < dimX )

          {

            if( ! interfaceMap[ blIdy * numThreadsPerBlock * dimX + numThreadsPerBlock * blIdx  + k*dimX + l ] )
@@ -252,8 +252,8 @@ updateBlocks( InterfaceMapType interfaceMap, 
      }

      

      

      for( int k = 0; k < numThreadsPerBlock; k++ ){ 

        for( int l = 0; l < numThreadsPerBlock; l++ ) {

      for( IndexType k = 0; k < numThreadsPerBlock; k++ ){ 

        for( IndexType l = 0; l < numThreadsPerBlock; l++ ) {

          if( blIdy * numThreadsPerBlock + k < dimY && blIdx * numThreadsPerBlock + l < dimX )      

            helpFunc[ blIdy * numThreadsPerBlock * dimX + numThreadsPerBlock * blIdx  + k*dimX + l ] = sArray[ (k + 1)* sizeSArray + l + 1 ];

          //std::cout<< sArray[k+1][l+1];
@@ -276,7 +276,7 @@ updateBlocks( const InterfaceMapType interfaceMap, 
        ArrayContainer BlockIterHost, int numThreadsPerBlock/*, Real **sArray*/ )

{  

//#pragma omp parallel for schedule( dynamic )

  for( int i = 0; i < BlockIterHost.getSize(); i++ )

  for( IndexType i = 0; i < BlockIterHost.getSize(); i++ )

  {

    if( BlockIterHost[ i ] )

    {
@@ -323,15 +323,15 @@ updateBlocks( const InterfaceMapType interfaceMap, 
      if( sArray == nullptr )

        std::cout << "Error while allocating memory for sArray." << std::endl;

      

      for( int k = 0; k < sizeSArray; k++ )

        for( int l = 0; l < sizeSArray; l++ )

          for( int m = 0; m < sizeSArray; m++ ){

      for( IndexType k = 0; k < sizeSArray; k++ )

        for( IndexType l = 0; l < sizeSArray; l++ )

          for( IndexType m = 0; m < sizeSArray; m++ ){

            sArray[ m * sizeSArray * sizeSArray + k * sizeSArray + l ] = std::numeric_limits< Real >::max();

          }

      

      

      for( int thrk = 0; thrk < numThreadsPerBlock; thrk++ )

        for( int thrj = 0; thrj < numThreadsPerBlock; thrj++ )

      for( IndexType thrk = 0; thrk < numThreadsPerBlock; thrk++ )

        for( IndexType thrj = 0; thrj < numThreadsPerBlock; thrj++ )

        {

          if( blIdx != 0 && thrj+1 < ykolik && thrk+1 < zkolik )

            sArray[(thrk+1 )* sizeSArray * sizeSArray + (thrj+1)*sizeSArray + 0] =
@@ -358,9 +358,9 @@ updateBlocks( const InterfaceMapType interfaceMap, 
                    aux[ (blIdz+1)*numThreadsPerBlock * dimX * dimY + blIdy * numThreadsPerBlock*dimX + blIdx*numThreadsPerBlock + thrj * dimX + thrk ];

        }

      

      for( int m = 0; m < numThreadsPerBlock; m++ ){

        for( int k = 0; k < numThreadsPerBlock; k++ ){

          for( int l = 0; l < numThreadsPerBlock; l++ ){

      for( IndexType m = 0; m < numThreadsPerBlock; m++ ){

        for( IndexType k = 0; k < numThreadsPerBlock; k++ ){

          for( IndexType l = 0; l < numThreadsPerBlock; l++ ){

            if( blIdy * numThreadsPerBlock + k < dimY && blIdx * numThreadsPerBlock + l < dimX && blIdz * numThreadsPerBlock + m < dimZ )

              sArray[(m+1) * sizeSArray * sizeSArray + (k+1) *sizeSArray + l+1] = 

                      aux[ blIdz * numThreadsPerBlock * dimX * dimY + blIdy * numThreadsPerBlock * dimX + blIdx*numThreadsPerBlock + m*dimX*dimY + k*dimX + l ];
@@ -379,9 +379,9 @@ updateBlocks( const InterfaceMapType interfaceMap, 
      s = "helpFunc-"+ std::to_string(numWhile) + ".tnl";

      helpFunc.save( s );*/

      

      for( int m = 0; m < numThreadsPerBlock; m++ ){

        for( int k = 0; k < numThreadsPerBlock; k++ ){ 

          for( int l = 0; l < numThreadsPerBlock; l++ ){

      for( IndexType m = 0; m < numThreadsPerBlock; m++ ){

        for( IndexType k = 0; k < numThreadsPerBlock; k++ ){ 

          for( IndexType l = 0; l < numThreadsPerBlock; l++ ){

            if( blIdy * numThreadsPerBlock + k < dimY && blIdx * numThreadsPerBlock + l < dimX && blIdz * numThreadsPerBlock + m < dimZ ){

              //std::cout << "proslo i = " << k * numThreadsPerBlock + l << std::endl;

              if( ! interfaceMap[ blIdz * numThreadsPerBlock * dimX * dimY + blIdy * numThreadsPerBlock * dimX + blIdx*numThreadsPerBlock + m*dimX*dimY + k*dimX + l  ] )
@@ -404,9 +404,9 @@ updateBlocks( const InterfaceMapType interfaceMap, 
      s = "helpFunc-"+ std::to_string(numWhile) + ".tnl";

      helpFunc.save( s );*/

      

      for( int m = numThreadsPerBlock-1; m >-1; m-- ){

        for( int k = 0; k < numThreadsPerBlock; k++ ){

          for( int l = 0; l <numThreadsPerBlock; l++ ){

      for( IndexType m = numThreadsPerBlock-1; m >-1; m-- ){

        for( IndexType k = 0; k < numThreadsPerBlock; k++ ){

          for( IndexType l = 0; l <numThreadsPerBlock; l++ ){

            if( blIdy * numThreadsPerBlock + k < dimY && blIdx * numThreadsPerBlock + l < dimX && blIdz * numThreadsPerBlock + m < dimZ )

            {

              if( ! interfaceMap[ blIdz * numThreadsPerBlock * dimX * dimY + blIdy * numThreadsPerBlock * dimX + blIdx*numThreadsPerBlock + m*dimX*dimY + k*dimX + l  ] )
@@ -427,9 +427,9 @@ updateBlocks( const InterfaceMapType interfaceMap, 
      s = "helpFunc-"+ std::to_string(numWhile) + ".tnl";

      helpFunc.save( s );*/

      

      for( int m = 0; m < numThreadsPerBlock; m++ ){

        for( int k = 0; k < numThreadsPerBlock; k++ ){

          for( int l = numThreadsPerBlock-1; l >-1; l-- ){

      for( IndexType m = 0; m < numThreadsPerBlock; m++ ){

        for( IndexType k = 0; k < numThreadsPerBlock; k++ ){

          for( IndexType l = numThreadsPerBlock-1; l >-1; l-- ){

            if( blIdy * numThreadsPerBlock + k < dimY && blIdx * numThreadsPerBlock + l < dimX && blIdz * numThreadsPerBlock + m < dimZ )

            {

              if( ! interfaceMap[ blIdz * numThreadsPerBlock * dimX * dimY + blIdy * numThreadsPerBlock * dimX + blIdx*numThreadsPerBlock + m*dimX*dimY + k*dimX + l  ] )
@@ -450,9 +450,9 @@ updateBlocks( const InterfaceMapType interfaceMap, 
      s = "helpFunc-"+ std::to_string(numWhile) + ".tnl";

      helpFunc.save( s );

      */

      for( int m = numThreadsPerBlock-1; m >-1; m-- ){

        for( int k = 0; k < numThreadsPerBlock; k++ ){

          for( int l = numThreadsPerBlock-1; l >-1; l-- ){

      for( IndexType m = numThreadsPerBlock-1; m >-1; m-- ){

        for( IndexType k = 0; k < numThreadsPerBlock; k++ ){

          for( IndexType l = numThreadsPerBlock-1; l >-1; l-- ){

            if( blIdy * numThreadsPerBlock + k < dimY && blIdx * numThreadsPerBlock + l < dimX && blIdz * numThreadsPerBlock + m < dimZ )

            {

              if( ! interfaceMap[ blIdz * numThreadsPerBlock * dimX * dimY + blIdy * numThreadsPerBlock * dimX + blIdx*numThreadsPerBlock + m*dimX*dimY + k*dimX + l  ] )
@@ -473,9 +473,9 @@ updateBlocks( const InterfaceMapType interfaceMap, 
      s = "helpFunc-"+ std::to_string(numWhile) + ".tnl";

      helpFunc.save( s );*/

      

      for( int m = 0; m < numThreadsPerBlock; m++ ){

        for( int k = numThreadsPerBlock-1; k > -1; k-- ){

          for( int l = 0; l <numThreadsPerBlock; l++ ){

      for( IndexType m = 0; m < numThreadsPerBlock; m++ ){

        for( IndexType k = numThreadsPerBlock-1; k > -1; k-- ){

          for( IndexType l = 0; l <numThreadsPerBlock; l++ ){

            if( blIdy * numThreadsPerBlock + k < dimY && blIdx * numThreadsPerBlock + l < dimX && blIdz * numThreadsPerBlock + m < dimZ )

            {

              if( ! interfaceMap[blIdz * numThreadsPerBlock * dimX * dimY + blIdy * numThreadsPerBlock * dimX + blIdx*numThreadsPerBlock + m*dimX*dimY + k*dimX + l  ] )
@@ -496,9 +496,9 @@ updateBlocks( const InterfaceMapType interfaceMap, 
      s = "helpFunc-"+ std::to_string(numWhile) + ".tnl";

      helpFunc.save( s );*/

      

      for( int m = numThreadsPerBlock-1; m >-1; m-- ){

        for( int k = numThreadsPerBlock-1; k > -1; k-- ){

          for( int l = 0; l <numThreadsPerBlock; l++ ){

      for( IndexType m = numThreadsPerBlock-1; m >-1; m-- ){

        for( IndexType k = numThreadsPerBlock-1; k > -1; k-- ){

          for( IndexType l = 0; l <numThreadsPerBlock; l++ ){

            if( blIdy * numThreadsPerBlock + k < dimY && blIdx * numThreadsPerBlock + l < dimX && blIdz * numThreadsPerBlock + m < dimZ )

            {

              if( ! interfaceMap[ blIdz * numThreadsPerBlock * dimX * dimY + blIdy * numThreadsPerBlock * dimX + blIdx*numThreadsPerBlock + m*dimX*dimY + k*dimX + l  ] )
@@ -519,9 +519,9 @@ updateBlocks( const InterfaceMapType interfaceMap, 
      s = "helpFunc-"+ std::to_string(numWhile) + ".tnl";

      helpFunc.save( s );*/

      

      for( int m = 0; m < numThreadsPerBlock; m++ ){

        for( int k = numThreadsPerBlock-1; k > -1; k-- ){

          for( int l = numThreadsPerBlock-1; l >-1; l-- ){

      for( IndexType m = 0; m < numThreadsPerBlock; m++ ){

        for( IndexType k = numThreadsPerBlock-1; k > -1; k-- ){

          for( IndexType l = numThreadsPerBlock-1; l >-1; l-- ){

            if( blIdy * numThreadsPerBlock + k < dimY && blIdx * numThreadsPerBlock + l < dimX && blIdz * numThreadsPerBlock + m < dimZ )

            {

              if( ! interfaceMap[ blIdz * numThreadsPerBlock * dimX * dimY + blIdy * numThreadsPerBlock * dimX + blIdx*numThreadsPerBlock + m*dimX*dimY + k*dimX + l  ] )
@@ -543,9 +543,9 @@ updateBlocks( const InterfaceMapType interfaceMap, 
      helpFunc.save( s );*/

      

      

      for( int m = numThreadsPerBlock-1; m >-1; m-- ){

        for( int k = numThreadsPerBlock-1; k > -1; k-- ){

          for( int l = numThreadsPerBlock-1; l >-1; l-- ){

      for( IndexType m = numThreadsPerBlock-1; m >-1; m-- ){

        for( IndexType k = numThreadsPerBlock-1; k > -1; k-- ){

          for( IndexType l = numThreadsPerBlock-1; l >-1; l-- ){

            if( blIdy * numThreadsPerBlock + k < dimY && blIdx * numThreadsPerBlock + l < dimX && blIdz * numThreadsPerBlock + m < dimZ )

            {

              if( ! interfaceMap[ blIdz * numThreadsPerBlock * dimX * dimY + blIdy * numThreadsPerBlock * dimX + blIdx*numThreadsPerBlock + m*dimX*dimY + k*dimX + l  ] )
@@ -571,9 +571,9 @@ updateBlocks( const InterfaceMapType interfaceMap, 
      }

      

      

      for( int k = 0; k < numThreadsPerBlock; k++ ){ 

        for( int l = 0; l < numThreadsPerBlock; l++ ) {

          for( int m = 0; m < numThreadsPerBlock; m++ ){

      for( IndexType k = 0; k < numThreadsPerBlock; k++ ){ 

        for( IndexType l = 0; l < numThreadsPerBlock; l++ ) {

          for( IndexType m = 0; m < numThreadsPerBlock; m++ ){

            if( blIdy * numThreadsPerBlock + k < dimY && blIdx * numThreadsPerBlock + l < dimX && blIdz * numThreadsPerBlock + m < dimZ ){      

              helpFunc[ blIdz * numThreadsPerBlock * dimX * dimY + blIdy * numThreadsPerBlock * dimX + blIdx*numThreadsPerBlock + m*dimX*dimY + k*dimX + l ] = 

                      sArray[ (m+1) * sizeSArray * sizeSArray + (k+1) *sizeSArray + l+1 ];
@@ -1224,26 +1224,6 @@ updateCell( MeshFunctionType& u, 
  }

}



template < typename T1, typename T2 >

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

{

  T1 tmp;

  if( fabs( a ) != std::numeric_limits< T1 >::max &&

          fabs( b ) != std::numeric_limits< T1 >::max &&

          fabs( a - b ) < ha/v )//TNL::sqrt( (ha * ha + hb * hb)/2 )/v )

  {

    tmp = ( ha * ha * b + hb * hb * a + 

            TNL::sign( value ) * ha * hb * TNL::sqrt( ( ha * ha + hb * hb )/( v * v ) - 

            ( a - b ) * ( a - b ) ) )/( ha * ha + hb * hb );

  }

  else

  {

    tmp = std::numeric_limits< T1 >::max;

  }

  

  return tmp;

}



template < typename T1 >

__cuda_callable__ void sortMinims( T1 pom[] )

{
@@ -1274,7 +1254,7 @@ __cuda_callable__ void sortMinims( T1 pom[] ) 
    tmp[3] = pom[5]; tmp[4] = pom[4]; tmp[5] = pom[3];

  }

  

  for( int i = 0; i < 6; i++ )

  for( unsigned int i = 0; i < 6; i++ )

  {

    pom[ i ] = tmp[ i ];

  }

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

+7 −1
Original line number Diff line number Diff line
@@ -18,6 +18,8 @@ 
#include <TNL/Pointers/SharedPointer.h>

#include "tnlFastSweepingMethod.h"



#include <TNL/Meshes/DistributedMeshes/DistributedGridIO.h>



template< typename Mesh,

          typename Communicator,

          typename Anisotropy,
@@ -46,6 +48,8 @@ class tnlDirectEikonalProblem 
      using MeshPointer = Pointers::SharedPointer< MeshType >;

      using DofVectorPointer = Pointers::SharedPointer< DofVectorType >;

      

      typedef Communicator CommunicatorType;

      

      static constexpr bool isTimeDependent() { return false; };



      static String getType();
@@ -79,6 +83,8 @@ class tnlDirectEikonalProblem 
         

         AnisotropyPointer anisotropy;



         Meshes::DistributedMeshes::DistrGridIOTypes distributedIOType;

         

   };



#include "tnlDirectEikonalProblem_impl.h"

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

+20 −5
Original line number Diff line number Diff line
@@ -113,10 +113,25 @@ tnlDirectEikonalProblem< Mesh, Communicator, Anisotropy, Real, Index >:: 
setInitialCondition( const Config::ParameterContainer& parameters,

                     DofVectorPointer& dofs )

{

  this->bindDofs( dofs );

  String inputFile = parameters.getParameter< String >( "input-file" );

  this->initialData->setMesh( this->getMesh() );

  std::cout<<"setInitialCondition" <<std::endl; 

  if(CommunicatorType::isDistributed())

  {

    std::cout<<"Nodes Distribution: " << u->getMesh().getDistributedMesh()->printProcessDistr() << std::endl;

    if(distributedIOType==Meshes::DistributedMeshes::MpiIO)

      Meshes::DistributedMeshes::DistributedGridIO<MeshFunctionType,Meshes::DistributedMeshes::MpiIO> ::load(inputFile, *u );

    if(distributedIOType==Meshes::DistributedMeshes::LocalCopy)

      Meshes::DistributedMeshes::DistributedGridIO<MeshFunctionType,Meshes::DistributedMeshes::LocalCopy> ::load(inputFile, *u );

    u->template synchronize<CommunicatorType>();

  }

  else

  {

    if( !this->initialData->boundLoad( inputFile ) )

      std::cerr << "I am not able to load the initial condition from the file " << inputFile << "." << std::endl;

    return false;

  }

  return true;

}

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

+25 −4
Original line number Diff line number Diff line
@@ -83,13 +83,34 @@ solve( const MeshPointer& mesh, 
  MeshFunctionType aux = *auxPtr;

  

  

//#ifdef HAVE_MPI

  bool a = Communicators::MpiCommunicator::IsInitialized();

#ifdef HAVE_MPI

  int i = Communicators::MpiCommunicator::GetRank( Communicators::MpiCommunicator::AllGroup );

  //printf( "Hello world from rank: %d ", i );

  //Communicators::MpiCommunicator::Request r = Communicators::MpiCommunicator::ISend( auxPtr, 0, 0, Communicators::MpiCommunicator::AllGroup );

  if( i == 1 )

    /*for( int k = 0; k < 16*16; k++ )

      aux[ k ] = 10;*/

    printf( "1: mesh x: %d\n", mesh->getDimensions().x() );

    printf( "1: mesh y: %d\n", mesh->getDimensions().y() );

    //aux.save("aux_proc1.tnl");

  if( i == 0 )

    printf( "0: mesh x: %d\n", mesh->getDimensions().x() );

    printf( "0: mesh y: %d\n", mesh->getDimensions().y() );

    //aux.save("aux_proc0.tnl");

    /*for( int k = 0; k < mesh->getDimensions().x()*mesh->getDimensions().y(); k++ )

      aux[ k ] = 10;

    for( int k = 0; k < mesh->getDimensions().x(); k++ ){

      for( int l = 0; l < mesh->getDimensions().y(); l++ )

        printf("%f.2\t",aux[ k * 16 + l ] );

    printf("\n");

    }*/

    

  /*bool a = Communicators::MpiCommunicator::IsInitialized();

  if( a )

    printf("Je Init\n");

  else

    printf("Neni Init\n");

//#endif

    printf("Neni Init\n");*/

#endif

  

  while( iteration < this->maxIterations )

  {

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

+8 −8

File changed.

Contains only whitespace changes.