cuIBM
A GPU-based Immersed Boundary Method code
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
NavierStokesSolver< memoryType > Class Template Reference

Solves the Navier-Stokes equations in a rectangular domain. More...

#include <NavierStokesSolver.h>

Inheritance diagram for NavierStokesSolver< memoryType >:
Inheritance graph
[legend]
Collaboration diagram for NavierStokesSolver< memoryType >:
Collaboration graph
[legend]

Public Member Functions

 NavierStokesSolver (parameterDB *pDB=NULL, domain *dInfo=NULL)
 Constructor. Copies the database and information about the computational grid. More...
 
virtual void initialise ()
 Initializes parameters, arrays and matrices required for the simulation. More...
 
void stepTime ()
 Calculates the variables at the next time step. More...
 
virtual void writeCommon ()
 Writes numerical solution at current time-step, as well as the number of iterations performed in each solver. More...
 
virtual void writeData ()
 Writes data into files. More...
 
bool finished ()
 Evaluates the condition required to stop the simulation. More...
 
virtual void shutDown ()
 Prints timing information and closes the different files. More...
 
virtual std::string name ()
 Returns the name of the current solver. More...
 

Protected Member Functions

void initialiseCommon ()
 Initializes parameters common to all Navier-Stokes solvers. More...
 
void initialiseArrays (int numQ, int numLambda)
 Initializes all arrays required to solve the Navier-Stokes equations. More...
 
virtual void initialiseFluxes ()
 
virtual void initialiseFluxes (real *q)
 Initializes velocity flux vectors. More...
 
void initialiseBoundaryArrays ()
 Initializes boundary velocity arrays with values stored in the database. More...
 
void assembleMatrices ()
 Assembles matrices of the intermediate flux solver and the Poisson solver. More...
 
void generateM ()
 
virtual void generateL ()
 
virtual void generateA (real alpha)
 
void generateBN ()
 Generates approximate inverse of the matrix resulting from implicit velocity terms. More...
 
virtual void generateQT ()
 
void updateQ (real gamma)
 Doing nothing. More...
 
virtual void generateC ()
 
void generateRN ()
 Generates explicit terms of the momentum equation. More...
 
void calculateExplicitQTerms ()
 
virtual void calculateExplicitLambdaTerms ()
 Doing nothing. Used in methods that use the explicit pressure term in the intermediate velocity solve step, such as FadlunEtAlSolver and DFModifiedSolver. More...
 
virtual void generateBC1 ()
 
virtual void generateBC2 ()
 
virtual void assembleRHS1 ()
 Assembles the right hand-side of the system for the intermediate flux. More...
 
void assembleRHS2 ()
 Assembles the right hand-side of the Poisson system. More...
 
virtual void solveIntermediateVelocity ()
 Solves for the intermediate flux velocity. More...
 
virtual void solvePoisson ()
 Solves the Poisson system. More...
 
virtual void projectionStep ()
 Projects the flux onto the divergence-free field. More...
 
void updateBoundaryConditions ()
 Doing nothing. More...
 
virtual void updateSolverState ()
 Doing nothing. Used in immersed boundary methods when the body moves. More...
 
template<>
void generateA (real alpha)
 Assembles the matrix from implicit terms of momentum equation. More...
 
template<>
void generateBC1 ()
 Generates inhomogeneous boundary conditions from the discrete Laplacian operator. More...
 
template<>
void generateBC2 ()
 Generates inhomogeneous boundary conditions from the discrete continuity equation. More...
 
template<>
void generateL ()
 Generates the discrete Laplacian matrix. More...
 
template<>
void generateM ()
 Generates the mass matrix and its inverse. More...
 
template<>
void generateQT ()
 Generates the discrete divergence matrix. More...
 
template<>
void calculateExplicitQTerms ()
 Generates explicit terms that arise from the velocity fluxes. More...
 
template<>
void initialiseFluxes ()
 Initializes velocity flux vectors (on the device). More...
 
template<>
void generateC ()
 Generates the matrix of the Poisson solver. More...
 

Protected Attributes

