Fixing the simple-problem example. (2be84baa) · Commits · TNL / tnl-dev

examples/navier-stokes/navierStokesSolver.h

+15 −0

Original line number	Diff line number	Diff line
		@@ -67,6 +67,21 @@ class navierStokesSolver

		bool setMeshGeometry( tnlLinearGridGeometry< 2, RealType, DeviceType, IndexType >& geometry ) const;

		template< typename InitMesh >
		bool initMesh( InitMesh& mesh, const tnlParameterContainer& parameters ) const;

		template< typename Real, typename Device, typename Index, template< int, typename, typename, typename > class Geometry >
		bool initMesh( tnlGrid< 1, Real, Device, Index, Geometry >& mesh,
		const tnlParameterContainer& parameters ) const;

		template< typename Real, typename Device, typename Index, template< int, typename, typename, typename > class Geometry >
		bool initMesh( tnlGrid< 2, Real, Device, Index, Geometry >& mesh,
		const tnlParameterContainer& parameters ) const;

		template< typename Real, typename Device, typename Index, template< int, typename, typename, typename > class Geometry >
		bool initMesh( tnlGrid< 3, Real, Device, Index, Geometry >& mesh,
		const tnlParameterContainer& parameters ) const;

		bool init( const tnlParameterContainer& parameters );

		bool setInitialCondition( const tnlParameterContainer& parameters );

examples/navier-stokes/navierStokesSolver_impl.h

+44 −19

