WaterLily
Introduction and Quickstart
WaterLily
— ModuleWaterLily.jl
Overview
WaterLily.jl is a simple and fast fluid simulator written in pure Julia. This is an experimental project to take advantage of the active scientific community in Julia to accelerate and enhance fluid simulations. Watch the JuliaCon2021 talk here:
Method/capabilities
WaterLily.jl solves the unsteady incompressible 2D or 3D Navier-Stokes equations on a Cartesian grid. The pressure Poisson equation is solved with a geometric multigrid method. Solid boundaries are modelled using the Boundary Data Immersion Method. The solver can run on serial CPU, multi-threaded CPU, or GPU backends.
Examples
The user can set the boundary conditions, the initial velocity field, the fluid viscosity (which determines the Reynolds number), and immerse solid obstacles using a signed distance function. These examples and others are found in the examples.
Flow over a circle
We define the size of the simulation domain as n
xm
cells. The circle has radius m/8
and is centered at (m/2,m/2)
. The flow boundary conditions are (U=1,0)
and Reynolds number is Re=U*radius/ν
where ν
(Greek "nu" U+03BD, not Latin lowercase "v") is the kinematic viscosity of the fluid.
using WaterLily
function circle(n,m;Re=250,U=1)
radius, center = m/8, m/2
body = AutoBody((x,t)->√sum(abs2, x .- center) - radius)
Simulation((n,m), (U,0), radius; ν=U*radius/Re, body)
end
The second to last line defines the circle geometry using a signed distance function. The AutoBody
function uses automatic differentiation to infer the other geometric parameter automatically. Replace the circle's distance function with any other, and now you have the flow around something else... such as a donut or the Julia logo. Finally, the last line defines the Simulation
by passing in parameters we've defined.
Now we can create a simulation (first line) and run it forward in time (third line)
circ = circle(3*2^6,2^7)
t_end = 10
sim_step!(circ,t_end)
Note we've set n,m
to be multiples of powers of 2, which is important when using the (very fast) Multi-Grid solver. We can now access and plot whatever variables we like. For example, we could print the velocity at I::CartesianIndex
using println(circ.flow.u[I])
or plot the whole pressure field using
using Plots
contour(circ.flow.p')
A set of flow metric functions have been implemented and the examples use these to make gifs such as the one above.
3D Taylor Green Vortex
The three-dimensional Taylor Green Vortex demonstrates many of the other available simulation options. First, you can simulate a nontrivial initial velocity field by passing in a vector function uλ(i,xyz)
where i ∈ (1,2,3)
indicates the velocity component uᵢ
and xyz=[x,y,z]
is the position vector.
function TGV(; pow=6, Re=1e5, T=Float64, mem=Array)
# Define vortex size, velocity, viscosity
L = 2^pow; U = 1; ν = U*L/Re
# Taylor-Green-Vortex initial velocity field
function uλ(i,xyz)
x,y,z = @. (xyz-1.5)*π/L # scaled coordinates
i==1 && return -U*sin(x)*cos(y)*cos(z) # u_x
i==2 && return U*cos(x)*sin(y)*cos(z) # u_y
return 0. # u_z
end
# Initialize simulation
return Simulation((L, L, L), (0, 0, 0), L; U, uλ, ν, T, mem)
end
This example also demonstrates the floating point type (T=Float64
) and array memory type (mem=Array
) options. For example, to run on an NVIDIA GPU we only need to import the CUDA.jl library and initialize the Simulation
memory on that device.
import CUDA
@assert CUDA.functional()
vortex = TGV(T=Float32,mem=CUDA.CuArray)
sim_step!(vortex,1)
For an AMD GPU, use import AMDGPU
and mem=AMDGPU.ROCArray
. Note that Julia 1.9 is required for AMD GPUs.
Moving bodies
You can simulate moving bodies in WaterLily by passing a coordinate map
to AutoBody
in addition to the sdf
.
using StaticArrays
function hover(L=2^5;Re=250,U=1,amp=π/4,ϵ=0.5,thk=2ϵ+√2)
# Line segment SDF
function sdf(x,t)
y = x .- SA[0,clamp(x[2],-L/2,L/2)]
√sum(abs2,y)-thk/2
end
# Oscillating motion and rotation
function map(x,t)
α = amp*cos(t*U/L); R = SA[cos(α) sin(α); -sin(α) cos(α)]
R * (x - SA[3L-L*sin(t*U/L),4L])
end
Simulation((6L,6L),(0,0),L;U,ν=U*L/Re,body=AutoBody(sdf,map),ϵ)
end
In this example, the sdf
function defines a line segment from -L/2 ≤ x[2] ≤ L/2
with a thickness thk
. To make the line segment move, we define a coordinate transformation function map(x,t)
. In this example, the coordinate x
is shifted by (3L,4L)
at time t=0
, which moves the center of the segment to this point. However, the horizontal shift varies harmonically in time, sweeping the segment left and right during the simulation. The example also rotates the segment using the rotation matrix R = [cos(α) sin(α); -sin(α) cos(α)]
where the angle α
is also varied harmonically. The combined result is a thin flapping line, similar to a cross-section of a hovering insect wing.
One important thing to note here is the use of StaticArrays
to define the sdf
and map
. This speeds up the simulation since it eliminates allocations at every grid cell and time step.
Circle inside an oscillating flow
This example demonstrates a 2D oscillating periodic flow over a circle.
function circle(n,m;Re=250,U=1)
# define a circle at the domain center
radius = m/8
body = AutoBody((x,t)->√sum(abs2, x .- (n/2,m/2)) - radius)
# define time-varying body force `g` and periodic direction `perdir`
accelScale, timeScale = U^2/2radius, radius/U
g(i,t) = i==1 ? -2accelScale*sin(t/timeScale) : 0
Simulation((n,m), (U,0), radius; ν=U*radius/Re, body, g, perdir=(1,))
end
The g
argument accepts a function with direction (i
) and time (t
) arguments. This allows you to create a spatially uniform body force with variations over time. In this example, the function adds a sinusoidal force in the "x" direction i=1
, and nothing to the other directions.
The perdir
argument is a tuple that specifies the directions to which periodic boundary conditions should be applied. Any number of directions may be defined as periodic, but in this example only the i=1
direction is used allowing the flow to accelerate freely in this direction.
WaterLily gives the posibility to set up a Simulation
using time-varying boundary conditions for the velocity field. This can be used to simulate a flow in an accelerating reference frame. The following example demonstrates how to set up a Simulation
with a time-varying velocity field.
using WaterLily
# define time-varying velocity boundary conditions
Ut(i,t::T;a0=0.5) where T = i==1 ? convert(T, a0*t) : zero(T)
# pass that to the function that creates the simulation
sim = Simulation((256,256), Ut, 32)
The Ut
function is used to define the time-varying velocity field. In this example, the velocity in the "x" direction is set to a0*t
where a0
is the acceleration of the reference frame. The Simulation
function is then called with the Ut
function as the second argument. The simulation will then run with the time-varying velocity field.
Periodic and convective boundary conditions
In addition to the standard free-slip (or reflective) boundary conditions, WaterLily also supports periodic boundary conditions. The following example demonstrates how to set up a Simulation
with periodic boundary conditions in the "y" direction.
using WaterLily,StaticArrays
# sdf an map for a moving circle in y-direction
function sdf(x,t)
norm2(SA[x[1]-192,mod(x[2]-384,384)-192])-32
end
function map(x,t)
x.-SA[0.,t/2]
end
# make a body
body = AutoBody(sdf, map)
# y-periodic boundary conditions
Simulation((512,384), (1,0), 32; body, perdir=(2,))
Additionally, the flag exitBC=true
can be passed to the Simulation
function to enable convective boundary conditions. This will apply a 1D convective exit in the x
direction (there is not way to change this at the moment). The exitBC
flag is set to false
by default. In this case, the boundary condition is set to the corresponding value of the u_BC
vector you specified when constructing the Simulation
.
using WaterLily
# make a body
body = AutoBody(sdf, map)
# y-periodic boundary conditions
Simulation((512,384), u_BC=(1,0), L=32; body, exitBC=true)
The following example demonstrates how to write simulation data to a .pvd
file using the WriteVTK
package and the WaterLily vtkwriter
function. The simplest writer can be instantiated with
using WaterLily,WriteVTK
# make a sim
sim = make_sim(...)
# make a writer
writer = vtkwriter("simple_writer")
# write the data
write!(writer,sim)
# don't forget to close the file
close(writer)
This would write the velocity and pressure fields to a file named simmple_writer.pvd
. The vtkwriter
function can also take a dictionary of custom attributes to write to the file. For example, to write the body (sdf) and λ₂ fields to the file, you could use the following code:
using WaterLily,WriteVTK
# make a writer with some attributes, need to output to CPU array to save file (|> Array)
velocity(a::Simulation) = a.flow.u |> Array;
pressure(a::Simulation) = a.flow.p |> Array;
_body(a::Simulation) = (measure_sdf!(a.flow.σ, a.body, WaterLily.time(a));
a.flow.σ |> Array;)
lamda(a::Simulation) = (@inside a.flow.σ[I] = WaterLily.λ₂(I, a.flow.u);
a.flow.σ |> Array;)
# this maps field names to values in the file
custom_attrib = Dict(
"Velocity" => velocity,
"Pressure" => pressure,
"Body" => _body,
"Lambda" => lamda
)
# make the writer
writer = vtkWriter("advanced_writer"; attrib=custom_attrib)
...
close(writer)
The functions that are passed to the attrib
(custom attributes) must follow the same structure as what is shown in this example, that is, given a Simulation
, return a N-dimensional (scalar or vector) field. The vtkwriter
function will automatically write the data to a .pvd
file, which can be read by Paraview. The prototype for the vtkwriter
function is:
# prototype vtk writer function
custom_vtk_function(a::Simulation) = ... |> Array
the ...
should be replaced with the code that generates the field you want to write to the file. The piping to a (CPU) Array
is necessary to ensure that the data is written to the CPU before being written to the file for GPU simulations.
This capability is very usefull to restart a simulation from a previous state. The ReadVTK
package is used to read simulation data from a .pvd
file. This .pvd
must have been writen with the vtkwriter
function and must contain at least the velocity
and pressure
fields. The following example demonstrates how to restart a simulation from a .pvd
file using the ReadVTK
package and the WaterLily vtkreader
function
using WaterLily,ReadVTK
sim = make_sim(...)
# restart the simulation
writer = restart_sim!(sim; fname="file_restart.pvd")
# this acctually append the data to the file used to restart
write!(writer, sim)
# don't forget to close the file
close(writer)
Internally, this function reads the last file in the .pvd
file and use that to set the velocity
and pressure
fields in the simulation. The sim_time
is also set to the last value saved in the .pvd
file. The function also returns a vtkwriter
that will append the new data to the file used to restart the simulation. Note the sim
that will be filled must be identical to the one saved to the file for this restart to work, that is, the same size, same body, etc.
Multi-threading and GPU backends
WaterLily uses KernelAbstractions.jl to multi-thread on CPU and run on GPU backends. The implementation method and speed-up are documented in our ParCFD abstract. In summary, a single macro WaterLily.@loop
is used for nearly every loop in the code base, and this uses KernelAbstractactions to generate optimized code for each back-end. The speed-up is more pronounce for large simulations, and we've benchmarked up to 23x-speed up on a Intel Core i7-10750H x6 processor, and 182x speed-up NVIDIA GeForce GTX 1650 Ti GPU card.
Note that multi-threading requires starting Julia with the --threads
argument, see the multi-threading section of the manual. If you are running Julia with multiple threads, KernelAbstractions will detect this and multi-thread the loops automatically. As in the Taylor-Green-Vortex examples above, running on a GPU requires initializing the Simulation
memory on the GPU, and care needs to be taken to move the data back to the CPU for visualization. See jelly fish for another non-trivial example.
Finally, KernelAbstractions does incur some CPU allocations for every loop, but other than this sim_step!
is completely non-allocating. This is one reason why the speed-up improves as the size of the simulation increases.
Development goals
- Immerse obstacles defined by 3D meshes using GeometryBasics.
- Multi-CPU/GPU simulations.
- Add free-surface physics with Volume-of-Fluid or Level-Set.
- Add external potential-flow domain boundary conditions.
If you have other suggestions or want to help, please raise an issue on github.
Types Methods and Functions
WaterLily
WaterLily.AbstractBody
WaterLily.AbstractPoisson
WaterLily.AutoBody
WaterLily.Bodies
WaterLily.Flow
WaterLily.MultiLevelPoisson
WaterLily.NoBody
WaterLily.Simulation
WaterLily.BC!
WaterLily.BC!
WaterLily.BCTuple
WaterLily.CIj
WaterLily.Jacobi!
WaterLily.L₂
WaterLily.accelerate!
WaterLily.apply!
WaterLily.curl
WaterLily.curvature
WaterLily.exitBC!
WaterLily.inside
WaterLily.inside_u
WaterLily.interp
WaterLily.ke
WaterLily.loc
WaterLily.measure
WaterLily.measure!
WaterLily.measure!
WaterLily.measure_sdf!
WaterLily.mom_step!
WaterLily.mult!
WaterLily.pcg!
WaterLily.reduce_sdf_map
WaterLily.residual!
WaterLily.sdf
WaterLily.sdf
WaterLily.sdf_map_d
WaterLily.sim_step!
WaterLily.sim_time
WaterLily.slice
WaterLily.solver!
WaterLily.δ
WaterLily.λ₂
WaterLily.ω
WaterLily.ω_mag
WaterLily.ω_θ
WaterLily.∂
WaterLily.∮nds
WaterLily.@inside
WaterLily.@loop
WaterLily.AbstractBody
— TypeAbstractBody
Immersed body Abstract Type. Any AbstractBody
subtype must implement
d = sdf(body::AbstractBody, x, t=0)
and
d,n,V = measure(body::AbstractBody, x, t=0)
where d
is the signed distance from x
to the body at time t
, and n
& V
are the normal and velocity vectors implied at x
.
WaterLily.AbstractPoisson
— TypePoisson{N,M}
Composite type for conservative variable coefficient Poisson equations:
∮ds β ∂x/∂n = σ
The resulting linear system is
Ax = [L+D+L']x = z
where A is symmetric, block-tridiagonal and extremely sparse. Moreover, D[I]=-∑ᵢ(L[I,i]+L'[I,i])
. This means matrix storage, multiplication, ect can be easily implemented and optimized without external libraries.
To help iteratively solve the system above, the Poisson structure holds helper arrays for inv(D)
, the error ϵ
, and residual r=z-Ax
. An iterative solution method then estimates the error ϵ=̃A⁻¹r
and increments x+=ϵ
, r-=Aϵ
.
WaterLily.AutoBody
— TypeAutoBody(sdf,map=(x,t)->x; compose=true) <: AbstractBody
sdf(x::AbstractVector,t::Real)::Real
: signed distance functionmap(x::AbstractVector,t::Real)::AbstractVector
: coordinate mapping functioncompose::Bool=true
: Flag for composingsdf=sdf∘map
Implicitly define a geometry by its sdf
and optional coordinate map
. Note: the map
is composed automatically if compose=true
, i.e. sdf(x,t) = sdf(map(x,t),t)
. Both parameters remain independent otherwise. It can be particularly heplful to set compose=false
when adding mulitple bodies together to create a more complex one.
WaterLily.Bodies
— TypeBodies(bodies, ops::AbstractVector)
bodies::Vector{AutoBody}
: Vector ofAutoBody
ops::Vector{Function}
: Vector of operators for the superposition of multipleAutoBody
s
Superposes multiple body::AutoBody
objects together according to the operators ops
. While this can be manually performed by the operators implemented for AutoBody
, adding too many bodies can yield a recursion problem of the sdf
and map
functions not fitting in the stack. This type implements the superposition of bodies by iteration instead of recursion, and the reduction of the sdf
and map
functions is done on the mesure
function, and not before. The operators vector ops
specifies the operation to call between two consecutive AutoBody
s in the bodies
vector. Note that +
(or the alias ∪
) is the only operation supported between Bodies
.
WaterLily.Flow
— TypeFlow{D::Int, T::Float, Sf<:AbstractArray{T,D}, Vf<:AbstractArray{T,D+1}, Tf<:AbstractArray{T,D+2}}
Composite type for a multidimensional immersed boundary flow simulation.
Flow solves the unsteady incompressible Navier-Stokes equations on a Cartesian grid. Solid boundaries are modelled using the Boundary Data Immersion Method. The primary variables are the scalar pressure p
(an array of dimension D
) and the velocity vector field u
(an array of dimension D+1
).
WaterLily.MultiLevelPoisson
— TypeMultiLevelPoisson{N,M}
Composite type used to solve the pressure Poisson equation with a geometric multigrid method. The only variable is levels
, a vector of nested Poisson
systems.
WaterLily.NoBody
— TypeNoBody
Use for a simulation without a body.
WaterLily.Simulation
— TypeSimulation(dims::NTuple, u_BC::Union{NTuple,Function}, L::Number;
U=norm2(u_BC), Δt=0.25, ν=0., ϵ=1, perdir=(1,)
uλ::nothing, g=nothing, exitBC=false,
body::AbstractBody=NoBody(),
T=Float32, mem=Array)
Constructor for a WaterLily.jl simulation:
dims
: Simulation domain dimensions.u_BC
: Simulation domain velocity boundary conditions, either a tupleu_BC[i]=uᵢ, i=eachindex(dims)
, or a time-varying functionf(i,t)
L
: Simulation length scale.U
: Simulation velocity scale.Δt
: Initial time step.ν
: Scaled viscosity (Re=UL/ν
).g
: Domain acceleration,g(i,t)=duᵢ/dt
ϵ
: BDIM kernel width.perdir
: Domain periodic boundary condition in the(i,)
direction.exitBC
: Convective exit boundary condition in thei=1
direction.uλ
: Function to generate the initial velocity field.body
: Immersed geometry.T
: Array element type.mem
: memory location.Array
,CuArray
,ROCm
to run on CPU, NVIDIA, or AMD devices, respectively.
See files in examples
folder for examples.
WaterLily.BC!
— FunctionBC!(a,A)
Apply boundary conditions to the ghost cells of a vector field. A Dirichlet condition a[I,i]=A[i]
is applied to the vector component normal to the domain boundary. For example aₓ(x)=Aₓ ∀ x ∈ minmax(X)
. A zero Neumann condition is applied to the tangential components.
WaterLily.BC!
— MethodBC!(a)
Apply zero Neumann boundary conditions to the ghost cells of a scalar field.
WaterLily.BCTuple
— MethodBCTuple(f,t,N)
Generate a tuple of N
values from either a boundary condition function f(i,t)
or the tuple of boundary conditions f=(fₓ,...).
WaterLily.CIj
— MethodCIj(j,I,jj)
Replace jᵗʰ component of CartesianIndex with k
WaterLily.Jacobi!
— MethodJacobi!(p::Poisson; it=1)
Jacobi smoother run it
times. Note: This runs for general backends, but is very slow to converge.
WaterLily.L₂
— MethodL₂(a)
L₂ norm of array a
excluding ghosts.
WaterLily.accelerate!
— Methodaccelerate!(r,t,g)
This function adds a uniform acceleration field g
at time t
to r
. If g ≠ nothing
, then g(i,t)=dUᵢ/dt
.
WaterLily.apply!
— Methodapply!(f, c)
Apply a vector function f(i,x)
to the faces of a uniform staggered array c
or a function f(x)
to the center of a uniform array c
.
WaterLily.curl
— Methodcurl(i,I,u)
Compute component i
of $𝛁×𝐮$ at the edge of cell I
. For example curl(3,CartesianIndex(2,2,2),u)
will compute ω₃(x=1.5,y=1.5,z=2)
as this edge produces the highest accuracy for this mix of cross derivatives on a staggered grid.
WaterLily.curvature
— Methodcurvature(A::AbstractMatrix)
Return H,K
the mean and Gaussian curvature from A=hessian(sdf)
. K=tr(minor(A))
in 3D and K=0
in 2D.
WaterLily.exitBC!
— MethodexitBC!(u,u⁰,U,Δt)
Apply a 1D convection scheme to fill the ghost cell on the exit of the domain.
WaterLily.inside
— Methodinside(a)
Return CartesianIndices range excluding a single layer of cells on all boundaries.
WaterLily.inside_u
— Methodinside_u(dims,j)
Return CartesianIndices range excluding the ghost-cells on the boundaries of a vector array on face j
with size dims
.
WaterLily.interp
— Methodinterp(x::SVector, arr::AbstractArray)
Linear interpolation from array `arr` at index-coordinate `x`.
Note: This routine works for any number of dimensions.
WaterLily.ke
— Methodke(I::CartesianIndex,u,U=0)
Compute $½∥𝐮-𝐔∥²$ at center of cell I
where U
can be used to subtract a background flow (by default, U=0
).
WaterLily.loc
— Methodloc(i,I) = loc(Ii)
Location in space of the cell at CartesianIndex I
at face i
. Using i=0
returns the cell center s.t. loc = I
.
WaterLily.measure!
— Functionmeasure!(sim::Simulation,t=timeNext(sim))
Measure a dynamic body
to update the flow
and pois
coefficients.
WaterLily.measure!
— Methodmeasure!(flow::Flow, body::AbstractBody; t=0, ϵ=1)
Queries the body geometry to fill the arrays:
flow.μ₀
, Zeroth kernel momentflow.μ₁
, First kernel moment scaled by the body normalflow.V
, Body velocity
at time t
using an immersion kernel of size ϵ
.
See Maertens & Weymouth, doi:10.1016/j.cma.2014.09.007.
WaterLily.measure
— Methodd,n,V = measure(body::AutoBody,x,t)
d,n,V = measure(body::Bodies,x,t)
Determine the implicit geometric properties from the sdf
and map
. The gradient of d=sdf(map(x,t))
is used to improve d
for pseudo-sdfs. The velocity is determined solely from the optional map
function.
WaterLily.measure_sdf!
— Functionmeasure_sdf!(a::AbstractArray, body::AbstractBody, t=0)
Uses sdf(body,x,t)
to fill a
.
WaterLily.mom_step!
— Methodmom_step!(a::Flow,b::AbstractPoisson)
Integrate the Flow
one time step using the Boundary Data Immersion Method and the AbstractPoisson
pressure solver to project the velocity onto an incompressible flow.
WaterLily.mult!
— Methodmult!(p::Poisson,x)
Efficient function for Poisson matrix-vector multiplication. Fills p.z = p.A x
with 0 in the ghost cells.
WaterLily.pcg!
— Methodpcg!(p::Poisson; it=6)
Conjugate-Gradient smoother with Jacobi preditioning. Runs at most it
iterations, but will exit early if the Gram-Schmidt update parameter |α| < 1%
or |r D⁻¹ r| < 1e-8
. Note: This runs for general backends and is the default smoother.
WaterLily.reduce_sdf_map
— Methodreduce_sdf_map(sdf_a,map_a,d_a,sdf_b,map_b,d_b,op,x,t)
Reduces two different sdf
and map
functions, and d
value.
WaterLily.residual!
— Methodresidual!(p::Poisson)
Computes the resiual r = z-Ax
and corrects it such that r = 0
if iD==0
which ensures local satisfiability and sum(r) = 0
which ensures global satisfiability.
The global correction is done by adjusting all points uniformly, minimizing the local effect. Other approaches are possible.
Note: These corrections mean x
is not strictly solving Ax=z
, but without the corrections, no solution exists.
WaterLily.sdf
— Methodd = sdf(body::AutoBody,x,t) = body.sdf(x,t)
WaterLily.sdf
— Methodd = sdf(a::Bodies,x,t)
Computes distance for Bodies
type.
WaterLily.sdf_map_d
— Methodsdf_map_d(ab::Bodies,x,t)
Returns the sdf
and map
functions, and the distance d
(d=sdf(x,t)
) for the Bodies
type.
WaterLily.sim_step!
— Methodsim_step!(sim::Simulation,t_end=sim(time)+Δt;max_steps=typemax(Int),remeasure=true,verbose=false)
Integrate the simulation sim
up to dimensionless time t_end
. If remeasure=true
, the body is remeasured at every time step. Can be set to false
for static geometries to speed up simulation.
WaterLily.sim_time
— Methodsim_time(sim::Simulation)
Return the current dimensionless time of the simulation tU/L
where t=sum(Δt)
, and U
,L
are the simulation velocity and length scales.
WaterLily.slice
— Methodslice(dims,i,j,low=1)
Return CartesianIndices
range slicing through an array of size dims
in dimension j
at index i
. low
optionally sets the lower extent of the range in the other dimensions.
WaterLily.solver!
— Methodsolver!(A::Poisson;log,tol,itmx)
Approximate iterative solver for the Poisson matrix equation Ax=b
.
A
: Poisson matrix with working arrays.A.x
: Solution vector. Can start with an initial guess.A.z
: Right-Hand-Side vector. Will be overwritten!A.n[end]
: stores the number of iterations performed.log
: Iftrue
, this function returns a vector holding theL₂
-norm of the residual at each iteration.tol
: Convergence tolerance on theL₂
-norm residual.itmx
: Maximum number of iterations.
WaterLily.δ
— Methodδ(i,N::Int)
δ(i,I::CartesianIndex{N}) where {N}
Return a CartesianIndex of dimension N
which is one at index i
and zero elsewhere.
WaterLily.λ₂
— Methodλ₂(I::CartesianIndex{3},u)
λ₂ is a deformation tensor metric to identify vortex cores. See https://en.wikipedia.org/wiki/Lambda2_method and Jeong, J., & Hussain, F., doi:10.1017/S0022112095000462
WaterLily.ω
— Methodω(I::CartesianIndex{3},u)
Compute 3-vector $𝛚=𝛁×𝐮$ at the center of cell I
.
WaterLily.ω_mag
— Methodω_mag(I::CartesianIndex{3},u)
Compute $∥𝛚∥$ at the center of cell I
.
WaterLily.ω_θ
— Methodω_θ(I::CartesianIndex{3},z,center,u)
Compute $𝛚⋅𝛉$ at the center of cell I
where $𝛉$ is the azimuth direction around vector z
passing through center
.
WaterLily.∂
— Method∂(i,j,I,u)
Compute $∂uᵢ/∂xⱼ$ at center of cell I
. Cross terms are computed less accurately than inline terms because of the staggered grid.
WaterLily.∮nds
— Method∮nds(p,body::AutoBody,t=0)
Surface normal integral of field p
over the body
.
WaterLily.@inside
— Macro@inside <expr>
Simple macro to automate efficient loops over cells excluding ghosts. For example,
@inside p[I] = sum(loc(0,I))
becomes
@loop p[I] = sum(loc(0,I)) over I ∈ inside(p)
See @loop
.
WaterLily.@loop
— Macro@loop <expr> over <I ∈ R>
Macro to automate fast CPU or GPU loops using KernelAbstractions.jl. The macro creates a kernel function from the expression <expr>
and evaluates that function over the CartesianIndices I ∈ R
.
For example
@loop a[I,i] += sum(loc(i,I)) over I ∈ R
becomes
@kernel function kern(a,i,@Const(I0))
I ∈ @index(Global,Cartesian)+I0
a[I,i] += sum(loc(i,I))
end
kern(get_backend(a),64)(a,i,R[1]-oneunit(R[1]),ndrange=size(R))
where get_backend
is used on the first variable in expr
(a
in this example).