parameterDBparamDB
 database that contains all the simulation parameters More...
 
domaindomInfo
 computational grid information More...
 
integrationScheme intgSchm
 instance of the class integrationScheme More...
 
real QCoeff
 
cusp::coo_matrix< int, real, memoryType > M
 diagonal mass matrix More...
 
cusp::coo_matrix< int, real, memoryType > Minv
 inverse of the mass matrix More...
 
cusp::coo_matrix< int, real, memoryType > L
 discrete Laplacian matrix More...
 
cusp::coo_matrix< int, real, memoryType > A
 combination of mass and Laplacian matrices More...
 
cusp::coo_matrix< int, real, memoryType > Q
 gradient matrix (and regularization matrix if immersed body) More...
 
cusp::coo_matrix< int, real, memoryType > QT
 transposed gradient matrix (and interpolation matrix if immersed body) More...
 
cusp::coo_matrix< int, real, memoryType > BN
 N-th order Taylor series expansion of the inverse of A. More...
 
cusp::coo_matrix< int, real, memoryType > C
 matrix of the Poisson system More...
 
cusp::array1d< real, memoryType > q
 velocity flux vector More...
 
cusp::array1d< real, memoryType > qStar
 intermediate velocity flux vector More...
 
cusp::array1d< real, memoryType > qOld
 velocity flux vector at the previous time-step More...
 
cusp::array1d< real, memoryType > lambda
 pressure vector (and body forces if immersed body) More...
 
cusp::array1d< real, memoryType > rn
 explicit terms of the momentum equation More...
 
cusp::array1d< real, memoryType > bc1
 inhomogeneous boundary conditions from the discrete Laplacian operator More...
 
cusp::array1d< real, memoryType > bc2
 inhomogeneous boundary conditions from the discrete continuity equation More...
 
cusp::array1d< real, memoryType > rhs1
 right hand-side of the system for the intermediate velocity flux vector More...
 
cusp::array1d< real, memoryType > rhs2
 right hand-side of the Poisson system More...
 
cusp::array1d< real, memoryType > H
 hold convective terms from previous time step More...
 
cusp::array1d< real, memoryType > temp1
 temporary 1D Cusp array More...
 
cusp::array1d< real, memoryType > temp2
 temporary 1D Cusp array More...
 
cusp::array1d< real, memoryType > bc [4]
 array that contains the boundary conditions of the rectangular domain More...
 
preconditioner< cusp::coo_matrix< int, real, memoryType > > * PC1
 preconditioner for the intermediate flux solver More...
 
preconditioner< cusp::coo_matrix< int, real, memoryType > > * PC2
 preconditioner for the Poisson solver More...
 
size_t timeStep
 iteration number More...
 
size_t subStep
 number of sub-iterations More...
 
size_t iterationCount1
 number of iteration to solve the intermediate fluxes More...
 
size_t iterationCount2
 number of iteration to solve the Poisson equation More...
 
Logger logger
 instance of the class Logger to track time of different tasks More...
 
std::ofstream iterationsFile
 file that contains the number of iterations More...
 

Detailed Description

template<typename memoryType>
class NavierStokesSolver< memoryType >

Solves the Navier-Stokes equations in a rectangular domain.

Methods are defined as virtual when they are redefined in a derived class with the same name.

Definition at line 26 of file NavierStokesSolver.h.

Constructor & Destructor Documentation

template<typename memoryType >
NavierStokesSolver< memoryType >::NavierStokesSolver ( parameterDB pDB = NULL,
domain dInfo = NULL 
)

Constructor. Copies the database and information about the computational grid.

Parameters
pDBdatabase that contains all the simulation parameters
dInfoinformation related to the computational grid

Definition at line 24 of file NavierStokesSolver.cu.

Member Function Documentation

template<typename memoryType >
void NavierStokesSolver< memoryType >::assembleMatrices ( )
protected

Assembles matrices of the intermediate flux solver and the Poisson solver.

Definition at line 338 of file NavierStokesSolver.cu.

References kernels::generateA(), and kernels::generateQT().

