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

Fixing finite differences.

parent 88cd0956
Loading
Loading
Loading
Loading
+1 −3
Original line number Diff line number Diff line
@@ -39,7 +39,7 @@ template< int Dimensions,
          int ZDifference,
          typename Real,
          typename Index >
class tnlBackwardFiniteDifference< tnlGrid< Dimensions, MeshReal, MeshDevice, MeshIndex >, Real, Index >
class tnlBackwardFiniteDifference< tnlGrid< Dimensions, MeshReal, MeshDevice, MeshIndex >, XDifference, YDifference, ZDifference, Real, Index >
: tnlDomain< Dimensions, MeshInteriorDomain >
{
   public:
@@ -49,8 +49,6 @@ class tnlBackwardFiniteDifference< tnlGrid< Dimensions, MeshReal, MeshDevice, Me
      typedef MeshDevice DeviceType;
      typedef Index IndexType;      
      
      static const int Dimensions = MeshType::meshDimensions;
      
      static constexpr int getMeshDimensions() { return Dimensions; }
      
      template< typename MeshFunction, typename MeshEntity >
+1 −3
Original line number Diff line number Diff line
@@ -39,7 +39,7 @@ template< int Dimensions,
          int ZDifference,
          typename Real,
          typename Index >
class tnlCentralFiniteDifference< tnlGrid< Dimensions, MeshReal, MeshDevice, MeshIndex >, Real, Index >
class tnlCentralFiniteDifference< tnlGrid< Dimensions, MeshReal, MeshDevice, MeshIndex >, XDifference, YDifference, ZDifference, Real, Index >
: tnlDomain< Dimensions, MeshInteriorDomain >
{
   public:
@@ -49,8 +49,6 @@ class tnlCentralFiniteDifference< tnlGrid< Dimensions, MeshReal, MeshDevice, Mes
      typedef MeshDevice DeviceType;
      typedef Index IndexType;      
      
      static const int Dimensions = MeshType::meshDimensions;
      
      static constexpr int getMeshDimensions() { return Dimensions; }
      
      template< typename MeshFunction, typename MeshEntity >
+8 −8
Original line number Diff line number Diff line
@@ -69,8 +69,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 1 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().getSpaceStepsProducts< -1 >();
         const typename MeshEntity::template NeighbourEntities< 1 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().template getSpaceStepsProducts< -1 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u[ neighbourEntities.template getEntityIndex< 1 >()] - u_c ) * hxDiv;
      }            
@@ -95,8 +95,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 1 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().getSpaceStepsProducts< -1 >();
         const typename MeshEntity::template NeighbourEntities< 1 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().template getSpaceStepsProducts< -1 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u_c - u[ neighbourEntities.template getEntityIndex< -1 >()] ) * hxDiv;
      }            
@@ -121,8 +121,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 1 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().getSpaceStepsProducts< -1 >();
         const typename MeshEntity::template NeighbourEntities< 1 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().template getSpaceStepsProducts< -1 >();
         return ( u[ neighbourEntities.template getEntityIndex< 1 >() ] -
                  u[ neighbourEntities.template getEntityIndex< -1 >() ] ) * ( 0.5 * hxDiv );
      }            
@@ -147,8 +147,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 1 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxSquareDiv = entity.getMesh().getSpaceStepsProducts< -2 >();
         const typename MeshEntity::template NeighbourEntities< 1 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxSquareDiv = entity.getMesh().template getSpaceStepsProducts< -2 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u[ neighbourEntities.template getEntityIndex< 1 >() ] -
                  2.0 * u_c +
+16 −16
Original line number Diff line number Diff line
@@ -69,8 +69,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().getSpaceStepsProducts< -1, 0 >();
         const typename MeshEntity::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().template getSpaceStepsProducts< -1, 0 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u[ neighbourEntities.template getEntityIndex< 1, 0 >()] - u_c ) * hxDiv;
      }            
@@ -92,8 +92,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hyDiv = entity.getMesh().getSpaceStepsProducts< 0, -1 >();
         const typename MeshEntity::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hyDiv = entity.getMesh().template getSpaceStepsProducts< 0, -1 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u[ neighbourEntities.template getEntityIndex< 0, 1 >()] - u_c ) * hyDiv;
      }            
@@ -118,8 +118,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().getSpaceStepsProducts< -1,  0 >();
         const typename MeshEntity::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().template getSpaceStepsProducts< -1,  0 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u_c - u[ neighbourEntities.template getEntityIndex< -1, 0 >()] ) * hxDiv;
      }            
@@ -141,8 +141,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hyDiv = entity.getMesh().getSpaceStepsProducts< 0, -1 >();
         const typename MeshEntity::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hyDiv = entity.getMesh().template getSpaceStepsProducts< 0, -1 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u_c - u[ neighbourEntities.template getEntityIndex< 0, -1 >()] ) * hyDiv;
      }            
@@ -167,8 +167,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().getSpaceStepsProducts< -1, 0 >();
         const typename MeshEntity::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().template getSpaceStepsProducts< -1, 0 >();
         return ( u[ neighbourEntities.template getEntityIndex< 1, 0 >() ] -
                  u[ neighbourEntities.template getEntityIndex< -1, 0 >() ] ) * ( 0.5 * hxDiv );
      }            
@@ -190,8 +190,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hyDiv = entity.getMesh().getSpaceStepsProducts< 0, -1 >();
         const typename MeshEntity::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hyDiv = entity.getMesh().template getSpaceStepsProducts< 0, -1 >();
         return ( u[ neighbourEntities.template getEntityIndex< 0,  1 >() ] -
                  u[ neighbourEntities.template getEntityIndex< 0, -1 >() ] ) * ( 0.5 * hyDiv );
      }            
