Fluid API reference

API Reference

namespace Fluid

Enums

enum VelocityOp

Operator when applying velocity to velocity field: add or set.

Values:

Add
Set

Functions

LinearSolver::Parameters FixedParams(unsigned iterations)

Create a linear solver parameters object with fixed solver type.

Return
parameters
Parameters
  • iterations: number of iterations to do

LinearSolver::Parameters IterativeParams(float errorTolerance)

Create a linear solver parameters object, solver will continue until error tolerance is reached.

Return
parameters
Parameters
  • errorTolerance: tolerance to reach before exiting

float DefaultParticleSize()

Variables

Renderer::ColorBlendState IntersectionBlend
Renderer::ColorBlendState UnionBlend
Renderer::Clear BoundariesClear
class Advection
#include <Advection.h>

Advects particles, velocity field or any field using a velocity field.

Public Functions

Advection(const Renderer::Device &device, const glm::ivec2 &size, float dt, Velocity &velocity, Velocity::InterpolationMode interpolationMode)

Initialize advection kernels and related object.

Parameters
  • device: vulkan device
  • size: size of velocity field
  • dt: delta time for integration
  • velocity: velocity field

void AdvectVelocity()

Self advect velocity.

void AdvectBind(Density &density)

Binds a density field to be advected.

Parameters
  • density: density field

void Advect()

Performs an advection of the density field. Asynchronous operation.

void AdvectParticleBind(Renderer::GenericBuffer &particles, Renderer::Texture &levelSet, Renderer::IndirectBuffer<Renderer::DispatchParams> &dispatchParams)

Binds praticles to be advected. Also use a level set to project out the particles if they enter it.

Parameters
  • particles: particles to be advected
  • levelSet: level set to project out particles
  • dispatchParams: contains number of particles

void AdvectParticles()

Advect particles. Asynchrounous operation.

class Cfl
#include <Cfl.h>

Calculates the CFL number of the velocity field. It’s an indication on how to choose your time step size. Ideally, the time step should be smaller than the CFL number.

Public Functions

void Compute()

Compute the CFL number. Non-blocking.

float Get()

Returns the CFL number. Blocking.

Return
cfl number

class Circle : public Vortex2D::Renderer::Transformable, public Vortex2D::Renderer::Drawable
#include <Boundaries.h>

Signed distance field of circle.

Public Functions

Circle(const Renderer::Device &device, float radius, float extent = 10.0f)

Initialize the circle with radius and extend of signed distance.

Parameters
  • device: vulkan device.
  • radius: radius of circle.
  • extent: extend how far from the circle the signed distance field is calculated.

void Initialize(const Renderer::RenderState &renderState)

Initialize the drawable for a particular state. This might include creating the correct pipeline. If it was already initialized, it will do nothing.

Parameters
  • renderState: the state to initialize with.

void Update(const glm::mat4 &projection, const glm::mat4 &view)

Update the MVP matrix of the drawable.

Parameters
  • projection: the projection matrix
  • view: the view matrix

void Draw(vk::CommandBuffer commandBuffer, const Renderer::RenderState &renderState)

Draw for the given render state. This has to be initialized before.

Parameters
  • commandBuffer: the command buffer to record into.
  • renderState: the render state to use.

class ConjugateGradient : public Vortex2D::Fluid::LinearSolver
#include <ConjugateGradient.h>

An iterative preconditioned conjugate linear solver. The preconditioner can be specified.

Public Functions

ConjugateGradient(const Renderer::Device &device, const glm::ivec2 &size, Preconditioner &preconditioner)

Initialize the solver with a size and preconditioner.

Parameters
  • device: vulkan device
  • size:
  • preconditioner:

void Bind(Renderer::GenericBuffer &d, Renderer::GenericBuffer &l, Renderer::GenericBuffer &b, Renderer::GenericBuffer &x)

Bind the buffers for the linear solver.

Parameters
  • d: the diagonal of the matrxi
  • l: the lower matrix
  • b: the right hand side
  • x: the unknowns

void BindRigidbody(float delta, Renderer::GenericBuffer &d, RigidBody &rigidBody)

Bind rigidbody with the linear solver’s matrix.