Referenced by TairaColoniusSolver< memoryType >::initialise(), DirectForcingSolver< memoryType >::initialise(), and DirectForcingSolver< memoryType >::updateSolverState().

template<typename memoryType >
void NavierStokesSolver< memoryType >::assembleRHS1 ( )
protectedvirtual

Assembles the right hand-side of the system for the intermediate flux.

Reimplemented in DirectForcingSolver< memoryType >.

Definition at line 413 of file NavierStokesSolver.cu.

Referenced by DirectForcingSolver< memoryType >::assembleRHS1().

template<typename memoryType >
void NavierStokesSolver< memoryType >::assembleRHS2 ( )
protected

Assembles the right hand-side of the Poisson system.

Definition at line 427 of file NavierStokesSolver.cu.

template<typename memoryType >
void NavierStokesSolver< memoryType >::calculateExplicitLambdaTerms ( )
protectedvirtual

Doing nothing. Used in methods that use the explicit pressure term in the intermediate velocity solve step, such as FadlunEtAlSolver and DFModifiedSolver.

Reimplemented in FadlunEtAlSolver< memoryType >, and DFModifiedSolver< memoryType >.

Definition at line 19 of file generateRN.inl.

template<>
void NavierStokesSolver< device_memory >::calculateExplicitQTerms ( )
protected

Generates explicit terms that arise from the velocity fluxes.

Includes the time derivative, convection and diffusion terms.

Definition at line 30 of file generateRN.inl.

References BSZ, XMINUS, XPLUS, YMINUS, and YPLUS.

template<typename memoryType>
void NavierStokesSolver< memoryType >::calculateExplicitQTerms ( )
protected
template<typename memoryType >
bool NavierStokesSolver< memoryType >::finished ( )

Evaluates the condition required to stop the simulation.

Returns
a Boolean to continue or stop the simulation

Definition at line 322 of file NavierStokesSolver.cu.

Referenced by main().

template<>
void NavierStokesSolver< device_memory >::generateA ( real  alpha)
protected

Assembles the matrix from implicit terms of momentum equation.

Generates the matrix A which solves the implicit velocity flux. Calculated from the matrices M (which contains the term from the time derivative) and the matrix L, which is the implicit diffusion matrix.

Parameters
alphaimplicit coefficient of the diffusive scheme

Definition at line 22 of file generateA.inl.

template<typename memoryType>
virtual void NavierStokesSolver< memoryType >::generateA ( real  alpha)
protectedvirtual
template<>
void NavierStokesSolver< device_memory >::generateBC1 ( )
protected

Generates inhomogeneous boundary conditions from the discrete Laplacian operator.

bottom

u

multiply by 0.5 for the Crank-Nicolson scheme and 2.0 for the non-uniform central difference

v

top

u

v

u

v

u

v

left

u

v

right

u

v

u

v

Definition at line 16 of file generateBC1.inl.

References CONVECTIVE, DIRICHLET, SPECIAL, XMINUS, XPLUS, YMINUS, and YPLUS.

template<typename memoryType>
virtual void NavierStokesSolver< memoryType >::generateBC1 ( )
protectedvirtual
template<>
void NavierStokesSolver< device_memory >::generateBC2 ( )
protected

Generates inhomogeneous boundary conditions from the discrete continuity equation.

Definition at line 16 of file generateBC2.inl.

References XMINUS, XPLUS, YMINUS, and YPLUS.

template<typename memoryType>
virtual void NavierStokesSolver< memoryType >::generateBC2 ( )
protectedvirtual
template<typename memoryType >
void NavierStokesSolver< memoryType >::generateBN ( )
protected

Generates approximate inverse of the matrix resulting from implicit velocity terms.

It computes the N-th order Taylor expansion of the inverse matrix. Currently, the order is N=1.

Definition at line 368 of file NavierStokesSolver.cu.

template<typename memoryType>
virtual void NavierStokesSolver< memoryType >::generateC ( )
protectedvirtual
template<>
void NavierStokesSolver< device_memory >::generateC ( )
protected

Generates the matrix of the Poisson solver.

Definition at line 392 of file NavierStokesSolver.cu.

