cuIBM
A GPU-based Immersed Boundary Method code
|
A fully discrete formulation of the direct forcing method. More...
#include <DirectForcingSolver.h>
Public Member Functions | |
DirectForcingSolver (parameterDB *pDB=NULL, domain *dInfo=NULL) | |
Constructor. Initializes the simulation parameters and the domain info. More... | |
virtual void | initialise () |
Initialize the vectors used in the simulation. More... | |
virtual void | writeData () |
Writes the velocity, pressure, force and mass flux data at every save point. More... | |
virtual std::string | name () |
Returns the name of the solver as a string. More... | |
Public Member Functions inherited from NSWithBody< memoryType > | |
virtual void | writeCommon () |
Writes flow variables and position of body points into files. More... | |
virtual void | shutDown () |
Closes iteration file and force file. More... | |
Public Member Functions inherited from NavierStokesSolver< memoryType > | |
NavierStokesSolver (parameterDB *pDB=NULL, domain *dInfo=NULL) | |
Constructor. Copies the database and information about the computational grid. More... | |
void | stepTime () |
Calculates the variables at the next time step. More... | |
bool | finished () |
Evaluates the condition required to stop the simulation. More... | |
Protected Member Functions | |
void | tagPoints () |
void | tagPoints (real *bx, real *by, real *uB, real *vB) |
Tags all the forcing nodes required for the type of linear interpolation explained in the paper by Fadlun et al. (2000). More... | |
virtual void | generateA (real alpha) |
virtual void | generateL () |
void | generateQT (int *QTRows, int *QTCols, real *QTVals) |
virtual void | generateQT () |
Compute the divergence operator (transpose of the gradient one). More... | |
void | updateQ () |
virtual void | generateC () |
Generates the right-hand side matrix in the Poisson step. More... | |
virtual void | assembleRHS1 () |
Assembles the matrix rhs1 for DirectForcingSolver. More... | |
void | updateRHS1 () |
virtual void | updateSolverState () |
Updates the matrices every time the body is moved. More... | |
virtual void | projectionStep () |
Projects the pressure gradient on to the intermediate velocity field to obtain the divergence-free velocity field at the next time step. More... | |
void | writeMassFluxInfo () |
Prints the min, max and sum of the divergences of the velocity field in every cell of the domain. More... | |
template<> | |
void | generateA (real alpha) |
Generates the matrix A on the device. More... | |
template<> | |
void | generateL () |
Sets up the Laplacian matrix for the implicit diffusion term. Uses second order central differences on the non-uniform grid. More... | |
template<> | |
void | updateQ () |
Update the gradient operator (device). More... | |
template<> | |
void | tagPoints () |
Tags the forcing nodes among the velocity nodes, i.e. the nodes at which the velocity interpolation is performed. More... | |
template<> | |
void | updateRHS1 () |
Update the RHS of the velocity system. More... | |
Protected Member Functions inherited from NSWithBody< memoryType > | |
virtual void | calculateForce () |
void | initialiseBodies () |
Stores the parameters of the simulation and initializes the location and motion of each immersed bodies. More... | |
void | updateBodies () |
Updates location and motion of each immersed body at current time. More... | |
template<> | |
void | calculateForce () |
Calculates forces acting on an immersed body. More... | |
Protected Member Functions inherited from NavierStokesSolver< memoryType > | |
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 () |
void | generateBN () |
Generates approximate inverse of the matrix resulting from implicit velocity terms. More... | |
void | updateQ (real gamma) |
Doing nothing. More... | |
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 () |
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... | |
void | updateBoundaryConditions () |
Doing nothing. 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 | |
cusp::array1d< int, host_memory > | tags |
vector used to differentiate forcing points from regular ones (host) More... | |
cusp::array1d< int, host_memory > | tags2 |
vector used to differentiate forcing points from regular ones (host) More... | |
cusp::array1d< int, device_memory > | tagsD |
vector used to differentiate forcing points from regular ones (device) More... | |
cusp::array1d< int, device_memory > | tags2D |
vector used to differentiate forcing points from regular ones (device) More... | |
vecH | coeffs |
coefficients of interpolation (host) More... | |
vecH | coeffs2 |
other coefficients of interpolation; quadratic interpolation (host) More... | |
vecD | coeffsD |
coefficients of interpolation (device) More... | |
vecD | coeffs2D |
other coefficients of interpolation; quadratic interpolation (device) More... | |
vecH | uv |
velocity field (host) More... | |
vecD | uvD |
velocity field (device) More... | |
cusp::array1d< real, memoryType > | pressure |
pressure field More... | |
Protected Attributes inherited from NSWithBody< memoryType > | |
bodies< memoryType > | B |
bodies in the flow More... | |
real | forceX |
force acting on each body in the x-direction More... | |
real | forceY |
force acting on each body in the y-direction More... | |
std::ofstream | forceFile |
file to write the forces More... | |
Protected Attributes inherited from NavierStokesSolver< memoryType > | |
parameterDB * | paramDB |
database that contains all the simulation parameters More... | |
domain * | domInfo |
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... | |
A fully discrete formulation of the direct forcing method.
Based on the method first proposed by Fadlun et al (2000) with modifications by Kim et al (2001). It does not follow the same equations that they used, but use a fractional step method starting with the discretized equations, based on the idea by Perot (1993).
This method does not use an explicit pressure term in the step where the intermediate velocity is calculated, and the pressure is directly obtained from the Poisson equation.
Definition at line 27 of file DirectForcingSolver.h.
DirectForcingSolver< memoryType >::DirectForcingSolver | ( | parameterDB * | pDB = NULL , |
domain * | dInfo = NULL |
||
) |
Constructor. Initializes the simulation parameters and the domain info.
Definition at line 18 of file DirectForcingSolver.cu.
Referenced by DirectForcingSolver< memoryType >::generateQT().
|
protectedvirtual |
Assembles the matrix rhs1 for DirectForcingSolver.
This function first calls the function assembleRHS1 from NavierStokesSolver. Then it called the function updateRHS1 to modify only the elements of the vector that correspond to the interpolation nodes.
Reimplemented from NavierStokesSolver< memoryType >.
Definition at line 99 of file DirectForcingSolver.cu.
References NavierStokesSolver< memoryType >::assembleRHS1(), and kernels::updateRHS1().
Referenced by DirectForcingSolver< memoryType >::generateQT().
|
protectedvirtual |
Generates the matrix A on the device.
alpha | Coefficient of the implicit part of the diffusion term. 1 for backward Euler, 0.5 for Crank-Nicolson and 0 for explicit Euler. |
Reimplemented from NavierStokesSolver< memoryType >.
Definition at line 18 of file generateA.inl.
|
protectedvirtual |
Reimplemented from NavierStokesSolver< memoryType >.
|
protectedvirtual |
Generates the right-hand side matrix in the Poisson step.
Because the fully discrete direct forcing method separates the domains inside and outside the immersed boundary, and Neumann boundary conditions are enforced at the immersed boundary, a point also needs to be fixed inside the immersed boundary, so that Poisson system does not have any zero eigenvalues.
In this function, phi at a point that corresponds to the center of the grid is fixed as zero. This is the cell with indices (nx/2, ny/2), and is the center in terms of the indices and not the physical location in space. Of course, this is invalid if the interior of the body does not include this point, and the solution obtained will be unphysical.
This needs to be done in a better way (i.e. by locating points that are inside the immersed boundaries, and fixing them to zero)
Reimplemented from NavierStokesSolver< memoryType >.
Reimplemented in FEAModifiedSolver< memoryType >, and FadlunEtAlSolver< memoryType >.
Definition at line 206 of file DirectForcingSolver.cu.
References NavierStokesSolver< memoryType >::generateC().
Referenced by DirectForcingSolver< memoryType >::generateQT().
|
protectedvirtual |
Sets up the Laplacian matrix for the implicit diffusion term. Uses second order central differences on the non-uniform grid.
DO NOT FORGET (IMPLEMENTED) A is generated as M-alpha*LHost But some values of L are the interpolation values of the direct forcing method These should not be multplied by alpha
The sign of the interpolation coeffs is the opposite of that in the toy python code since A is obtained by subtracting alpha*L from M
(NOT IMPLEMENTED YET) Also remember to take alpha into account in BC1 for moving bodies
x-component
y-component
Reimplemented from NavierStokesSolver< memoryType >.
Definition at line 13 of file generateL.inl.
|
protectedvirtual |
Reimplemented from NavierStokesSolver< memoryType >.
|
inlineprotected |
Definition at line 57 of file DirectForcingSolver.h.
References DirectForcingSolver< memoryType >::assembleRHS1(), DirectForcingSolver< memoryType >::DirectForcingSolver(), DirectForcingSolver< memoryType >::generateC(), DirectForcingSolver< memoryType >::generateQT(), DirectForcingSolver< memoryType >::initialise(), DirectForcingSolver< memoryType >::projectionStep(), DirectForcingSolver< memoryType >::updateQ(), DirectForcingSolver< memoryType >::updateRHS1(), DirectForcingSolver< memoryType >::updateSolverState(), DirectForcingSolver< memoryType >::writeData(), and DirectForcingSolver< memoryType >::writeMassFluxInfo().
|
protectedvirtual |
Compute the divergence operator (transpose of the gradient one).
Reimplemented from NavierStokesSolver< memoryType >.
Reimplemented in FEAModifiedSolver< memoryType >, FadlunEtAlSolver< memoryType >, and DFImprovedSolver< memoryType >.
Definition at line 46 of file generateQT.inl.
References NavierStokesSolver< memoryType >::generateQT(), and kernels::updateQ().
Referenced by DirectForcingSolver< memoryType >::generateQT().
|
virtual |
Initialize the vectors used in the simulation.
Reimplemented from NavierStokesSolver< memoryType >.
Reimplemented in DiffusionSolver< memoryType >.
Definition at line 29 of file DirectForcingSolver.cu.
References NavierStokesSolver< memoryType >::assembleMatrices(), NavierStokesSolver< memoryType >::initialiseArrays(), NSWithBody< memoryType >::initialiseBodies(), and NavierStokesSolver< memoryType >::initialiseCommon().
Referenced by DirectForcingSolver< memoryType >::generateQT(), and DiffusionSolver< memoryType >::initialise().
|
inlinevirtual |
Returns the name of the solver as a string.
Reimplemented from NavierStokesSolver< memoryType >.
Reimplemented in FEAModifiedSolver< memoryType >, FadlunEtAlSolver< memoryType >, DiffusionSolver< memoryType >, DFModifiedSolver< memoryType >, and DFImprovedSolver< memoryType >.
Definition at line 93 of file DirectForcingSolver.h.
|
protectedvirtual |
Projects the pressure gradient on to the intermediate velocity field to obtain the divergence-free velocity field at the next time step.
Reimplemented from NavierStokesSolver< memoryType >.
Reimplemented in DFModifiedSolver< memoryType >, and DiffusionSolver< memoryType >.
Definition at line 153 of file DirectForcingSolver.cu.
References NavierStokesSolver< memoryType >::projectionStep().
Referenced by DirectForcingSolver< memoryType >::generateQT().
|
protected |
Tags the forcing nodes among the velocity nodes, i.e. the nodes at which the velocity interpolation is performed.
Definition at line 13 of file tagPoints.inl.
|
protected |
|
protected |
Tags all the forcing nodes required for the type of linear interpolation explained in the paper by Fadlun et al. (2000).
It uses a raytracing algorithm to detect points that are near the boundary, and just outside it. For more information about the algorithm, read the section on ray-crossings in the Search and Intersection chapter of the book Computational Geometry in C by Joseph O'Rourke.
bx | host array of the x-coordinates of the boundary points |
by | host array of the y-coordinates of the boundary points |
uB | host array of the x-components of the boundary velocities |
vB | host array of the y-components of the boundary velocities |
if the ray intersects the boundary segment top endpoint must be strictly above the ray bottom can be on or below the ray
Definition at line 60 of file tagPoints.inl.
|
protected |
Update the gradient operator (device).
After the Q matrix has been set up for the entire grid, this function removes the non-zeros that correspond to the velocity nodes where the interpolation is performed, i.e. it sets the fluxes there to zero. This routine is specific to DirectForcingSolver.
Definition at line 20 of file generateQT.inl.
|
protected |
|
protected |
Update the RHS of the velocity system.
The vector rhs1 is first set up as if it would have been in the absence of an immersed boundary. This function then changes values in the vector that correspond to the forcing nodes on the grid, replacing them with the rhs values from the interpolation relations.
Definition at line 20 of file updateRHS1.inl.
|
protected |
Referenced by DirectForcingSolver< memoryType >::generateQT().
|
protectedvirtual |
Updates the matrices every time the body is moved.
Change the location of the body points, re-tags all the points on the velocity grid to locate the new forcing nodes. Reassembles the matrices. Moving bodies have not been tested for DirectForcingSolver.
Reimplemented from NavierStokesSolver< memoryType >.
Definition at line 75 of file DirectForcingSolver.cu.
References NavierStokesSolver< memoryType >::assembleMatrices(), and NSWithBody< memoryType >::updateBodies().
Referenced by DirectForcingSolver< memoryType >::generateQT().
|
virtual |
Writes the velocity, pressure, force and mass flux data at every save point.
Reimplemented from NavierStokesSolver< memoryType >.
Definition at line 167 of file DirectForcingSolver.cu.
References NSWithBody< memoryType >::calculateForce(), kernels::forceX(), and NSWithBody< memoryType >::writeCommon().
Referenced by DirectForcingSolver< memoryType >::generateQT().
|
protected |
Prints the min, max and sum of the divergences of the velocity field in every cell of the domain.
The divergence is calculated as QT*q, which is technically the sum of the mass fluxes in every cell. This QT is also differs depending on the Direct Forcing method used.
Definition at line 118 of file DirectForcingSolver.cu.
Referenced by DirectForcingSolver< memoryType >::generateQT().
|
protected |
coefficients of interpolation (host)
Definition at line 37 of file DirectForcingSolver.h.
|
protected |
other coefficients of interpolation; quadratic interpolation (host)
Definition at line 37 of file DirectForcingSolver.h.
|
protected |
other coefficients of interpolation; quadratic interpolation (device)
Definition at line 39 of file DirectForcingSolver.h.
|
protected |
coefficients of interpolation (device)
Definition at line 39 of file DirectForcingSolver.h.
|
protected |
pressure field
Definition at line 44 of file DirectForcingSolver.h.
|
protected |
vector used to differentiate forcing points from regular ones (host)
Definition at line 32 of file DirectForcingSolver.h.
|
protected |
vector used to differentiate forcing points from regular ones (host)
Definition at line 32 of file DirectForcingSolver.h.
|
protected |
vector used to differentiate forcing points from regular ones (device)
Definition at line 35 of file DirectForcingSolver.h.
|
protected |
vector used to differentiate forcing points from regular ones (device)
Definition at line 35 of file DirectForcingSolver.h.
|
protected |
velocity field (host)
Definition at line 42 of file DirectForcingSolver.h.
|
protected |
velocity field (device)
Definition at line 43 of file DirectForcingSolver.h.