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

Refactoring mean-curvature flow.

parent e96dbb79
Loading
Loading
Loading
Loading
+1 −1
Original line number Diff line number Diff line
@@ -155,7 +155,7 @@ class tnlNeighbourGridEntityGetter<
              cerr << "entity.getCoordinates()  + CoordinatesType( stepX, stepY ) = " << entity.getCoordinates()  + CoordinatesType( stepX, stepY )
                   << " entity.getGrid().getDimensions() = " << entity.getGrid().getDimensions()
                   << " EntityDimensions = " << EntityDimensions );
            return NeighbourGridEntityType( this->grid,
            return NeighbourGridEntityType( this->entity.getGrid(),
                                            CoordinatesType( entity.getCoordinates().x() + stepX,
                                                             entity.getCoordinates().y() + stepY ) );
      }
+4 −4
Original line number Diff line number Diff line
@@ -153,12 +153,12 @@ class tnlNeighbourGridEntityGetter<
              cerr << "entity.getCoordinates() = " << entity.getCoordinates()
                   << " entity.getGrid().getDimensions() = " << entity.getGrid().getDimensions()
                   << " EntityDimensions = " << EntityDimensions );
         tnlAssert( entity.getCoordinates() + CoordinatesType( stepX, stepY ) >= CoordinatesType( 0, 0, 0 ) &&
                    entity.getCoordinates() + CoordinatesType( stepX, stepY ) < entity.getGrid().getDimensions(),
              cerr << "entity.getCoordinates()  + CoordinatesType( stepX, stepY ) = " << entity.getCoordinates()  + CoordinatesType( stepX, stepY )
         tnlAssert( entity.getCoordinates() + CoordinatesType( stepX, stepY, stepZ ) >= CoordinatesType( 0, 0, 0 ) &&
                    entity.getCoordinates() + CoordinatesType( stepX, stepY, stepZ ) < entity.getGrid().getDimensions(),
              cerr << "entity.getCoordinates()  + CoordinatesType( stepX, stepY ) = " << entity.getCoordinates()  + CoordinatesType( stepX, stepY, stepZ )
                   << " entity.getGrid().getDimensions() = " << entity.getGrid().getDimensions()
                   << " EntityDimensions = " << EntityDimensions );
         return NeighbourGridEntity( CoordinatesType( entity.getCoordinates().x() + stepX,
         return NeighbourGridEntityType( this->entity.getGrid(), CoordinatesType( entity.getCoordinates().x() + stepX,
                                                      entity.getCoordinates().y() + stepY,
                                                      entity.getCoordinates().z() + stepZ ) );
      }
+0 −1
Original line number Diff line number Diff line
@@ -143,7 +143,6 @@ class tnlOneSideDiffNonlinearOperator< tnlGrid< 3, MeshReal, Device, MeshIndex >
             typename MeshEntity >
   __cuda_callable__
   Real getValue( const MeshType& mesh,
                  const IndexType cellIndex,
                  const MeshEntity& entity,
                  const Vector& u,
                  const RealType& time) const;
+1 −1
Original line number Diff line number Diff line
@@ -236,12 +236,12 @@ __cuda_callable__
Real
tnlOneSideDiffNonlinearOperator< tnlGrid< 3, MeshReal, Device, MeshIndex >, OperatorQ, Real, Index >::
getValue( const MeshType& mesh,
          const IndexType cellIndex,
          const MeshEntity& entity,
          const Vector& u,
          const Real& time ) const
{
   const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
   const IndexType& cellIndex = entity.getIndex();
   return operatorQ.getValueStriped( mesh, entity, u, time ) *
      ( ( (  u[ neighbourEntities.template getEntityIndex<  1,  0, 0 >() ] - u[ cellIndex ] ) * mesh.getHxInverse() / operatorQ.getValue( mesh, entity, u, time ) -
          ( -u[ neighbourEntities.template getEntityIndex< -1,  0, 0 >() ] + u[ cellIndex ] ) * mesh.getHxInverse() / operatorQ.getValue( mesh, neighbourEntities.template getEntity< -1,  0,  0 >(), u, time ) ) * mesh.getHxInverse() + 
+34 −33
Original line number Diff line number Diff line
@@ -114,8 +114,9 @@ getValue( const MeshType& mesh,
   const IndexType& cellIndex = entity.getIndex();
   const typename MeshEntity::template NeighbourEntities< 1 >& neighbourEntities = entity.getNeighbourEntities();      
   // TODO: fix this
   return sqrt( this->eps + ( u[ neighbourEntities.template getEntityIndex< 1 >( cellIndex ) ] - u[ cellIndex ]) * 
          ( u[ neighbourEntities.template getEntityIndex< 1 >( cellIndex ) ] - u[ cellIndex ]) *
   return sqrt( this->eps + 
      ( u[ neighbourEntities.template getEntityIndex< 1 >() ] - u[ cellIndex ]) * 
      ( u[ neighbourEntities.template getEntityIndex< 1 >() ] - u[ cellIndex ]) *
          mesh.template getSpaceStepsProducts< -2 >() );
}

@@ -152,10 +153,10 @@ getValueStriped( const MeshType& mesh,
{
   const IndexType& cellIndex = entity.getIndex();
   const typename MeshEntity::template NeighbourEntities< 1 >& neighbourEntities = entity.getNeighbourEntities();      
   return sqrt( this->eps + 0.5*( ( u[ neighbourEntities.template getEntityIndex< 1 >( cellIndex ) ] - u[ cellIndex ]) * 
          ( u[ neighbourEntities.template getEntityIndex< 1 >( cellIndex ) ] - u[ cellIndex ]) *
          mesh.template getSpaceStepsProducts< -1, 0, 0 >() * mesh.template getSpaceStepsProducts< -1, 0, 0 >() + ( - u[ neighbourEntities.template getEntityIndex< -1 >( cellIndex ) ] + u[ cellIndex ] ) 
          * ( - u[ neighbourEntities.template getEntityIndex< -1 >( cellIndex ) ] + u[ cellIndex ] ) * mesh.template getSpaceStepsProducts< -1, 0, 0 >() * mesh.template getSpaceStepsProducts< -1, 0, 0 >() ) );
   return sqrt( this->eps + 0.5*( ( u[ neighbourEntities.template getEntityIndex< 1 >() ] - u[ cellIndex ]) * 
          (  u[ neighbourEntities.template getEntityIndex<  1 >() ] - u[ cellIndex ]) * mesh.template getSpaceStepsProducts< -1 >() * mesh.template getSpaceStepsProducts< -1 >() + 
          ( -u[ neighbourEntities.template getEntityIndex< -1 >() ] + u[ cellIndex ] ) 
          * ( - u[ neighbourEntities.template getEntityIndex< -1 >() ] + u[ cellIndex ] ) * mesh.template getSpaceStepsProducts< -1 >() * mesh.template getSpaceStepsProducts< -1 >() ) );
}

template< typename MeshReal,
@@ -286,10 +287,10 @@ getValue( const MeshType& mesh,
{
   const IndexType& cellIndex = entity.getIndex();
   const typename MeshEntity::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();      
   return sqrt( this->eps + ( u[ neighbourEntities.template getEntityIndex< 0,1 >( cellIndex ) ] - u[ cellIndex ] ) * 
          ( u[ neighbourEntities.template getEntityIndex< 0,1 >( cellIndex ) ] - u[ cellIndex ] )
          * mesh.template getSpaceStepsProducts< 0, -1, 0 >() * mesh.template getSpaceStepsProducts< 0, -1, 0 >() + ( u[ neighbourEntities.template getEntityIndex< 1,0 >( cellIndex ) ] - u[ cellIndex ] ) 
          * ( u[ neighbourEntities.template getEntityIndex< 1,0 >( cellIndex ) ] - u[ cellIndex ] ) * mesh.template getSpaceStepsProducts< -1, 0, 0 >() * mesh.template getSpaceStepsProducts< -1, 0, 0 >() );
   return sqrt( this->eps + ( u[ neighbourEntities.template getEntityIndex< 0,1 >() ] - u[ cellIndex ] ) * 
          ( u[ neighbourEntities.template getEntityIndex< 0, 1 >() ] - u[ cellIndex ] )
          * mesh.template getSpaceStepsProducts< 0, -1 >() * mesh.template getSpaceStepsProducts< 0, -1 >() + ( u[ neighbourEntities.template getEntityIndex< 1,0 >() ] - u[ cellIndex ] ) 
          * ( u[ neighbourEntities.template getEntityIndex< 1, 0 >() ] - u[ cellIndex ] ) * mesh.template getSpaceStepsProducts< -1, 0 >() * mesh.template getSpaceStepsProducts< -1, 0 >() );
}

template< typename MeshReal,
@@ -326,14 +327,14 @@ getValueStriped( const MeshType& mesh,
{
   const IndexType& cellIndex = entity.getIndex();
   const typename MeshEntity::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();      
   return sqrt( this->eps + 0.5*( ( u[ neighbourEntities.template getEntityIndex< 0,1 >( cellIndex ) ] - u[ cellIndex ] ) * 
          ( u[ neighbourEntities.template getEntityIndex< 0,1 >( cellIndex ) ] - u[ cellIndex ] )
          * mesh.template getSpaceStepsProducts< 0, -1, 0 >() * mesh.template getSpaceStepsProducts< 0, -1, 0 >() + ( u[ neighbourEntities.template getEntityIndex< 1,0 >( cellIndex ) ] - u[ cellIndex ] ) 
          * ( u[ neighbourEntities.template getEntityIndex< 1,0 >( cellIndex ) ] - u[ cellIndex ] ) * mesh.template getSpaceStepsProducts< -1, 0, 0 >() * mesh.template getSpaceStepsProducts< -1, 0, 0 >()
          + ( - u[ neighbourEntities.template getEntityIndex< -1,0 >( cellIndex ) ] + u[ cellIndex ]) * ( - u[ neighbourEntities.template getEntityIndex< -1,0 >( cellIndex ) ] + u[ cellIndex ]) 
          * mesh.template getSpaceStepsProducts< -1, 0, 0 >() * mesh.template getSpaceStepsProducts< -1, 0, 0 >()
          + ( - u[ neighbourEntities.template getEntityIndex< 0,-1 >( cellIndex ) ] + u[ cellIndex ]) * ( - u[ neighbourEntities.template getEntityIndex< 0,-1 >( cellIndex ) ] + u[ cellIndex ]) * 
          mesh.template getSpaceStepsProducts< 0, -1, 0 >() * mesh.template getSpaceStepsProducts< 0, -1, 0 >() ) );
   return sqrt( this->eps + 0.5*( ( u[ neighbourEntities.template getEntityIndex< 0,1 >() ] - u[ cellIndex ] ) * 
          ( u[ neighbourEntities.template getEntityIndex< 0, 1 >() ] - u[ cellIndex ] )
          * mesh.template getSpaceStepsProducts< 0, -1 >() * mesh.template getSpaceStepsProducts< 0, -1 >() + ( u[ neighbourEntities.template getEntityIndex< 1,0 >() ] - u[ cellIndex ] ) 
          * ( u[ neighbourEntities.template getEntityIndex< 1, 0 >() ] - u[ cellIndex ] ) * mesh.template getSpaceStepsProducts< -1, 0 >() * mesh.template getSpaceStepsProducts< -1, 0 >()
          + ( - u[ neighbourEntities.template getEntityIndex< -1, 0 >() ] + u[ cellIndex ]) * ( - u[ neighbourEntities.template getEntityIndex< -1,0 >() ] + u[ cellIndex ]) 
          * mesh.template getSpaceStepsProducts< -1, 0 >() * mesh.template getSpaceStepsProducts< -1, 0 >()
          + ( - u[ neighbourEntities.template getEntityIndex< 0,-1 >() ] + u[ cellIndex ]) * ( - u[ neighbourEntities.template getEntityIndex< 0,-1 >() ] + u[ cellIndex ]) * 
          mesh.template getSpaceStepsProducts< 0, -1 >() * mesh.template getSpaceStepsProducts< 0, -1 >() ) );
}

template< typename MeshReal,
@@ -462,11 +463,11 @@ getValue( const MeshType& mesh,
{
   const IndexType& cellIndex = entity.getIndex();
   const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();      
   return sqrt( 1.0 + ( u[ neighbourEntities.template getEntityIndex< 0,1,0 >( cellIndex ) ] - u[ cellIndex ] ) * 
          ( u[ neighbourEntities.template getEntityIndex< 0,1,0 >( cellIndex ) ] - u[ cellIndex ] )
          * mesh.template getSpaceStepsProducts< 0, -1, 0 >() * mesh.template getSpaceStepsProducts< 0, -1, 0 >() + ( u[ neighbourEntities.template getEntityIndex< 1,0,0 >( cellIndex ) ] - u[ cellIndex ] ) 
          * ( u[ neighbourEntities.template getEntityIndex< 1,0,0 >( cellIndex ) ] - u[ cellIndex ] ) * mesh.template getSpaceStepsProducts< -1, 0, 0 >() * mesh.template getSpaceStepsProducts< -1, 0, 0 >() 
          + ( u[ neighbourEntities.template getEntityIndex< 0,0,1 >( cellIndex ) ] - u[ cellIndex ] ) * ( u[ neighbourEntities.template getEntityIndex< 0,0,1 >( cellIndex ) ] - u[ cellIndex ] )
   return sqrt( 1.0 + ( u[ neighbourEntities.template getEntityIndex< 0,1,0 >() ] - u[ cellIndex ] ) * 
          ( u[ neighbourEntities.template getEntityIndex< 0,1,0 >() ] - u[ cellIndex ] )
          * mesh.template getSpaceStepsProducts< 0, -1, 0 >() * mesh.template getSpaceStepsProducts< 0, -1, 0 >() + ( u[ neighbourEntities.template getEntityIndex< 1,0,0 >() ] - u[ cellIndex ] ) 
          * ( u[ neighbourEntities.template getEntityIndex< 1,0,0 >() ] - u[ cellIndex ] ) * mesh.template getSpaceStepsProducts< -1, 0, 0 >() * mesh.template getSpaceStepsProducts< -1, 0, 0 >() 
          + ( u[ neighbourEntities.template getEntityIndex< 0,0,1 >() ] - u[ cellIndex ] ) * ( u[ neighbourEntities.template getEntityIndex< 0,0,1 >() ] - u[ cellIndex ] )
            * mesh.template getSpaceStepsProducts< 0, 0, -1 >() * mesh.template getSpaceStepsProducts< 0, 0, -1 >() );
}
   
@@ -486,16 +487,16 @@ getValueStriped( const MeshType& mesh,
{
   const IndexType& cellIndex = entity.getIndex();
   const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();      
   return sqrt( this->eps + 0.5*( ( u[ neighbourEntities.template getEntityIndex< 0,1,0 >( cellIndex ) ] - u[ cellIndex ] ) * 
           ( u[ neighbourEntities.template getEntityIndex< 0,1,0 >( cellIndex ) ] - u[ cellIndex ] ) * mesh.template getSpaceStepsProducts< 0, -1, 0 >() * mesh.template getSpaceStepsProducts< 0, -1, 0 >() 
           + ( u[ neighbourEntities.template getEntityIndex< 1,0,0 >( cellIndex ) ] - u[ cellIndex ] ) * ( u[ neighbourEntities.template getEntityIndex< 1,0,0 >( cellIndex ) ] - u[ cellIndex ] ) 
           * mesh.template getSpaceStepsProducts< -1, 0, 0 >() * mesh.template getSpaceStepsProducts< -1, 0, 0 >() + ( - u[ neighbourEntities.template getEntityIndex< -1,0,0 >( cellIndex ) ] + u[ cellIndex ]) 
           * ( - u[ neighbourEntities.template getEntityIndex< -1,0,0 >( cellIndex ) ] + u[ cellIndex ]) * mesh.template getSpaceStepsProducts< -1, 0, 0 >() * mesh.template getSpaceStepsProducts< -1, 0, 0 >()
           + ( - u[ neighbourEntities.template getEntityIndex< 0,-1,0 >( cellIndex ) ] + u[ cellIndex ]) * 
           ( - u[ neighbourEntities.template getEntityIndex< 0,-1,0 >( cellIndex ) ] + u[ cellIndex ]) * mesh.template getSpaceStepsProducts< 0, -1, 0 >() * mesh.template getSpaceStepsProducts< 0, -1, 0 >() 
           + ( u[ neighbourEntities.template getEntityIndex< 0,0,1 >( cellIndex ) ] - u[ cellIndex ] ) * ( u[ neighbourEntities.template getEntityIndex< 0,0,1 >( cellIndex ) ] - u[ cellIndex ] )
           * mesh.template getSpaceStepsProducts< 0, 0, -1 >() * mesh.template getSpaceStepsProducts< 0, 0, -1 >() + ( - u[ neighbourEntities.template getEntityIndex< 0,0,-1 >( cellIndex ) ] + u[ cellIndex ]) * 
           ( - u[ neighbourEntities.template getEntityIndex< 0,0,-1 >( cellIndex ) ] + u[ cellIndex ]) * mesh.template getSpaceStepsProducts< 0, 0, -1 >() * mesh.template getSpaceStepsProducts< 0, 0, -1 >()
   return sqrt( this->eps + 0.5*( ( u[ neighbourEntities.template getEntityIndex< 0,1,0 >() ] - u[ cellIndex ] ) * 
           ( u[ neighbourEntities.template getEntityIndex< 0,1,0 >() ] - u[ cellIndex ] ) * mesh.template getSpaceStepsProducts< 0, -1, 0 >() * mesh.template getSpaceStepsProducts< 0, -1, 0 >() 
           + ( u[ neighbourEntities.template getEntityIndex< 1,0,0 >() ] - u[ cellIndex ] ) * ( u[ neighbourEntities.template getEntityIndex< 1,0,0 >() ] - u[ cellIndex ] ) 
           * mesh.template getSpaceStepsProducts< -1, 0, 0 >() * mesh.template getSpaceStepsProducts< -1, 0, 0 >() + ( - u[ neighbourEntities.template getEntityIndex< -1,0,0 >() ] + u[ cellIndex ]) 
           * ( - u[ neighbourEntities.template getEntityIndex< -1,0,0 >() ] + u[ cellIndex ]) * mesh.template getSpaceStepsProducts< -1, 0, 0 >() * mesh.template getSpaceStepsProducts< -1, 0, 0 >()
           + ( - u[ neighbourEntities.template getEntityIndex< 0,-1,0 >() ] + u[ cellIndex ]) * 
           ( - u[ neighbourEntities.template getEntityIndex< 0,-1,0 >() ] + u[ cellIndex ]) * mesh.template getSpaceStepsProducts< 0, -1, 0 >() * mesh.template getSpaceStepsProducts< 0, -1, 0 >() 
           + ( u[ neighbourEntities.template getEntityIndex< 0,0,1 >() ] - u[ cellIndex ] ) * ( u[ neighbourEntities.template getEntityIndex< 0,0,1 >() ] - u[ cellIndex ] )
           * mesh.template getSpaceStepsProducts< 0, 0, -1 >() * mesh.template getSpaceStepsProducts< 0, 0, -1 >() + ( - u[ neighbourEntities.template getEntityIndex< 0,0,-1 >() ] + u[ cellIndex ]) * 
           ( - u[ neighbourEntities.template getEntityIndex< 0,0,-1 >() ] + u[ cellIndex ]) * mesh.template getSpaceStepsProducts< 0, 0, -1 >() * mesh.template getSpaceStepsProducts< 0, 0, -1 >()
           ) );
}   
#endif	/* TNLONESIDEDIFFOPERATORQ_IMPL_H */
Loading