template<>
void NavierStokesSolver< device_memory >::generateL ( )
protected

Generates the discrete Laplacian matrix.

The kinematic viscosity is included in the matrix. The matrix is stored in a COO format.

x-component

< calculate the row of the matrix

< scaling factor (to obtain the normalised matrix)

< no south coefficient for the bottom row

< no west coefficient for the leftmost column

y-component

Definition at line 15 of file generateL.inl.

template<typename memoryType>
virtual void NavierStokesSolver< memoryType >::generateL ( )
protectedvirtual
template<>
void NavierStokesSolver< device_memory >::generateM ( )
protected

Generates the mass matrix and its inverse.

The mass diagonal matrix is stored in a COO format. The time-increment value is included in the matrix.

Definition at line 18 of file generateM.inl.

template<typename memoryType>
void NavierStokesSolver< memoryType >::generateM ( )
protected
template<>
void NavierStokesSolver< device_memory >::generateQT ( )
protected

Generates the discrete divergence matrix.

The gradient matrix is also computed as the transpose of the divergence one.

Definition at line 16 of file generateQT.inl.

References kernels::generateQT().

template<typename memoryType>
virtual void NavierStokesSolver< memoryType >::generateQT ( )
protectedvirtual
template<typename memoryType >
void NavierStokesSolver< memoryType >::generateRN ( )
protected

Generates explicit terms of the momentum equation.

Definition at line 89 of file generateRN.inl.

template<typename memoryType >
void NavierStokesSolver< memoryType >::initialise ( )
virtual

Initializes parameters, arrays and matrices required for the simulation.

Reimplemented in DirectForcingSolver< memoryType >, TairaColoniusSolver< memoryType >, and DiffusionSolver< memoryType >.

Definition at line 35 of file NavierStokesSolver.cu.

Referenced by NSWithBody< memoryType >::initialiseBodies(), and main().

template<typename memoryType >
void NavierStokesSolver< memoryType >::initialiseArrays ( int  numQ,
int  numLambda 
)
protected

Initializes all arrays required to solve the Navier-Stokes equations.

Parameters
numQtotal number velocity (or flux) unknowns (x- and y- directions)
numLambdanumber of pressure unknowns (plus number of body force unknowns)

Definition at line 98 of file NavierStokesSolver.cu.

Referenced by TairaColoniusSolver< memoryType >::initialise(), and DirectForcingSolver< memoryType >::initialise().

template<typename memoryType >
void NavierStokesSolver< memoryType >::initialiseBoundaryArrays ( )
protected

Initializes boundary velocity arrays with values stored in the database.

Top and Bottom

Left and Right

Definition at line 211 of file NavierStokesSolver.cu.

References boundaryCondition::value, XMINUS, XPLUS, YMINUS, and YPLUS.

template<typename memoryType >
void NavierStokesSolver< memoryType >::initialiseCommon ( )
protected

Initializes parameters common to all Navier-Stokes solvers.

Definition at line 55 of file NavierStokesSolver.cu.

References io::writeGrid().

Referenced by TairaColoniusSolver< memoryType >::initialise(), and DirectForcingSolver< memoryType >::initialise().

template<typename memoryType>
virtual void NavierStokesSolver< memoryType >::initialiseFluxes ( )
protectedvirtual
template<typename memoryType >
void NavierStokesSolver< memoryType >::initialiseFluxes ( real q)
protectedvirtual

Initializes velocity flux vectors.

Parameters
qthe velocity flux vector

Definition at line 165 of file NavierStokesSolver.cu.

References io::readData().

template<>
void NavierStokesSolver< device_memory >::initialiseFluxes ( )
protected

Initializes velocity flux vectors (on the device).

It creates a raw pointer before calling a method to initialize the flux vector.

Definition at line 145 of file NavierStokesSolver.cu.

template<typename memoryType>
virtual std::string NavierStokesSolver< memoryType >::name ( )
inlinevirtual
template<typename memoryType >
void NavierStokesSolver< memoryType >::projectionStep ( )
protectedvirtual
template<typename memoryType >
void NavierStokesSolver< memoryType >::shutDown ( )
virtual