@@ -217,8 +217,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxSquareDiv = entity.getMesh().getSpaceStepsProducts< -2, 0 >();
         const typename MeshEntity::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxSquareDiv = entity.getMesh().template getSpaceStepsProducts< -2, 0 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u[ neighbourEntities.template getEntityIndex<  1, 0 >() ] -
                  2.0 * u_c +
@@ -242,8 +242,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hySquareDiv = entity.getMesh().getSpaceStepsProducts< 0, -2 >();
         const typename MeshEntity::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hySquareDiv = entity.getMesh().template getSpaceStepsProducts< 0, -2 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u[ neighbourEntities.template getEntityIndex< 0,  1 >() ] -
                  2.0 * u_c +
+24 −24
Original line number Diff line number Diff line
@@ -37,8 +37,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().getSpaceStepsProducts< -1, 0, 0 >();
         const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().template getSpaceStepsProducts< -1, 0, 0 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u[ neighbourEntities.template getEntityIndex< 1, 0, 0 >()] - u_c ) * hxDiv;
      }            
@@ -60,8 +60,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hyDiv = entity.getMesh().getSpaceStepsProducts< 0, -1, 0 >();
         const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hyDiv = entity.getMesh().template getSpaceStepsProducts< 0, -1, 0 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u[ neighbourEntities.template getEntityIndex< 0, 1, 0 >()] - u_c ) * hyDiv;
      }            
@@ -83,8 +83,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hzDiv = entity.getMesh().getSpaceStepsProducts< 0, 0, -1 >();
         const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hzDiv = entity.getMesh().template getSpaceStepsProducts< 0, 0, -1 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u[ neighbourEntities.template getEntityIndex< 0, 0, 1 >()] - u_c ) * hzDiv;
      }            
@@ -109,8 +109,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().getSpaceStepsProducts< -1,  0 >();
         const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().template getSpaceStepsProducts< -1,  0 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u_c - u[ neighbourEntities.template getEntityIndex< -1, 0, 0 >()] ) * hxDiv;
      }            
@@ -132,8 +132,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hyDiv = entity.getMesh().getSpaceStepsProducts< 0, -1, 0 >();
         const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hyDiv = entity.getMesh().template getSpaceStepsProducts< 0, -1, 0 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u_c - u[ neighbourEntities.template getEntityIndex< 0, -1, 0 >()] ) * hyDiv;
      }            
@@ -155,8 +155,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hzDiv = entity.getMesh().getSpaceStepsProducts< 0, 0, -1 >();
         const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hzDiv = entity.getMesh().template getSpaceStepsProducts< 0, 0, -1 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u_c - u[ neighbourEntities.template getEntityIndex< 0, 0, -1 >()] ) * hzDiv;
      }            
@@ -181,8 +181,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().getSpaceStepsProducts< -1, 0, 0 >();
         const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxDiv = entity.getMesh().template getSpaceStepsProducts< -1, 0, 0 >();
         return ( u[ neighbourEntities.template getEntityIndex< 1, 0, 0 >() ] -
                  u[ neighbourEntities.template getEntityIndex< -1, 0, 0 >() ] ) * ( 0.5 * hxDiv );
      }            
@@ -204,8 +204,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hyDiv = entity.getMesh().getSpaceStepsProducts< 0, -1, 0 >();
         const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hyDiv = entity.getMesh().template getSpaceStepsProducts< 0, -1, 0 >();
         return ( u[ neighbourEntities.template getEntityIndex< 0, 1, 0 >() ] -
                  u[ neighbourEntities.template getEntityIndex< 0, -1, 0 >() ] ) * ( 0.5 * hyDiv );
      }            
@@ -227,8 +227,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hzDiv = entity.getMesh().getSpaceStepsProducts< 0, 0, -1 >();
         const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hzDiv = entity.getMesh().template getSpaceStepsProducts< 0, 0, -1 >();
         return ( u[ neighbourEntities.template getEntityIndex< 0, 0, 1 >() ] -
                  u[ neighbourEntities.template getEntityIndex< 0, 0, -1 >() ] ) * ( 0.5 * hzDiv );
      }            
@@ -253,8 +253,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxSquareDiv = entity.getMesh().getSpaceStepsProducts< -2, 0, 0 >();
         const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hxSquareDiv = entity.getMesh().template getSpaceStepsProducts< -2, 0, 0 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u[ neighbourEntities.template getEntityIndex<  1, 0, 0 >() ] -
                  2.0 * u_c +
@@ -278,8 +278,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hySquareDiv = entity.getMesh().getSpaceStepsProducts< 0, -2, 0 >();
         const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hySquareDiv = entity.getMesh().template getSpaceStepsProducts< 0, -2, 0 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u[ neighbourEntities.template getEntityIndex< 0,  1, 0 >() ] -
                  2.0 * u_c +
@@ -303,8 +303,8 @@ class tnlFiniteDifferences<
                            const MeshEntity& entity,
                            const Real& time = 0 )
      {         
         const typename EntityType::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hzSquareDiv = entity.getMesh().getSpaceStepsProducts< 0, 0, -2 >();
         const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
         const Real& hzSquareDiv = entity.getMesh().template getSpaceStepsProducts< 0, 0, -2 >();
         const Real& u_c = u[ entity.getIndex() ];
         return ( u[ neighbourEntities.template getEntityIndex< 0, 0,  1 >() ] -
                  2.0 * u_c +
Loading