Original line number	Diff line number	Diff line
		@@ -72,6 +72,46 @@ bool navierStokesSolver< Mesh, EulerScheme > :: setMeshGeometry( tnlLinearGridGe
		return true;
		}

		template< typename Mesh, typename EulerScheme >
		template< typename Real, typename Device, typename Index,
		template< int, typename, typename, typename > class Geometry >
		bool navierStokesSolver< Mesh, EulerScheme >::initMesh( tnlGrid< 1, Real, Device, Index, Geometry >& mesh,
		const tnlParameterContainer& parameters ) const
		{

		}

		template< typename Mesh, typename EulerScheme >
		template< typename Real, typename Device, typename Index,
		template< int, typename, typename, typename > class Geometry >
		bool navierStokesSolver< Mesh, EulerScheme >::initMesh( tnlGrid< 2, Real, Device, Index, Geometry >& mesh,
		const tnlParameterContainer& parameters ) const
		{
		tnlTuple< 2, IndexType > meshes;
		meshes.x() = parameters.GetParameter< int >( "x-size" );
		meshes.y() = parameters.GetParameter< int >( "y-size" );
		if( meshes.x() <= 0 )
		{
		cerr << "Error: x-size must be positive integer number! It is " << meshes. x() << " now." << endl;
		return false;
		}
		if( meshes.y() <= 0 )
		{
		cerr << "Error: y-size must be positive integer number! It is " << meshes. y() << " now." << endl;
		return false;
		}
		mesh.setDimensions( meshes. x(), meshes. y() );
		return this->setMeshGeometry( mesh.getGeometry() );
		}

		template< typename Mesh, typename EulerScheme >
		template< typename Real, typename Device, typename Index,
		template< int, typename, typename, typename > class Geometry >
		bool navierStokesSolver< Mesh, EulerScheme >::initMesh( tnlGrid< 3, Real, Device, Index, Geometry >& mesh,
		const tnlParameterContainer& parameters ) const
		{

		}

		template< typename Mesh, typename EulerScheme >
		bool navierStokesSolver< Mesh, EulerScheme >::init( const tnlParameterContainer& parameters )
		@@ -109,23 +149,8 @@ bool navierStokesSolver< Mesh, EulerScheme > :: init( const tnlParameterContaine
		/****
		* Set-up the space discretization
		*/
		tnlTuple< 2, IndexType > meshes;
		meshes. x() = parameters. GetParameter< int >( "x-size" );
		meshes. y() = parameters. GetParameter< int >( "y-size" );
		if( meshes. x() <= 0 )
		{
		cerr << "Error: x-size must be positive integer number! It is " << meshes. x() << " now." << endl;
		return false;
		}
		if( meshes. y() <= 0 )
		{
		cerr << "Error: y-size must be positive integer number! It is " << meshes. y() << " now." << endl;
		if( ! this->initMesh( this->mesh, parameters ) )
		return false;
		}
		this -> mesh. setDimensions( meshes. x(), meshes. y() );
		this -> setMeshGeometry( this -> mesh. getGeometry() );
		RealType hx = this -> mesh. getParametricStep(). x();
		RealType hy = this -> mesh. getParametricStep(). y();
		mesh.refresh();
		mesh.save( tnlString( "mesh.tnl" ) );
		nsSolver.setMesh( this->mesh );

examples/simple-solver/simpleProblemSolver.h

+16 −1

Original line number	Diff line number	Diff line
		@@ -23,7 +23,7 @@
		#include <solvers/tnlSolverMonitor.h>
		#include <core/tnlLogger.h>
		#include <core/vectors/tnlVector.h>
		#include <core/tnlSharedVector.h>
		#include <core/vectors/tnlSharedVector.h>

		template< typename Mesh >
		class simpleProblemSolver
		@@ -45,6 +45,21 @@ class simpleProblemSolver
		void writeProlog( tnlLogger& logger,
		const tnlParameterContainer& parameters ) const;

		template< typename InitMesh >
		bool initMesh( InitMesh& mesh, const tnlParameterContainer& parameters ) const;

		template< typename Real, typename Device, typename Index, template< int, typename, typename, typename > class Geometry >
		bool initMesh( tnlGrid< 1, Real, Device, Index, Geometry >& mesh,
		const tnlParameterContainer& parameters ) const;

		template< typename Real, typename Device, typename Index, template< int, typename, typename, typename > class Geometry >
		bool initMesh( tnlGrid< 2, Real, Device, Index, Geometry >& mesh,
		const tnlParameterContainer& parameters ) const;

		template< typename Real, typename Device, typename Index, template< int, typename, typename, typename > class Geometry >
		bool initMesh( tnlGrid< 3, Real, Device, Index, Geometry >& mesh,
		const tnlParameterContainer& parameters ) const;

		bool init( const tnlParameterContainer& parameters );

		bool setInitialCondition( const tnlParameterContainer& parameters );

examples/simple-solver/simpleProblemSolver_impl.h

+83 −11

Original line number	Diff line number	Diff line
		@@ -59,6 +59,76 @@ void simpleProblemSolver< Mesh> :: writeProlog( tnlLogger& logger,
		logger. WriteParameter< int >( "Simple parameter:", 1 );
		}

		template< typename Mesh >
		template< typename Real, typename Device, typename Index,
		template< int, typename, typename, typename > class Geometry >
		bool simpleProblemSolver< Mesh >::initMesh( tnlGrid< 1, Real, Device, Index, Geometry >& mesh,
		const tnlParameterContainer& parameters ) const
		{
		tnlTuple< 1, IndexType > meshes;
		meshes.x() = parameters.GetParameter< int >( "x-size" );
		if( meshes.x() <= 0 )
		{
		cerr << "Error: x-size must be positive integer number! It is " << meshes. x() << " now." << endl;
		return false;
		}
		mesh.setDimensions( meshes. x() );
		return true;
		}

		template< typename Mesh >
		template< typename Real, typename Device, typename Index,
		template< int, typename, typename, typename > class Geometry >
		bool simpleProblemSolver< Mesh >::initMesh( tnlGrid< 2, Real, Device, Index, Geometry >& mesh,
		const tnlParameterContainer& parameters ) const
		{
		tnlTuple< 2, IndexType > meshes;
		meshes.x() = parameters.GetParameter< int >( "x-size" );
		meshes.y() = parameters.GetParameter< int >( "y-size" );
		if( meshes.x() <= 0 )
		{
		cerr << "Error: x-size must be positive integer number! It is " << meshes. x() << " now." << endl;
		return false;
		}
		if( meshes.y() <= 0 )
		{
		cerr << "Error: y-size must be positive integer number! It is " << meshes. y() << " now." << endl;
		return false;
		}
		mesh.setDimensions( meshes. x(), meshes. y() );
		return true;
		}

		template< typename Mesh >
		template< typename Real, typename Device, typename Index,
		template< int, typename, typename, typename > class Geometry >
		bool simpleProblemSolver< Mesh >::initMesh( tnlGrid< 3, Real, Device, Index, Geometry >& mesh,
		const tnlParameterContainer& parameters ) const
		{
		tnlTuple< 3, IndexType > meshes;
		meshes.x() = parameters.GetParameter< int >( "x-size" );
		meshes.y() = parameters.GetParameter< int >( "y-size" );
		meshes.z() = parameters.GetParameter< int >( "z-size" );
		if( meshes.x() <= 0 )
		{
		cerr << "Error: x-size must be positive integer number! It is " << meshes. x() << " now." << endl;
		return false;
		}
		if( meshes.y() <= 0 )
		{
		cerr << "Error: y-size must be positive integer number! It is " << meshes. y() << " now." << endl;
		return false;
		}
		if( meshes.z() <= 0 )
		{
		cerr << "Error: z-size must be positive integer number! It is " << meshes. z() << " now." << endl;
		return false;
		}

		mesh.setDimensions( meshes.x(), meshes.y(). meshses.z() );
		return true;
		}

		template< typename Mesh >
		bool simpleProblemSolver< Mesh> :: init( const tnlParameterContainer& parameters )
		{
		@@ -70,25 +140,27 @@ bool simpleProblemSolver< Mesh> :: init( const tnlParameterContainer& parameters

		/****
		* 2. Set-up geometry of the problem domain using some mesh like tnlGrid.
		* Implement additional template specializations of the method initMesh
		* if necessary.
		*/
		if( ! this->initMesh( this->mesh, parameters ) )
		return false;
		if( ! this->mesh.save( "mesh.tnl" ) )
		{
		cerr << "I am not able to save the mesh into a file mesh.tnl." << endl;
		return false;
		}

		/****
		* 3. Set-up DOFs and supporting grid functions
		*/
		dofVector. setSize( 100 );
		const IndexType& dofs = this->mesh.getDofs();
		dofVector. setSize( 2*dofs );

		/****
		* You may use tnlSharedVector if you need to split the dofVector into more
		* grid functions like the following example:
		*/

		const IndexType& dofs = 50;
		/****
		* Usually you will replace this by:
		*
		* const IndexType& dofs = mesh. getDofs();
		*
		*/

		this -> u. bind( & dofVector. getData()[ 0 * dofs ], dofs );
		this -> v. bind( & dofVector. getData()[ 1 * dofs ], dofs );
		/****

src/implementation/solvers/cfd/navier-stokes/tnlNavierStokesSolver_impl.h

+59 −10

Original line number	Diff line number	Diff line
		@@ -27,7 +27,7 @@ tnlNavierStokesSolver< AdvectionScheme, DiffusionScheme, BoundaryConditions >::t
		: advection( 0 ),
		u1Viscosity( 0 ),
		u2Viscosity( 0 ),
		eViscosity( 0 ),
		temperatureViscosity( 0 ),
		mu( 0.0 ),
		gravity( 0.0 ),
		R( 0.0 ),
		@@ -63,11 +63,11 @@ template< typename AdvectionScheme,
		typename BoundaryConditions >
		void tnlNavierStokesSolver< AdvectionScheme, DiffusionScheme, BoundaryConditions >::setDiffusionScheme( DiffusionSchemeType& u1Viscosity,
		DiffusionSchemeType& u2Viscosity,
		DiffusionSchemeType& eViscosity )
		DiffusionSchemeType& temperatureViscosity )
		{
		this->u1Viscosity = &u1Viscosity;
		this->u2Viscosity = &u2Viscosity;
		this->eViscosity = &eViscosity;
		this->temperatureViscosity = &temperatureViscosity;
		}

		template< typename AdvectionScheme,
		@@ -89,6 +89,7 @@ void tnlNavierStokesSolver< AdvectionScheme, DiffusionScheme, BoundaryConditions
		this->u2.setSize( this->mesh->getDofs() );
		this->p.setSize( this->mesh->getDofs() );
		this->temperature.setSize( this->mesh->getDofs() );
		this->rhsDofVector.setSize( this->getDofs() );
		}

		template< typename AdvectionScheme,
		@@ -353,7 +354,7 @@ void tnlNavierStokesSolver< AdvectionScheme,
		this->advection->setRhoU2( dofs_rho_u2 );
		this->advection->setE( dofs_e );
		this->advection->setP( this->p );
		this->eViscosity->setFunction( dofs_e );
		this->temperatureViscosity->setFunction( this->temperature );

		rho_t.bind( & fu. getData()[ 0 ], dofs );
		rho_u1_t.bind( & fu. getData()[ dofs ], dofs );
		@@ -451,15 +452,32 @@ void tnlNavierStokesSolver< AdvectionScheme,
		/***
		* Add the viscosity term
		*/
		rho_u1_t[ c ] += this->mu * u1Viscosity->getDiffusion( c );
		rho_u1_t[ c ] += this->mu*u1Viscosity->getDiffusion( c );;
		rho_u2_t[ c ] += this->mu*u2Viscosity->getDiffusion( c );
		e_t[ c ] += eViscosity->getDiffusion( c );

		IndexType e = this->mesh->getElementIndex( i+1, j );
		IndexType w = this->mesh->getElementIndex( i-1, j );
		IndexType n = this->mesh->getElementIndex( i, j+1 );
		IndexType s = this->mesh->getElementIndex( i, j-1 );
		RealType u1_e = 0.5*( this->u1[ c ] + this->u1[ e ] );
		RealType u1_w = 0.5*( this->u1[ c ] + this->u1[ w ] );
		RealType u2_n = 0.5*( this->u2[ c ] + this->u2[ n ] );
		RealType u2_s = 0.5*( this->u2[ c ] + this->u2[ s ] );
		RealType hx = mesh->getParametricStep().x();
		RealType hy = mesh->getParametricStep().y();
		RealType u1_x_e = ( this->u1[ e ] - this->u1[ c ] ) / hx;
		RealType u1_x_w = ( this->u1[ c ] - this->u1[ w ] ) / hx;
		RealType u2_y_n = ( this->u2[ n ] - this->u2[ c ] ) / hy;
		RealType u2_y_s = ( this->u2[ c ] - this->u2[ s ] ) / hy;


		e_t[ c ] += this->mu( ( u1_e u1_x_e - u1_w * u1_x_w ) / hx +
		( u2_n * u2_y_n - u2_s * u2_y_s ) / hy );
		e_t[ c ] = 0.0;
		}

		/*rhsDofVector = fu;
		makeSnapshot( 0.0, 1 );
		getchar();*/
		writeExplicitRhs( time, -1 );
		getchar();
		}

		template< typename AdvectionScheme,
		@@ -548,5 +566,36 @@ typename tnlNavierStokesSolver< AdvectionScheme,
		0.5 * ( u1u1 + u2u2 ) );
		}

		template< typename AdvectionScheme,
		typename DiffusionScheme,
		typename BoundaryConditions >
		bool tnlNavierStokesSolver< AdvectionScheme,
		DiffusionScheme,
		BoundaryConditions >::writeExplicitRhs( const RealType& t,
		const IndexType step )
		{
		tnlSharedVector< RealType, DeviceType, IndexType > dofs_rho, dofs_rho_u1, dofs_rho_u2;

		const IndexType& dofs = mesh->getDofs();
		dofs_rho. bind( & this->rhsDofVector.getData()[ 0 ], dofs );
		dofs_rho_u1. bind( & this->rhsDofVector.getData()[ dofs ], dofs );
		dofs_rho_u2. bind( & this->rhsDofVector.getData()[ 2 * dofs ], dofs );

		tnlString fileName;
		FileNameBaseNumberEnding( "rho-t-", step, 5, ".tnl", fileName );
		if( ! dofs_rho. save( fileName ) )
		return false;

		FileNameBaseNumberEnding( "rho-u1-t-", step, 5, ".tnl", fileName );
		if( ! dofs_rho_u1. save( fileName ) )
		return false;

		FileNameBaseNumberEnding( "rho-u2-t-", step, 5, ".tnl", fileName );
		if( ! dofs_rho_u2. save( fileName ) )
		return false;

		return true;
		}


		#endif /* tnlNavierStokesSolver_IMPL_H_ */