Prints timing information and closes the different files.

Reimplemented in NSWithBody< memoryType >.

Definition at line 591 of file NavierStokesSolver.cu.

References io::printTimingInfo().

Referenced by main().

template<typename memoryType >
void NavierStokesSolver< memoryType >::solveIntermediateVelocity ( )
protectedvirtual

Solves for the intermediate flux velocity.

Definition at line 446 of file NavierStokesSolver.cu.

template<typename memoryType >
void NavierStokesSolver< memoryType >::solvePoisson ( )
protectedvirtual

Solves the Poisson system.

Reimplemented in DiffusionSolver< memoryType >.

Definition at line 491 of file NavierStokesSolver.cu.

template<typename memoryType >
void NavierStokesSolver< memoryType >::stepTime ( )

Calculates the variables at the next time step.

Definition at line 256 of file NavierStokesSolver.cu.

Referenced by main().

template<typename memoryType >
void NavierStokesSolver< memoryType >::updateBoundaryConditions ( )
protected

Doing nothing.

Definition at line 311 of file NavierStokesSolver.cu.

template<typename memoryType >
void NavierStokesSolver< memoryType >::updateQ ( real  gamma)
protected

Doing nothing.

Parameters
gammacoefficient of the convection term at the current time step

Definition at line 300 of file NavierStokesSolver.cu.

template<typename memoryType >
void NavierStokesSolver< memoryType >::updateSolverState ( )
protectedvirtual

Doing nothing. Used in immersed boundary methods when the body moves.

Reimplemented in DirectForcingSolver< memoryType >, and TairaColoniusSolver< memoryType >.

Definition at line 286 of file NavierStokesSolver.cu.

template<typename memoryType >
void NavierStokesSolver< memoryType >::writeCommon ( )
virtual

Writes numerical solution at current time-step, as well as the number of iterations performed in each solver.

Reimplemented in NSWithBody< memoryType >.

Definition at line 557 of file NavierStokesSolver.cu.

References io::writeData().

Referenced by NSWithBody< memoryType >::writeCommon().

template<typename memoryType >
void NavierStokesSolver< memoryType >::writeData ( )
virtual

Writes data into files.

Reimplemented in DirectForcingSolver< memoryType >, and TairaColoniusSolver< memoryType >.

Definition at line 577 of file NavierStokesSolver.cu.

Referenced by main().

Member Data Documentation

template<typename memoryType>
cusp::coo_matrix<int, real, memoryType> NavierStokesSolver< memoryType >::A
protected

combination of mass and Laplacian matrices

Definition at line 36 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::array1d<real, memoryType> NavierStokesSolver< memoryType >::bc[4]
protected

array that contains the boundary conditions of the rectangular domain

Definition at line 46 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::array1d<real, memoryType> NavierStokesSolver< memoryType >::bc1
protected

inhomogeneous boundary conditions from the discrete Laplacian operator

Definition at line 46 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::array1d<real, memoryType> NavierStokesSolver< memoryType >::bc2
protected

inhomogeneous boundary conditions from the discrete continuity equation

Definition at line 46 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::coo_matrix<int, real, memoryType> NavierStokesSolver< memoryType >::BN
protected

N-th order Taylor series expansion of the inverse of A.

Definition at line 36 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::coo_matrix<int, real, memoryType> NavierStokesSolver< memoryType >::C
protected

matrix of the Poisson system

Definition at line 36 of file NavierStokesSolver.h.

template<typename memoryType>
domain* NavierStokesSolver< memoryType >::domInfo
protected

computational grid information

Definition at line 30 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::array1d<real, memoryType> NavierStokesSolver< memoryType >::H
protected

hold convective terms from previous time step

Definition at line 46 of file NavierStokesSolver.h.

template<typename memoryType>
integrationScheme NavierStokesSolver< memoryType >::intgSchm
protected

instance of the class integrationScheme

Definition at line 31 of file NavierStokesSolver.h.

template<typename memoryType>
size_t NavierStokesSolver< memoryType >::iterationCount1
protected

