diff --git a/examples/advection/LaxFridrichs.h b/examples/advection/LaxFridrichs.h
index 4f38d390f61e9c3d15dd9c533ad4c621cff9f13d..d27679548e3a21f4237e2885e99f0d7b4ffb1fb9 100644
--- a/examples/advection/LaxFridrichs.h
+++ b/examples/advection/LaxFridrichs.h
@@ -36,12 +36,12 @@ class LaxFridrichs< Meshes::Grid< 1, MeshReal, Device, MeshIndex >, Real, Index,
typedef Functions::MeshFunction< MeshType > MeshFunctionType;
typedef VelocityFunction VelocityFunctionType;
typedef Functions::VectorField< Dimensions, VelocityFunctionType > VelocityFieldType;
- typedef SharedPointer< VelocityFieldType > VelocityFieldPointer;
+ typedef SharedPointer< VelocityFieldType, DeviceType > VelocityFieldPointer;
static void configSetup( Config::ConfigDescription& config,
const String& prefix = "" )
{
- config.addEntry< double >( "viscosity", "Value of artificial (numerical) viscosity in the Lax-Fridrichs shceme", 1.0 );
+ config.addEntry< double >( "viscosity", "Value of artificial (numerical) viscosity in the Lax-Fridrichs scheme", 1.0 );
}
LaxFridrichs()
@@ -54,7 +54,8 @@ class LaxFridrichs< Meshes::Grid< 1, MeshReal, Device, MeshIndex >, Real, Index,
const Config::ParameterContainer& parameters,
const String& prefix = "" )
{
- return this->velocityField->setup( meshPointer, parameters, prefix + "velocity-field-" );
+ this->artificialViscosity = parameters.getParameter< double >( "viscosity" );
+ return true;
}
static String getType();
@@ -127,12 +128,12 @@ class LaxFridrichs< Meshes::Grid< 2, MeshReal, Device, MeshIndex >, Real, Index,
typedef Functions::MeshFunction< MeshType > MeshFunctionType;
typedef VelocityFunction VelocityFunctionType;
typedef Functions::VectorField< Dimensions, VelocityFunctionType > VelocityFieldType;
- typedef SharedPointer< VelocityFieldType > VelocityFieldPointer;
+ typedef SharedPointer< VelocityFieldType, DeviceType > VelocityFieldPointer;
static void configSetup( Config::ConfigDescription& config,
const String& prefix = "" )
{
- config.addEntry< double >( "viscosity", "Value of artificial (numerical) viscosity in the Lax-Fridrichs shceme", 1.0 );
+ config.addEntry< double >( "viscosity", "Value of artificial (numerical) viscosity in the Lax-Fridrichs scheme", 1.0 );
}
LaxFridrichs()
@@ -145,7 +146,8 @@ class LaxFridrichs< Meshes::Grid< 2, MeshReal, Device, MeshIndex >, Real, Index,
const Config::ParameterContainer& parameters,
const String& prefix = "" )
{
- return this->velocityField->setup( meshPointer, parameters, prefix + "velocity-field-" );
+ this->artificialViscosity = parameters.getParameter< double >( "viscosity" );
+ return true;
}
static String getType();
@@ -224,12 +226,12 @@ class LaxFridrichs< Meshes::Grid< 3, MeshReal, Device, MeshIndex >, Real, Index,
typedef Functions::MeshFunction< MeshType > MeshFunctionType;
typedef VelocityFunction VelocityFunctionType;
typedef Functions::VectorField< Dimensions, VelocityFunctionType > VelocityFieldType;
- typedef SharedPointer< VelocityFieldType > VelocityFieldPointer;
+ typedef SharedPointer< VelocityFieldType, DeviceType > VelocityFieldPointer;
static void configSetup( Config::ConfigDescription& config,
const String& prefix = "" )
{
- config.addEntry< double >( "viscosity", "Value of artificial (numerical) viscosity in the Lax-Fridrichs shceme", 1.0 );
+ config.addEntry< double >( "viscosity", "Value of artificial (numerical) viscosity in the Lax-Fridrichs scheme", 1.0 );
}
LaxFridrichs()
@@ -242,7 +244,8 @@ class LaxFridrichs< Meshes::Grid< 3, MeshReal, Device, MeshIndex >, Real, Index,
const Config::ParameterContainer& parameters,
const String& prefix = "" )
{
- return this->velocityField->setup( meshPointer, parameters, prefix + "velocity-field-" );
+ this->artificialViscosity = parameters.getParameter< double >( "viscosity" );
+ return true;
}
static String getType();
diff --git a/examples/advection/LaxFridrichs_impl.h b/examples/advection/LaxFridrichs_impl.h
deleted file mode 100644
index dc65664b3b5b8ef599d61988924b3a705b6f94ca..0000000000000000000000000000000000000000
--- a/examples/advection/LaxFridrichs_impl.h
+++ /dev/null
@@ -1,361 +0,0 @@
-#ifndef LaxFridrichs_IMPL_H
-#define LaxFridrichs_IMPL_H
-
-namespace TNL {
-
-/****
- * 1D problem
- */
-template< typename MeshReal,
- typename Device,
- typename MeshIndex,
- typename Real,
- typename Index >
-String
-LaxFridrichs< Meshes::Grid< 1, MeshReal, Device, MeshIndex >, Real, Index >::
-getType()
-{
- return String( "LaxFridrichs< " ) +
- MeshType::getType() + ", " +
- TNL::getType< Real >() + ", " +
- TNL::getType< Index >() + " >";
-}
-
-template< typename MeshReal,
- typename Device,
- typename MeshIndex,
- typename Real,
- typename Index >
-template< typename MeshFunction, typename MeshEntity >
-__cuda_callable__
-Real
-LaxFridrichs< Meshes::Grid< 1, MeshReal, Device, MeshIndex >, Real, Index >::
-operator()( const MeshFunction& u,
- const MeshEntity& entity,
- const Real& time ) const
-{
- /****
- * Implement your explicit form of the differential operator here.
- * The following example is the Laplace operator approximated
- * by the Finite difference method.
- */
- static_assert( MeshEntity::entityDimensions == 1, "Wrong mesh entity dimensions." );
- static_assert( MeshFunction::getEntitiesDimensions() == 1, "Wrong preimage function" );
- const typename MeshEntity::template NeighbourEntities< 1 >& neighbourEntities = entity.getNeighbourEntities();
-
- const RealType& hxInverse = entity.getMesh().template getSpaceStepsProducts< -1 >();
- const IndexType& center = entity.getIndex();
- const IndexType& east = neighbourEntities.template getEntityIndex< 1 >();
- const IndexType& west = neighbourEntities.template getEntityIndex< -1 >();
- double a;
- a = this->advectionSpeedX;
- return (0.5 / this->tau ) * this->artificalViscosity *
- ( u[ west ]- 2.0 * u[ center ] + u[ east ] )
- - (a = this->advectionSpeedX * ( u[ east ] - u[west] ) ) * hxInverse * 0.5;
-
-}
-
-template< typename MeshReal,
- typename Device,
- typename MeshIndex,
- typename Real,
- typename Index >
-template< typename MeshEntity >
-__cuda_callable__
-Index
-LaxFridrichs< Meshes::Grid< 1, MeshReal, Device, MeshIndex >, Real, Index >::
-getLinearSystemRowLength( const MeshType& mesh,
- const IndexType& index,
- const MeshEntity& entity ) const
-{
- /****
- * Return a number of non-zero elements in a line (associated with given grid element) of
- * the linear system.
- * The following example is the Laplace operator approximated
- * by the Finite difference method.
- */
-
- return 2*Dimensions + 1;
-}
-
-template< typename MeshReal,
- typename Device,
- typename MeshIndex,
- typename Real,
- typename Index >
- template< typename MeshEntity, typename Vector, typename MatrixRow >
-__cuda_callable__
-void
-LaxFridrichs< Meshes::Grid< 1, MeshReal, Device, MeshIndex >, Real, Index >::
-updateLinearSystem( const RealType& time,
- const RealType& tau,
- const MeshType& mesh,
- const IndexType& index,
- const MeshEntity& entity,
- const MeshFunctionType& u,
- Vector& b,
- MatrixRow& matrixRow ) const
-{
- /****
- * Setup the non-zero elements of the linear system here.
- * The following example is the Laplace operator appriximated
- * by the Finite difference method.
- */
-
- const typename MeshEntity::template NeighbourEntities< 1 >& neighbourEntities = entity.getNeighbourEntities();
- const RealType& lambdaX = tau * entity.getMesh().template getSpaceStepsProducts< -2 >();
- const IndexType& center = entity.getIndex();
- const IndexType& east = neighbourEntities.template getEntityIndex< 1 >();
- const IndexType& west = neighbourEntities.template getEntityIndex< -1 >();
- matrixRow.setElement( 0, west, - lambdaX );
- matrixRow.setElement( 1, center, 2.0 * lambdaX );
- matrixRow.setElement( 2, east, - lambdaX );
-}
-
-/****
- * 2D problem
- */
-template< typename MeshReal,
- typename Device,
- typename MeshIndex,
- typename Real,
- typename Index >
-String
-LaxFridrichs< Meshes::Grid< 2, MeshReal, Device, MeshIndex >, Real, Index >::
-getType()
-{
- return String( "LaxFridrichs< " ) +
- MeshType::getType() + ", " +
- TNL::getType< Real >() + ", " +
- TNL::getType< Index >() + " >";
-}
-
-template< typename MeshReal,
- typename Device,
- typename MeshIndex,
- typename Real,
- typename Index >
-template< typename MeshFunction, typename MeshEntity >
-__cuda_callable__
-Real
-LaxFridrichs< Meshes::Grid< 2, MeshReal, Device, MeshIndex >, Real, Index >::
-operator()( const MeshFunction& u,
- const MeshEntity& entity,
- const Real& time ) const
-{
- /****
- * Implement your explicit form of the differential operator here.
- * The following example is the Laplace operator approximated
- * by the Finite difference method.
- */
- static_assert( MeshEntity::entityDimensions == 2, "Wrong mesh entity dimensions." );
- static_assert( MeshFunction::getEntitiesDimensions() == 2, "Wrong preimage function" );
- const typename MeshEntity::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
-
- const RealType& hxInverse = entity.getMesh().template getSpaceStepsProducts< -1, 0 >();
- const RealType& hyInverse = entity.getMesh().template getSpaceStepsProducts< 0, -1 >();
- const IndexType& center = entity.getIndex();
- const IndexType& east = neighbourEntities.template getEntityIndex< 1, 0 >();
- const IndexType& west = neighbourEntities.template getEntityIndex< -1, 0 >();
- const IndexType& north = neighbourEntities.template getEntityIndex< 0, 1 >();
- const IndexType& south = neighbourEntities.template getEntityIndex< 0, -1 >();
- double a;
- double b;
- a = this->advectionSpeedX;
- b = this->advectionSpeedY;
- return ( 0.25 / this->tau ) * this->artificalViscosity *
- ( u[ west ] + u[ east ] + u[ south ] + u[ north ] - 4 * u[ center ] ) -
- (a * ( u[ east ] - u[west] ) ) * hxInverse * 0.5 -
- (b * ( u[ north ] - u[ south ] ) ) * hyInverse * 0.5;
-
-
-}
-
-template< typename MeshReal,
- typename Device,
- typename MeshIndex,
- typename Real,
- typename Index >
-template< typename MeshEntity >
-__cuda_callable__
-Index
-LaxFridrichs< Meshes::Grid< 2, MeshReal, Device, MeshIndex >, Real, Index >::
-getLinearSystemRowLength( const MeshType& mesh,
- const IndexType& index,
- const MeshEntity& entity ) const
-{
- /****
- * Return a number of non-zero elements in a line (associated with given grid element) of
- * the linear system.
- * The following example is the Laplace operator approximated
- * by the Finite difference method.
- */
-
- return 2*Dimensions + 1;
-}
-
-template< typename MeshReal,
- typename Device,
- typename MeshIndex,
- typename Real,
- typename Index >
- template< typename MeshEntity, typename Vector, typename MatrixRow >
-__cuda_callable__
-void
-LaxFridrichs< Meshes::Grid< 2, MeshReal, Device, MeshIndex >, Real, Index >::
-updateLinearSystem( const RealType& time,
- const RealType& tau,
- const MeshType& mesh,
- const IndexType& index,
- const MeshEntity& entity,
- const MeshFunctionType& u,
- Vector& b,
- MatrixRow& matrixRow ) const
-{
- /****
- * Setup the non-zero elements of the linear system here.
- * The following example is the Laplace operator appriximated
- * by the Finite difference method.
- */
-
- const typename MeshEntity::template NeighbourEntities< 2 >& neighbourEntities = entity.getNeighbourEntities();
- const RealType& lambdaX = tau * entity.getMesh().template getSpaceStepsProducts< -2, 0 >();
- const RealType& lambdaY = tau * entity.getMesh().template getSpaceStepsProducts< 0, -2 >();
- const IndexType& center = entity.getIndex();
- const IndexType& east = neighbourEntities.template getEntityIndex< 1, 0 >();
- const IndexType& west = neighbourEntities.template getEntityIndex< -1, 0 >();
- const IndexType& north = neighbourEntities.template getEntityIndex< 0, 1 >();
- const IndexType& south = neighbourEntities.template getEntityIndex< 0, -1 >();
- matrixRow.setElement( 0, south, -lambdaY );
- matrixRow.setElement( 1, west, -lambdaX );
- matrixRow.setElement( 2, center, 2.0 * ( lambdaX + lambdaY ) );
- matrixRow.setElement( 3, east, -lambdaX );
- matrixRow.setElement( 4, north, -lambdaY );
-}
-
-/****
- * 3D problem
- */
-template< typename MeshReal,
- typename Device,
- typename MeshIndex,
- typename Real,
- typename Index >
-String
-LaxFridrichs< Meshes::Grid< 3, MeshReal, Device, MeshIndex >, Real, Index >::
-getType()
-{
- return String( "LaxFridrichs< " ) +
- MeshType::getType() + ", " +
- TNL::getType< Real >() + ", " +
- TNL::getType< Index >() + " >";
-}
-
-template< typename MeshReal,
- typename Device,
- typename MeshIndex,
- typename Real,
- typename Index >
-template< typename MeshFunction, typename MeshEntity >
-__cuda_callable__
-Real
-LaxFridrichs< Meshes::Grid< 3, MeshReal, Device, MeshIndex >, Real, Index >::
-operator()( const MeshFunction& u,
- const MeshEntity& entity,
- const Real& time ) const
-{
- /****
- * Implement your explicit form of the differential operator here.
- * The following example is the Laplace operator approximated
- * by the Finite difference method.
- */
- static_assert( MeshEntity::entityDimensions == 3, "Wrong mesh entity dimensions." );
- static_assert( MeshFunction::getEntitiesDimensions() == 3, "Wrong preimage function" );
- const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
-
- const RealType& hxSquareInverse = entity.getMesh().template getSpaceStepsProducts< -2, 0, 0 >();
- const RealType& hySquareInverse = entity.getMesh().template getSpaceStepsProducts< 0, -2, 0 >();
- const RealType& hzSquareInverse = entity.getMesh().template getSpaceStepsProducts< 0, 0, -2 >();
- const IndexType& center = entity.getIndex();
- const IndexType& east = neighbourEntities.template getEntityIndex< 1, 0, 0 >();
- const IndexType& west = neighbourEntities.template getEntityIndex< -1, 0, 0 >();
- const IndexType& north = neighbourEntities.template getEntityIndex< 0, 1, 0 >();
- const IndexType& south = neighbourEntities.template getEntityIndex< 0, -1, 0 >();
- const IndexType& up = neighbourEntities.template getEntityIndex< 0, 0, 1 >();
- const IndexType& down = neighbourEntities.template getEntityIndex< 0, 0, -1 >();
- return ( u[ west ] - 2.0 * u[ center ] + u[ east ] ) * hxSquareInverse +
- ( u[ south ] - 2.0 * u[ center ] + u[ north ] ) * hySquareInverse +
- ( u[ up ] - 2.0 * u[ center ] + u[ down ] ) * hzSquareInverse;
-}
-
-template< typename MeshReal,
- typename Device,
- typename MeshIndex,
- typename Real,
- typename Index >
-template< typename MeshEntity >
-__cuda_callable__
-Index
-LaxFridrichs< Meshes::Grid< 3, MeshReal, Device, MeshIndex >, Real, Index >::
-getLinearSystemRowLength( const MeshType& mesh,
- const IndexType& index,
- const MeshEntity& entity ) const
-{
- /****
- * Return a number of non-zero elements in a line (associated with given grid element) of
- * the linear system.
- * The following example is the Laplace operator approximated
- * by the Finite difference method.
- */
-
- //return 2*Dimensions + 1;
-}
-
-template< typename MeshReal,
- typename Device,
- typename MeshIndex,
- typename Real,
- typename Index >
- template< typename MeshEntity, typename Vector, typename MatrixRow >
-__cuda_callable__
-void
-LaxFridrichs< Meshes::Grid< 3, MeshReal, Device, MeshIndex >, Real, Index >::
-updateLinearSystem( const RealType& time,
- const RealType& tau,
- const MeshType& mesh,
- const IndexType& index,
- const MeshEntity& entity,
- const MeshFunctionType& u,
- Vector& b,
- MatrixRow& matrixRow ) const
-{
- /****
- * Setup the non-zero elements of the linear system here.
- * The following example is the Laplace operator appriximated
- * by the Finite difference method.
- */
-
- /* const typename MeshEntity::template NeighbourEntities< 3 >& neighbourEntities = entity.getNeighbourEntities();
- const RealType& lambdaX = tau * entity.getMesh().template getSpaceStepsProducts< -2, 0, 0 >();
- const RealType& lambdaY = tau * entity.getMesh().template getSpaceStepsProducts< 0, -2, 0 >();
- const RealType& lambdaZ = tau * entity.getMesh().template getSpaceStepsProducts< 0, 0, -2 >();
- const IndexType& center = entity.getIndex();
- const IndexType& east = neighbourEntities.template getEntityIndex< 1, 0, 0 >();
- const IndexType& west = neighbourEntities.template getEntityIndex< -1, 0, 0 >();
- const IndexType& north = neighbourEntities.template getEntityIndex< 0, 1, 0 >();
- const IndexType& south = neighbourEntities.template getEntityIndex< 0, -1, 0 >();
- const IndexType& up = neighbourEntities.template getEntityIndex< 0, 0, 1 >();
- const IndexType& down = neighbourEntities.template getEntityIndex< 0, 0, -1 >();
- matrixRow.setElement( 0, down, -lambdaZ );
- matrixRow.setElement( 1, south, -lambdaY );
- matrixRow.setElement( 2, west, -lambdaX );
- matrixRow.setElement( 3, center, 2.0 * ( lambdaX + lambdaY + lambdaZ ) );
- matrixRow.setElement( 4, east, -lambdaX );
- matrixRow.setElement( 5, north, -lambdaY );
- matrixRow.setElement( 6, up, -lambdaZ );*/
-}
-
-} // namespace TNL
-
-#endif /* LaxFridrichsIMPL_H */
-
diff --git a/examples/advection/advection.h b/examples/advection/advection.h
index 7e9c6c64296c2f4975729f36eb18d59ea23ba2e2..62d7bf0e7ba8a8ac9e7a6cb98f6d4d87d6615d09 100644
--- a/examples/advection/advection.h
+++ b/examples/advection/advection.h
@@ -37,7 +37,7 @@ template< typename ConfigTag >class advectionConfig
Functions::VectorField< 3, Functions::Analytic::Constant< 3 > >::configSetup( config, "velocity-field-" );
typedef Meshes::Grid< 3 > MeshType;
- LaxFridrichs< MeshType >::ConfigSetup( config, "lax-fridrichs" );
+ LaxFridrichs< MeshType >::configSetup( config, "lax-fridrichs" );
config.addEntry< String >( "boundary-conditions-type", "Choose the boundary conditions type.", "dirichlet");
config.addEntryEnum< String >( "dirichlet" );
@@ -62,8 +62,10 @@ class advectionSetter
static bool run( const Config::ParameterContainer & parameters )
{
- enum { Dimensions = MeshType::getMeshDimensions() };
- typedef LaxFridrichs< MeshType, Real, Index > ApproximateOperator;
+ static const int Dimensions = MeshType::getMeshDimensions();
+ typedef Functions::Analytic::Constant< Dimensions, RealType > ConstantFunctionType;
+ typedef Functions::VectorField< Dimensions, ConstantFunctionType > VelocityFieldType;
+ typedef LaxFridrichs< MeshType, Real, Index, VelocityFieldType > ApproximateOperator;
typedef advectionRhs< MeshType, Real > RightHandSide;
typedef Containers::StaticVector < MeshType::getMeshDimensions(), Real > Vertex;
diff --git a/examples/advection/advectionProblem.h b/examples/advection/advectionProblem.h
index 73d34dcf2e050c00efc73d669485fcc3997c3df3..3b9c3a3f35f116d5f9e2266a52c7a17e7a72c71f 100644
--- a/examples/advection/advectionProblem.h
+++ b/examples/advection/advectionProblem.h
@@ -12,7 +12,8 @@ namespace TNL {
template< typename Mesh,
typename BoundaryCondition,
typename RightHandSide,
- typename DifferentialOperator >
+ typename DifferentialOperator,
+ typename VelocityField = Functions::MeshFunction< Mesh > >
class advectionProblem:
public PDEProblem< Mesh,
typename DifferentialOperator::RealType,
@@ -30,6 +31,8 @@ public PDEProblem< Mesh,
typedef SharedPointer< DifferentialOperator > DifferentialOperatorPointer;
typedef SharedPointer< BoundaryCondition > BoundaryConditionPointer;
typedef SharedPointer< RightHandSide, DeviceType > RightHandSidePointer;
+ typedef VelocityField VelocityFieldType;
+ typedef SharedPointer< VelocityFieldType, DeviceType > VelocityFieldPointer;
using typename BaseType::MeshType;
using typename BaseType::MeshPointer;
@@ -93,7 +96,9 @@ public PDEProblem< Mesh,
BoundaryConditionPointer boundaryConditionPointer;
- RightHandSidePointer rightHandSidePointer;
+ RightHandSidePointer rightHandSidePointer;
+
+ VelocityFieldPointer velocityField;
String velocityType;
};
diff --git a/examples/advection/advectionProblem_impl.h b/examples/advection/advectionProblem_impl.h
index 9b0f7ea672328a7661ac3c96dd8d309f0b4a3a6f..f2de7ea03ee13be4b4b41775c89f974b7c41dc4f 100644
--- a/examples/advection/advectionProblem_impl.h
+++ b/examples/advection/advectionProblem_impl.h
@@ -12,9 +12,10 @@ namespace TNL {
template< typename Mesh,
typename BoundaryCondition,
typename RightHandSide,
- typename DifferentialOperator >
+ typename DifferentialOperator,
+ typename VelocityField >
String
-advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator >::
+advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator, VelocityField >::
getTypeStatic()
{
return String( "advectionProblem< " ) + Mesh :: getTypeStatic() + " >";
@@ -23,9 +24,10 @@ getTypeStatic()
template< typename Mesh,
typename BoundaryCondition,
typename RightHandSide,
- typename DifferentialOperator >
+ typename DifferentialOperator,
+ typename VelocityField >
String
-advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator >::
+advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator, VelocityField >::
getPrologHeader() const
{
return String( "advection" );
@@ -34,9 +36,10 @@ getPrologHeader() const
template< typename Mesh,
typename BoundaryCondition,
typename RightHandSide,
- typename DifferentialOperator >
+ typename DifferentialOperator,
+ typename VelocityField >
void
-advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator >::
+advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator, VelocityField >::
writeProlog( Logger& logger, const Config::ParameterContainer& parameters ) const
{
/****
@@ -48,9 +51,10 @@ writeProlog( Logger& logger, const Config::ParameterContainer& parameters ) cons
template< typename Mesh,
typename BoundaryCondition,
typename RightHandSide,
- typename DifferentialOperator >
+ typename DifferentialOperator,
+ typename VelocityField >
bool
-advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator >::
+advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator, VelocityField >::
setup( const MeshPointer& meshPointer,
const Config::ParameterContainer& parameters,
const String& prefix )
@@ -62,8 +66,9 @@ setup( const MeshPointer& meshPointer,
differentialOperatorPointer->setAdvectionSpeedX(advectionSpeedX);
const double advectionSpeedY = parameters.getParameter< double >( "advection-speedY" );
differentialOperatorPointer->setAdvectionSpeedY(advectionSpeedY);*/
-
- if( ! this->differentialOperatorPointer->setup( meshPointer, parameters, prefix + "advection-" ) ||
+
+ if( ! this->velocityField->setup( meshPointer, parameters, prefix + "velocity-field-" ) ||
+ ! this->differentialOperatorPointer->setup( meshPointer, parameters, prefix + "advection-" ) ||
! this->boundaryConditionPointer->setup( meshPointer, parameters, prefix + "boundary-conditions-" ) ||
! this->rightHandSidePointer->setup( parameters, prefix + "right-hand-side-" ) )
return false;
@@ -73,9 +78,10 @@ setup( const MeshPointer& meshPointer,
template< typename Mesh,
typename BoundaryCondition,
typename RightHandSide,
- typename DifferentialOperator >
-typename advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator >::IndexType
-advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator >::
+ typename DifferentialOperator,
+ typename VelocityField >
+typename advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator, VelocityField >::IndexType
+advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator, VelocityField >::
getDofs( const MeshPointer& mesh ) const
{
/****
@@ -88,9 +94,10 @@ getDofs( const MeshPointer& mesh ) const
template< typename Mesh,
typename BoundaryCondition,
typename RightHandSide,
- typename DifferentialOperator >
+ typename DifferentialOperator,
+ typename VelocityField >
void
-advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator >::
+advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator, VelocityField >::
bindDofs( const MeshPointer& meshPointer,
DofVectorPointer& dofVector )
{
@@ -101,9 +108,10 @@ bindDofs( const MeshPointer& meshPointer,
template< typename Mesh,
typename BoundaryCondition,
typename RightHandSide,
- typename DifferentialOperator >
+ typename DifferentialOperator,
+ typename VelocityField >
bool
-advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator >::
+advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator, VelocityField >::
setInitialCondition( const Config::ParameterContainer& parameters,
const MeshPointer& meshPointer,
DofVectorPointer& dofs,
@@ -122,10 +130,11 @@ setInitialCondition( const Config::ParameterContainer& parameters,
template< typename Mesh,
typename BoundaryCondition,
typename RightHandSide,
- typename DifferentialOperator >
+ typename DifferentialOperator,
+ typename VelocityField >
template< typename Matrix >
bool
-advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator >::
+advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator, VelocityField >::
setupLinearSystem( const MeshPointer& mesh,
Matrix& matrix )
{
@@ -148,9 +157,10 @@ setupLinearSystem( const MeshPointer& mesh,
template< typename Mesh,
typename BoundaryCondition,
typename RightHandSide,
- typename DifferentialOperator >
+ typename DifferentialOperator,
+ typename VelocityField >
bool
-advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator >::
+advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator, VelocityField >::
makeSnapshot( const RealType& time,
const IndexType& step,
const MeshPointer& mesh,
@@ -171,9 +181,10 @@ makeSnapshot( const RealType& time,
template< typename Mesh,
typename BoundaryCondition,
typename RightHandSide,
- typename DifferentialOperator >
+ typename DifferentialOperator,
+ typename VelocityField >
void
-advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator >::
+advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator, VelocityField >::
getExplicitRHS( const RealType& time,
const RealType& tau,
const MeshPointer& mesh,
@@ -212,10 +223,11 @@ getExplicitRHS( const RealType& time,
template< typename Mesh,
typename BoundaryCondition,
typename RightHandSide,
- typename DifferentialOperator >
+ typename DifferentialOperator,
+ typename VelocityField >
template< typename Matrix >
void
-advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator >::
+advectionProblem< Mesh, BoundaryCondition, RightHandSide, DifferentialOperator, VelocityField >::
assemblyLinearSystem( const RealType& time,
const RealType& tau,
const MeshPointer& mesh,