Parameters
  • delta: solver delta
  • d: diagonal matrix
  • rigidBody: rigidbody to bind to

void Solve(Parameters &params, const std::vector<RigidBody *> &rigidbodies = {})

Solve iteratively solve the linear equations in data.

float GetError()

Return
the max error

class Density : public Vortex2D::Renderer::RenderTexture, public Vortex2D::Renderer::Sprite
#include <Density.h>

Density field, used to represent smoke swirling.

class Depth
#include <Multigrid.h>

Contains the sizes of the multigrid hierarchy.

Public Functions

Depth(const glm::ivec2 &size)

Initialize with the finest size.

Parameters
  • size: the base size.

int GetMaxDepth() const

The calculated depth of the multigrid.

Return
the depth.

glm::ivec2 GetDepthSize(std::size_t i) const

Gets the depth for a given level.

Return
the size
Parameters
  • i: the level

class Diagonal : public Vortex2D::Fluid::Preconditioner
#include <Diagonal.h>

Diagonal preconditioner. Simplest of preconditioner, useful to verify if the preconditioned conjugate gradient works.

Public Functions

void Bind(Renderer::GenericBuffer &d, Renderer::GenericBuffer &l, Renderer::GenericBuffer &b, Renderer::GenericBuffer &x)

Bind the linear equation buffers.

Parameters
  • d: the diagonal of the matrix
  • l: the lower matrix
  • b: the right hand side
  • x: the unknown buffer

void Record(vk::CommandBuffer commandBuffer)

Record the preconditioner.

Parameters
  • commandBuffer: the command buffer to record into.

class DistanceField : public Vortex2D::Renderer::AbstractSprite
#include <Boundaries.h>

Sprite of a distance field.

Public Functions

DistanceField(const Renderer::Device &device, Renderer::RenderTexture &levelSet, float scale = 1.0f)

Initialize the price with the level set and scale.

Parameters
  • device: vulkan device
  • levelSet: level set to use as sprite
  • scale: scale of the level set

void Draw(vk::CommandBuffer commandBuffer, const Renderer::RenderState &renderState)

Draw for the given render state. This has to be initialized before.

Parameters
  • commandBuffer: the command buffer to record into.
  • renderState: the render state to use.

class Extrapolation
#include <Extrapolation.h>

Class to extrapolate values into the neumann and/or dirichlet boundaries.

Public Functions

void Extrapolate()

Will extrapolate values from buffer into the dirichlet and neumann boundaries.

void ConstrainBind(Renderer::Texture &solidPhi)

Binds a solid level set to use later and constrain the velocity against.

Parameters
  • solidPhi: solid level set

void ConstrainVelocity()

Constrain the velocity, i.e. ensure that the velocity normal to the solid level set is 0.

class GaussSeidel : public Vortex2D::Fluid::LinearSolver, public Vortex2D::Fluid::Preconditioner
#include <GaussSeidel.h>

An iterative black and red successive over relaxation linear solver.

Public Functions

void Bind(Renderer::GenericBuffer &d, Renderer::GenericBuffer &l, Renderer::GenericBuffer &b, Renderer::GenericBuffer &x)

Bind the buffers for the linear solver.

Parameters
  • d: the diagonal of the matrxi
  • l: the lower matrix
  • b: the right hand side
  • x: the unknowns

void BindRigidbody(float delta, Renderer::GenericBuffer &d, RigidBody &rigidBody)

Bind rigidbody with the linear solver’s matrix.

Parameters
  • delta: solver delta
  • d: diagonal matrix
  • rigidBody: rigidbody to bind to

void Solve(Parameters &params, const std::vector<RigidBody *> &rigidbodies = {})

Iterative solving of the linear equations in data.

float GetError()

Return
the max error

void Record(vk::CommandBuffer commandBuffer)

Record the preconditioner.

Parameters
  • commandBuffer: the command buffer to record into.

void Record(vk::CommandBuffer commandBuffer, int iterations)

Record a determined number of iterations.

Parameters
  • commandBuffer:
  • iterations:

void SetW(float w)

Set the w factor of the GS iterations : x_new = w * x_new + (1-w) * x_old.

Parameters
  • w:

void SetPreconditionerIterations(int iterations)

set number of iterations to be used when GS is a preconditioner

Parameters
  • iterations:

class IncompletePoisson : public Vortex2D::Fluid::Preconditioner
#include <IncompletePoisson.h>

Incomplete poisson preconditioner. Slightly better than a simple diagonal preconditioner.

Public Functions

void Bind(Renderer::GenericBuffer &d, Renderer::GenericBuffer &l, Renderer::GenericBuffer &b, Renderer::GenericBuffer &x)

Bind the linear equation buffers.

Parameters
  • d: the diagonal of the matrix
  • l: the lower matrix
  • b: the right hand side
  • x: the unknown buffer

void Record(vk::CommandBuffer commandBuffer)

Record the preconditioner.

Parameters
  • commandBuffer: the command buffer to record into.

class Jacobi : public Vortex2D::Fluid::Preconditioner
#include <Jacobi.h>

An iterative jacobi linear solver.

Public Functions

void Bind(Renderer::GenericBuffer &d, Renderer::GenericBuffer &l, Renderer::GenericBuffer &b, Renderer::GenericBuffer &x)

Bind the linear equation buffers.

Parameters
  • d: the diagonal of the matrix
  • l: the lower matrix
  • b: the right hand side
  • x: the unknown buffer

void Record(vk::CommandBuffer commandBuffer)

Record the preconditioner.

Parameters
  • commandBuffer: the command buffer to record into.

void SetW(float w)

Set the w factor of the GS iterations : x_new = w * x_new + (1-w) * x_old.

Parameters
  • w:

void SetPreconditionerIterations(int iterations)

set number of iterations to be used when GS is a preconditioner

Parameters
  • iterations:

class LevelSet : public Vortex2D::Renderer::RenderTexture
#include <LevelSet.h>

A signed distance field, which can be re-initialized. In other words, a level set.

Public Functions

void Reinitialise()

Reinitialise the level set, i.e. ensure it is a correct signed distance field.

void ShrinkWrap()

Shrink wrap wholes.

void ExtrapolateBind(Renderer::Texture &solidPhi)

Bind a solid level set, which will be used to extrapolate into this level set.

Parameters
  • solidPhi:

void Extrapolate()

Extrapolate this level set into the solid level set it was attached to. This only performs a single cell extrapolation.

struct LinearSolver
#include <LinearSolver.h>

An interface to represent a linear solver.

Subclassed by Vortex2D::Fluid::ConjugateGradient, Vortex2D::Fluid::GaussSeidel, Vortex2D::Fluid::Multigrid

Public Functions

virtual void Bind(Renderer::GenericBuffer &d, Renderer::GenericBuffer &l, Renderer::GenericBuffer &b, Renderer::GenericBuffer &x) = 0

Bind the buffers for the linear solver.

Parameters
  • d: the diagonal of the matrxi
  • l: the lower matrix
  • b: the right hand side
  • x: the unknowns

virtual void BindRigidbody(float delta, Renderer::GenericBuffer &d, RigidBody &rigidBody) = 0

Bind rigidbody with the linear solver’s matrix.

Parameters
  • delta: solver delta
  • d: diagonal matrix
  • rigidBody: rigidbody to bind to

virtual void Solve(Parameters &params, const std::vector<RigidBody *> &rigidBodies = {}) = 0

Solves the linear equations.

Parameters
  • params: solver iteration/error parameters
  • rigidBodies: rigidbody to include in solver’s matrix

virtual float GetError() = 0

Return
the max error

struct Data
#include <LinearSolver.h>

The various parts of linear equations.

struct DebugCopy
#include <LinearSolver.h>

Copies the linear solver data in the debug linear solver data.

Public Functions

void Copy()

Copies the linear solver data in the debug linear solver data.

struct DebugData
#include <LinearSolver.h>

Contains the linear equations as texture, so it can easily be visualised in RenderDoc.

class Error
#include <LinearSolver.h>

Calculates the max residual error of the linear system.

Public Functions

void Bind(Renderer::GenericBuffer &d, Renderer::GenericBuffer &l, Renderer::GenericBuffer &div, Renderer::GenericBuffer &pressure)

Bind the linear system.

Parameters
  • d: the diagonal of the matrxi
  • l: the lower matrix
  • b: the right hand side
  • x: the unknowns

Error &Submit()

Submit the error calculation.

Return
this.

Error &Wait()

Wait for error to be calculated.

Return
this.

float GetError()

Get the maximum error.

Return
The error.

struct Parameters
#include <LinearSolver.h>

Parameters for an iterative linear solvers.

Public Types

enum SolverType

Run the solver a fixed number of step or until we reached a minimum error.

Values:

Fixed
Iterative

Public Functions

Parameters(SolverType type, unsigned iterations, float errorTolerance = 0.0f)

Construct parameters with max iterations and max error.

Parameters
  • type: fixed or iterative type of solver
  • iterations: max number of iterations to perform
  • errorTolerance: solver stops when the error is smaller than this.

bool IsFinished(float initialError) const

Checks if we’ve reacched the parameters.

Return
if we can stop the linear solver.
Parameters
  • initialError: the initial error

void Reset()

Sets the out error and out iterations to 0.

class LocalGaussSeidel : public Vortex2D::Fluid::Preconditioner
#include <GaussSeidel.h>

A version of the gauss seidel that can only be applied on sizes (16,16) or smaller.

Public Functions

void Bind(Renderer::GenericBuffer &d, Renderer::GenericBuffer &l, Renderer::GenericBuffer &b, Renderer::GenericBuffer &x)

Bind the linear equation buffers.

Parameters
  • d: the diagonal of the matrix
  • l: the lower matrix
  • b: the right hand side
  • x: the unknown buffer

void Record(vk::CommandBuffer commandBuffer)

Record the preconditioner.

Parameters
  • commandBuffer: the command buffer to record into.

class Multigrid : public Vortex2D::Fluid::LinearSolver, public Vortex2D::Fluid::Preconditioner
#include <Multigrid.h>

Multigrid preconditioner. It creates a hierarchy of twice as small set of linear equations. It applies a few iterations of jacobi on each level and transfers the error on the level above. It then copies the error down, adds to the current solution and apply a few more iterations of jacobi.

Public Functions

Multigrid(const Renderer::Device &device, const glm::ivec2 &size, float delta, int numSmoothingIterations = 3, SmootherSolver smoother = SmootherSolver::Jacobi)

Initialize multigrid for given size and delta.

Parameters
  • device: vulkan device
  • size: of the linear equations
  • delta: timestep delta

void Bind(Renderer::GenericBuffer &d, Renderer::GenericBuffer &l, Renderer::GenericBuffer &b, Renderer::GenericBuffer &x)

Bind the buffers for the linear solver.

Parameters
  • d: the diagonal of the matrxi
  • l: the lower matrix
  • b: the right hand side
  • x: the unknowns

void BuildHierarchiesBind(Pressure &pressure, Renderer::Texture &solidPhi, Renderer::Texture &liquidPhi)

Bind the level sets from which the hierarchy is built.

Parameters
  • pressure: The current linear equations
  • solidPhi: the solid level set
  • liquidPhi: the liquid level set

void BuildHierarchies()

Computes the hierarchy to be used by the multigrid. Asynchronous operation.

void Record(vk::CommandBuffer commandBuffer)

Record the preconditioner.

Parameters
  • commandBuffer: the command buffer to record into.

void BindRigidbody(float delta, Renderer::GenericBuffer &d, RigidBody &rigidBody)

Bind rigidbody with the linear solver’s matrix.

Parameters
  • delta: solver delta
  • d: diagonal matrix
  • rigidBody: rigidbody to bind to

void Solve(Parameters &params, const std::vector<RigidBody *> &rigidBodies = {})

Solves the linear equations.

Parameters
  • params: solver iteration/error parameters
  • rigidBodies: rigidbody to include in solver’s matrix

float GetError()

Return
the max error

class ParticleCount : public Vortex2D::Renderer::RenderTexture
#include <Particles.h>

Container for particles used in the advection of the fluid simulation. Also a level set that is built from the particles.

Public Functions

void Scan()

Count the number of particles and update the internal data structures.

int GetTotalCount()

Calculate the total number of particles and return it.

Return

Renderer::IndirectBuffer<Renderer::DispatchParams> &GetDispatchParams()

Calculate the dispatch parameters to use on the particle buffer.

Return

void LevelSetBind(LevelSet &levelSet)

Bind a solid level set, which will be used to interpolate the particles out of.

Parameters
  • levelSet:

void Phi()

Calculate the level set from the particles.

void VelocitiesBind(Velocity &velocity, Renderer::GenericBuffer &valid)

Bind the velocities, used for advection of the particles.

Parameters
  • velocity:
  • valid:

void TransferToGrid()

Interpolate the velocities of the particles to the velocities field.

void TransferFromGrid()

Interpolate the velocities field in to the particles’ velocity.

class Polygon : public Vortex2D::Renderer::Transformable, public Vortex2D::Renderer::Drawable
#include <Boundaries.h>

Signed distance field of a poylgon.

Subclassed by Vortex2D::Fluid::Rectangle

Public Functions

Polygon(const Renderer::Device &device, std::vector<glm::vec2> points, bool inverse = false, float extent = 10.0f)

Initialize polygon with set of points and extent of signed distance.

Parameters
  • device: vulkan device
  • points: clockwise oriented set of points (mininum 3).
  • inverse: flag if the distance field should be inversed.
  • extent: extend how far from the poylon the signed distance field is calculated.

void Initialize(const Renderer::RenderState &renderState)

Initialize the drawable for a particular state. This might include creating the correct pipeline. If it was already initialized, it will do nothing.

Parameters
  • renderState: the state to initialize with.

void Update(const glm::mat4 &projection, const glm::mat4 &view)

Update the MVP matrix of the drawable.

Parameters
  • projection: the projection matrix
  • view: the view matrix

void Draw(vk::CommandBuffer commandBuffer, const Renderer::RenderState &renderState)

Draw for the given render state. This has to be initialized before.

Parameters
  • commandBuffer: the command buffer to record into.
  • renderState: the render state to use.

struct Preconditioner
#include <Preconditioner.h>

An interface to represent a linear solver preconditioner.

Subclassed by Vortex2D::Fluid::Diagonal, Vortex2D::Fluid::GaussSeidel, Vortex2D::Fluid::IncompletePoisson, Vortex2D::Fluid::Jacobi, Vortex2D::Fluid::LocalGaussSeidel, Vortex2D::Fluid::Multigrid

Public Functions

virtual void Bind(Renderer::GenericBuffer &d, Renderer::GenericBuffer &l, Renderer::GenericBuffer &b, Renderer::GenericBuffer &x) = 0

Bind the linear equation buffers.

Parameters
  • d: the diagonal of the matrix
  • l: the lower matrix
  • b: the right hand side
  • x: the unknown buffer

virtual void Record(vk::CommandBuffer commandBuffer) = 0

Record the preconditioner.

Parameters
  • commandBuffer: the command buffer to record into.

class PrefixScan
#include <PrefixScan.h>

The prefix sum operator.

void PrefixSym(int input[], int n, int output[])
{
    output[0] = input[0];

    for (int i = 1; i < n; i++)
        output[i] = output[i-1] + input[i];
}

class Bound
#include <PrefixScan.h>

A prefix scan object bound with input/output buffers, ready to be dispatched.

class Pressure
#include <Pressure.h>

build the linear equation and compute the divergence from the resulting solution.

Public Functions

Renderer::Work::Bound BindMatrixBuild(const glm::ivec2 &size, Renderer::GenericBuffer &diagonal, Renderer::GenericBuffer &lower, Renderer::Texture &liquidPhi, Renderer::Texture &solidPhi)

Bind the various buffes for the linear system Ax = b.

Return
Parameters
  • size: size of the linear system
  • diagonal: diagonal of A
  • lower: lower matrix of A
  • liquidPhi: liquid level set
  • solidPhi: solid level set

void BuildLinearEquation()

Build the matrix A and right hand side b.

void ApplyPressure()

Apply the solution of the equation Ax = b, i.e. the pressure to the velocity to make it non-divergent.

class Rectangle : public Vortex2D::Fluid::Polygon
#include <Boundaries.h>

Signed distance field of a rectangle.

Public Functions

Rectangle(const Renderer::Device &device, const glm::vec2 &size, bool inverse = false, float extent = 10.0f)

Initialize rectangle with size and extend of signed distance.

Parameters
  • device: vulkan device.
  • size: rectangle size
  • inverse: flag if the distance field should be inverted.
  • extent: extent how far from the rectangle the signed distance field is calculated.

void Initialize(const Renderer::RenderState &renderState)

Initialize the drawable for a particular state. This might include creating the correct pipeline. If it was already initialized, it will do nothing.

Parameters
  • renderState: the state to initialize with.

void Update(const glm::mat4 &projection, const glm::mat4 &view)

Update the MVP matrix of the drawable.

Parameters
  • projection: the projection matrix
  • view: the view matrix

void Draw(vk::CommandBuffer commandBuffer, const Renderer::RenderState &renderState)

Draw for the given render state. This has to be initialized before.

Parameters
  • commandBuffer: the command buffer to record into.
  • renderState: the render state to use.

class Reduce
#include <Reduce.h>

Parallel reduction of a buffer into one value. The operator and type of data is specified by inheriting the class.

Subclassed by Vortex2D::Fluid::ReduceJ, Vortex2D::Fluid::ReduceMax, Vortex2D::Fluid::ReduceSum

Public Functions

Reduce::Bound Bind(Renderer::GenericBuffer &input, Renderer::GenericBuffer &output)

Bind the reduce operation.

Return
a bound object that can be recorded in a command buffer.
Parameters
  • input: input buffer
  • output: output buffer

class Bound
#include <Reduce.h>

Bound input and output buffer for a reduce operation.

Public Functions

void Record(vk::CommandBuffer commandBuffer)

Record the reduce operation.

Parameters
  • commandBuffer: the command buffer to record into.

class ReduceJ : public Vortex2D::Fluid::Reduce
#include <Reduce.h>

Reduce operation on a struct with a 2d vector and 1 float (i.e. 3 floats) with addition.

Public Functions

ReduceJ(const Renderer::Device &device, const glm::ivec2 &size)

Initialize reduce with device and 2d size.

Parameters
  • device:
  • size:

class ReduceMax : public Vortex2D::Fluid::Reduce
#include <Reduce.h>

Reduce operation on float with max of absolute.

Public Functions

ReduceMax(const Renderer::Device &device, const glm::ivec2 &size)

Initialize reduce with device and 2d size.

Parameters
  • device:
  • size:

class ReduceSum : public Vortex2D::Fluid::Reduce
#include <Reduce.h>

Reduce operation on float with addition.

Public Functions

ReduceSum(const Renderer::Device &device, const glm::ivec2 &size)

Initialize reduce with device and 2d size.

Parameters
  • device:
  • size:

class RigidBody : public Vortex2D::Renderer::Transformable
#include <Rigidbody.h>

Rigidbody that can interact with the fluid: either be push by it, or influence it, or both.

Public Functions

virtual void ApplyForces()

function to override and apply forces from this rigidbody to the external rigidbody

virtual void ApplyVelocities()

Override and apply velocities from the external rigidbody to the this rigidbody.

void SetMassData(float mass, float inertia)

Sets the mass and inertia of the rigidbody.

Parameters
  • mass: of the body
  • inertia: of the body

void SetVelocities(const glm::vec2 &velocity, float angularVelocity)

sets the velocities and angular velocities of the body

Parameters
  • velocity:
  • angularVelocity:

void UpdatePosition()

Upload the transform matrix to the GPU.

void RenderPhi()

Render the current object orientation in an internal texture and the external one.

void BindPhi(Renderer::RenderTexture &phi)

Bind the rendertexture where this rigidbodies shape will be rendered.

Parameters
  • phi: render texture of the world

void BindDiv(Renderer::GenericBuffer &div, Renderer::GenericBuffer &diagonal)

Bind a the right hand side and diagonal of the linear system Ax = b. This is to apply the rigid body influence to the system.

Parameters
  • div: right hand side of the linear system Ax=b
  • diagonal: diagonal of matrix A

void BindVelocityConstrain(Fluid::Velocity &velocity)

Bind velocities to constrain based on the body’s velocity.

Parameters
  • velocity:

void BindForce(Renderer::GenericBuffer &d, Renderer::GenericBuffer &pressure)

Bind pressure, to have the pressure update the body’s forces.

Parameters
  • d: diagonal of matrix A
  • pressure: solved pressure buffer

void BindPressure(float delta, Renderer::GenericBuffer &d, Renderer::GenericBuffer &s, Renderer::GenericBuffer &z)

Bind pressure, to have the pressure update the body’s forces.

Parameters
  • delta:
  • d:
  • s:
  • z:

void Div()

Apply the body’s velocities to the linear equations matrix A and right hand side b.

void Force()

Apply the pressure to body, updating its forces.

void Pressure()

Reduce the force for pressure update.

void VelocityConstrain()

Constrain the velocities field based on the body’s velocity.

Velocity GetForces()

Download the forces from the GPU and return them.

Return

vk::Flags<Type> GetType()

Type of this body.

Return

void SetType(Type type)

Set the type of the body.

Parameters
  • type:

Renderer::RenderTexture &Phi()

the local level set of the body

Return

class RigidBodySolver
#include <Rigidbody.h>

Interface to call the external rigidbody solver.

Public Functions

virtual void Step(float delta) = 0

perfoms a single step of the solver.

Parameters
  • delta: of simulation

class SmokeWorld : public Vortex2D::Fluid::World
#include <World.h>

A concrete implementation of World to simulate ‘smoke’, or more accurately dye in a liquid. The liquid cannot change location or size.

Public Functions

void FieldBind(Density &density)

Bind a density field to be moved around with the fluid.

Parameters
  • density: the density field

class Transfer
#include <Transfer.h>

Prolongates or restrict a level set on a finer or coarser level set.

Public Functions

Transfer(const Renderer::Device &device)

Initialize prolongate and restrict compute pipelines.

Parameters
  • device:

void ProlongateBind(std::size_t level, const glm::ivec2 &fineSize, Renderer::GenericBuffer &fine, Renderer::GenericBuffer &fineDiagonal, Renderer::GenericBuffer &coarse, Renderer::GenericBuffer &coarseDiagonal)

Prolongate a level set on a finer level set. Setting the 4 cells to the value of the coarser grid. Multiple level sets can be bound and indexed.

Parameters
  • level: the index of the bound level set to prolongate
  • fineSize: size of the finer level set
  • fine: the finer level set
  • fineDiagonal: the diagonal of the linear equation matrix at size fineSize
  • coarse: the coarse level set
  • coarseDiagonal: the diagonal of the linear equation matrix at size half of fineSize

void RestrictBind(std::size_t level, const glm::ivec2 &fineSize, Renderer::GenericBuffer &fine, Renderer::GenericBuffer &fineDiagonal, Renderer::GenericBuffer &coarse, Renderer::GenericBuffer &coarseDiagonal)

Restricing the level set on a coarser level set. Averages 4 cells into one. Multiple level sets can be bound and indexed.

Parameters
  • level: the index of the bound level set to prolongate
  • fineSize: size of the finer level set
  • fine: the finer level set
  • fineDiagonal: the diagonal of the linear equation matrix at size fineSize
  • coarse: the coarse level set
  • coarseDiagonal: the diagonal of the linear equation matrix at size half of fineSize

void Prolongate(vk::CommandBuffer commandBuffer, std::size_t level)

Prolongate the level set, using the bound level sets at the specified index.

Parameters
  • commandBuffer: command buffer to record into.
  • level: index of bound level sets.

void Restrict(vk::CommandBuffer commandBuffer, std::size_t level)

Restrict the level set, using the bound level sets at the specified index.

Parameters
  • commandBuffer: command buffer to record into.
  • level: index of bound level sets.

class Velocity : public Vortex2D::Renderer::RenderTexture
#include <Velocity.h>

The Velocity field. Can be used to calculate a difference between different states. Contains three fields: intput and output, used for ping-pong algorithms, and d, the difference between two velocity fields.

Public Types

enum InterpolationMode

Velocity interpolation when querying in the shader with non-integer locations.

Values:

Linear = 0
Cubic = 1

Public Functions

Renderer::Texture &Output()

An output texture used for algorithms that used the velocity as input and need to create a new velocity field.

Return

Renderer::Texture &D()

A difference velocity field, calculated with the difference between this velocity field, and the output velocity field.

Return

void CopyBack(vk::CommandBuffer commandBuffer)

Copy the output field to the main field.

Parameters
  • commandBuffer:

void Clear(vk::CommandBuffer commandBuffer)

Clear the velocity field.

Parameters
  • commandBuffer:

void SaveCopy()

Copy to the difference field.

void VelocityDiff()

Calculate the difference between the difference field and this velocity field, store it in the diference field.

class WaterWorld : public Vortex2D::Fluid::World
#include <World.h>

A concrete implementation of World to simulate water.

Public Functions

Renderer::RenderCommand RecordParticleCount(Renderer::RenderTarget::DrawableList drawables)

The water simulation uses particles to define the water area. In fact, the level set is built from the particles. This means to be able to set an area, we can’t use RecordLiquidPhi. To define the particle area, simply draw a regular shape. The colour r is used to determine if we add or remove particles, use r = 4 to add and r = -4 to remove.

Return
render command
Parameters
  • drawables: list of drawables object with colour 8 or -8

void ParticlePhi()

Using the particles, create a level set (phi) encompassing all the particles. This can be viewed with LiquidDistanceField.

class World
#include <World.h>

The main class of the framework. Each instance manages a grid and this class is used to set forces, define boundaries, solve the incompressbility equations and do the advection.

Subclassed by Vortex2D::Fluid::SmokeWorld, Vortex2D::Fluid::WaterWorld

Public Functions

World(const Renderer::Device &device, const glm::ivec2 &size, float dt, int numSubSteps = 1, Velocity::InterpolationMode interpolationMode = Velocity::InterpolationMode::Linear)

Construct an Engine with a size and time step.

Parameters
  • device: vulkan device
  • size: dimensions of the simulation
  • dt: timestamp of the simulation, e.g. 0.016 for 60FPS simulations.
  • numSubSteps: the number of sub-steps to perform per step call. Reduces loss of fluid.

void Step(LinearSolver::Parameters &params)

Perform one step of the simulation.

Renderer::RenderCommand RecordVelocity(Renderer::RenderTarget::DrawableList drawables, VelocityOp op)

Record drawables to the velocity field. The colour (r,g) will be used as the velocity (x, y)

Return
render command
Parameters
  • drawables: a list of drawable field
  • op: operation of the render: add velocity or set velocity

void SubmitVelocity(Renderer::RenderCommand &renderCommand)

submit the render command created with RecordVelocity

Parameters
  • renderCommand: the render command

Renderer::RenderCommand RecordLiquidPhi(Renderer::RenderTarget::DrawableList drawables)

Record drawables to the liquid level set, i.e. to define the fluid area. The drawables need to make a signed distance field, if not the result is undefined.

Return
render command
Parameters
  • drawables: a list of signed distance field drawables

Renderer::RenderCommand RecordStaticSolidPhi(Renderer::RenderTarget::DrawableList drawables)

Record drawables to the solid level set, i.e. to define the boundary area. The drawables need to make a signed distance field, if not the result is undefined.

Return
render command
Parameters
  • drawables: a list of signed distance field drawables

DistanceField LiquidDistanceField()

Create sprite that can be rendered to visualize the liquid level set.

Return
a sprite

DistanceField SolidDistanceField()

Create sprite that can be rendered to visualize the solid level set.

Return
a sprite

void AddRigidbody(RigidBody &rigidbody)

Add a rigibody to the solver.

Parameters
  • rigidbody:

void RemoveRigidBody(RigidBody &rigidbody)

Remove a rigidbody from the solver.

Parameters
  • rigidbody:

void AttachRigidBodySolver(RigidBodySolver &rigidbodySolver)

Attach a rigidbody solver, e.g. box2d.

Parameters
  • rigidbodySolver:

float GetCFL()

Calculate the CFL number, i.e. the width divided by the max velocity.

Return
CFL number

Renderer::Texture &GetVelocity()

Get the velocity, can be used to display it.

Return
velocity field reference