number of iteration to solve the intermediate fluxes

Definition at line 64 of file NavierStokesSolver.h.

template<typename memoryType>
size_t NavierStokesSolver< memoryType >::iterationCount2
protected

number of iteration to solve the Poisson equation

Definition at line 64 of file NavierStokesSolver.h.

template<typename memoryType>
std::ofstream NavierStokesSolver< memoryType >::iterationsFile
protected

file that contains the number of iterations

Definition at line 71 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::coo_matrix<int, real, memoryType> NavierStokesSolver< memoryType >::L
protected

discrete Laplacian matrix

Definition at line 36 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::array1d<real, memoryType> NavierStokesSolver< memoryType >::lambda
protected

pressure vector (and body forces if immersed body)

Definition at line 46 of file NavierStokesSolver.h.

template<typename memoryType>
Logger NavierStokesSolver< memoryType >::logger
protected

instance of the class Logger to track time of different tasks

Definition at line 69 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::coo_matrix<int, real, memoryType> NavierStokesSolver< memoryType >::M
protected

diagonal mass matrix

Definition at line 36 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::coo_matrix<int, real, memoryType> NavierStokesSolver< memoryType >::Minv
protected

inverse of the mass matrix

Definition at line 36 of file NavierStokesSolver.h.

template<typename memoryType>
parameterDB* NavierStokesSolver< memoryType >::paramDB
protected

database that contains all the simulation parameters

Definition at line 29 of file NavierStokesSolver.h.

template<typename memoryType>
preconditioner< cusp::coo_matrix<int, real, memoryType> >* NavierStokesSolver< memoryType >::PC1
protected

preconditioner for the intermediate flux solver

Definition at line 61 of file NavierStokesSolver.h.

template<typename memoryType>
preconditioner< cusp::coo_matrix<int, real, memoryType> > * NavierStokesSolver< memoryType >::PC2
protected

preconditioner for the Poisson solver

Definition at line 61 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::coo_matrix<int, real, memoryType> NavierStokesSolver< memoryType >::Q
protected

gradient matrix (and regularization matrix if immersed body)

Definition at line 36 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::array1d<real, memoryType> NavierStokesSolver< memoryType >::q
protected

velocity flux vector

Definition at line 46 of file NavierStokesSolver.h.

template<typename memoryType>
real NavierStokesSolver< memoryType >::QCoeff
protected

Definition at line 33 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::array1d<real, memoryType> NavierStokesSolver< memoryType >::qOld
protected

velocity flux vector at the previous time-step

Definition at line 46 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::array1d<real, memoryType> NavierStokesSolver< memoryType >::qStar
protected

intermediate velocity flux vector

Definition at line 46 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::coo_matrix<int, real, memoryType> NavierStokesSolver< memoryType >::QT
protected

transposed gradient matrix (and interpolation matrix if immersed body)

Definition at line 36 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::array1d<real, memoryType> NavierStokesSolver< memoryType >::rhs1
protected

right hand-side of the system for the intermediate velocity flux vector

Definition at line 46 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::array1d<real, memoryType> NavierStokesSolver< memoryType >::rhs2
protected

right hand-side of the Poisson system

Definition at line 46 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::array1d<real, memoryType> NavierStokesSolver< memoryType >::rn
protected

explicit terms of the momentum equation

Definition at line 46 of file NavierStokesSolver.h.

template<typename memoryType>
size_t NavierStokesSolver< memoryType >::subStep
protected

number of sub-iterations

Definition at line 64 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::array1d<real, memoryType> NavierStokesSolver< memoryType >::temp1
protected

temporary 1D Cusp array

Definition at line 46 of file NavierStokesSolver.h.

template<typename memoryType>
cusp::array1d<real, memoryType> NavierStokesSolver< memoryType >::temp2
protected

temporary 1D Cusp array

Definition at line 46 of file NavierStokesSolver.h.

template<typename memoryType>
size_t NavierStokesSolver< memoryType >::timeStep
protected

iteration number

Definition at line 64 of file NavierStokesSolver.h.


The documentation for this class was generated from the following files: