Added grid render.
Added particle toggle. Added boundary shape adjustment.
This commit is contained in:
171
fluid/fluid.go
171
fluid/fluid.go
@@ -4,6 +4,7 @@
|
||||
package fluid
|
||||
|
||||
import (
|
||||
"fluids/utils"
|
||||
"math"
|
||||
)
|
||||
|
||||
@@ -30,15 +31,16 @@ type Fluid struct {
|
||||
pNumX, pNumY int //SPATIAL HASHING: number of particle cells x,y
|
||||
pNumCells int //SPATIAL HASHING: number of particle cells in total
|
||||
pInvSpacing float32 //SPATIAL HASHING: one over the particle spacing
|
||||
particleDensity []float32 //density of particles within each field cell
|
||||
particleRestDensity float32
|
||||
ParticleDensity []float32 //density of particles within each field cell
|
||||
ParticleRestDensity float32
|
||||
fluidDensity float32
|
||||
angle float32
|
||||
dt float32 //delta time step
|
||||
gravitynormal FieldVector //unit vector for gravity
|
||||
stepacceleration FieldVector //step acceleration due to gravity
|
||||
flipPicRatio float32
|
||||
numSubSteps int //number of simulation substeps
|
||||
numSubSteps int //number of simulation substeps
|
||||
Block bool //rectangular or circular field container
|
||||
}
|
||||
|
||||
func NewFluid(dimensions FieldVector, spacing float32) *Fluid {
|
||||
@@ -46,9 +48,10 @@ func NewFluid(dimensions FieldVector, spacing float32) *Fluid {
|
||||
Field: NewVelocityField(dimensions, spacing),
|
||||
dt: FluidDefaultTimeStep,
|
||||
fluidDensity: FluidDefaultDensity,
|
||||
particleRestDensity: FluidDefaultRestDensity,
|
||||
ParticleRestDensity: FluidDefaultRestDensity,
|
||||
flipPicRatio: FluidDefaultFlipPicRatio,
|
||||
numSubSteps: FluidDefaultSubSteps,
|
||||
Block: true,
|
||||
}
|
||||
|
||||
f.Initialize()
|
||||
@@ -59,7 +62,7 @@ func (f *Fluid) Initialize() {
|
||||
f.InitializeParticles()
|
||||
|
||||
f.SetAngle(FluidDefaultAngle)
|
||||
f.particleDensity = make([]float32, f.Field.Nx*f.Field.Ny)
|
||||
f.ParticleDensity = make([]float32, f.Field.Nx*f.Field.Ny)
|
||||
}
|
||||
|
||||
func (f *Fluid) InitializeParticles() {
|
||||
@@ -158,32 +161,6 @@ func (f *Fluid) IntegrateParticles() {
|
||||
}
|
||||
}
|
||||
|
||||
// returns adjusted value between start and end, whichever is closer
|
||||
// unless value is between, in which case it returns value
|
||||
func clamp(value, start, end int) int {
|
||||
if value < start {
|
||||
return start
|
||||
}
|
||||
|
||||
if value < end {
|
||||
return value
|
||||
}
|
||||
|
||||
return end
|
||||
}
|
||||
|
||||
func clamp32(value, start, end float32) float32 {
|
||||
if value < start {
|
||||
return start
|
||||
}
|
||||
|
||||
if value < end {
|
||||
return value
|
||||
}
|
||||
|
||||
return end
|
||||
}
|
||||
|
||||
func (f *Fluid) PerformParticleSpatialHash() {
|
||||
|
||||
clear(f.numCellParticles)
|
||||
@@ -191,8 +168,8 @@ func (f *Fluid) PerformParticleSpatialHash() {
|
||||
//find spatcial hash cell for this particle's position, increment count there
|
||||
for i := range f.Particles {
|
||||
p := &f.Particles[i]
|
||||
xi := clamp(int(p.Position.X*f.pInvSpacing), 0, f.pNumX-1)
|
||||
yi := clamp(int(p.Position.Y*f.pInvSpacing), 0, f.pNumY-1)
|
||||
xi := utils.Clamp(int(p.Position.X*f.pInvSpacing), 0, f.pNumX-1)
|
||||
yi := utils.Clamp(int(p.Position.Y*f.pInvSpacing), 0, f.pNumY-1)
|
||||
cellnumber := xi*f.pNumY + yi
|
||||
f.numCellParticles[cellnumber]++
|
||||
}
|
||||
@@ -209,8 +186,8 @@ func (f *Fluid) PerformParticleSpatialHash() {
|
||||
//firstCellParticle bucket as we go
|
||||
for i := range f.Particles {
|
||||
p := &f.Particles[i]
|
||||
xi := clamp(int(p.Position.X*f.pInvSpacing), 0, f.pNumX-1)
|
||||
yi := clamp(int(p.Position.Y*f.pInvSpacing), 0, f.pNumY-1)
|
||||
xi := utils.Clamp(int(p.Position.X*f.pInvSpacing), 0, f.pNumX-1)
|
||||
yi := utils.Clamp(int(p.Position.Y*f.pInvSpacing), 0, f.pNumY-1)
|
||||
cellnumber := xi*f.pNumY + yi
|
||||
f.firstCellParticle[cellnumber]--
|
||||
f.cellParticleIds[f.firstCellParticle[cellnumber]] = i
|
||||
@@ -291,32 +268,78 @@ func (f *Fluid) HandleParticleCollisions() {
|
||||
minDist2 := minDist * minDist
|
||||
*/
|
||||
|
||||
minX := f.particleRadius
|
||||
maxX := float32(f.Field.Nx-1)*f.Field.H - f.particleRadius
|
||||
minY := f.particleRadius
|
||||
maxY := float32(f.Field.Ny-1)*f.Field.H - f.particleRadius
|
||||
if f.Block {
|
||||
minX := f.Field.H + f.particleRadius
|
||||
maxX := float32(f.Field.Nx-1)*f.Field.H - f.particleRadius
|
||||
minY := f.Field.H + f.particleRadius
|
||||
maxY := float32(f.Field.Ny-1)*f.Field.H - f.particleRadius
|
||||
|
||||
for i := range f.Particles {
|
||||
p := &f.Particles[i]
|
||||
for i := range f.Particles {
|
||||
p := &f.Particles[i]
|
||||
|
||||
if p.Position.X < minX {
|
||||
p.Position.X = minX
|
||||
if p.Position.X < minX {
|
||||
p.Position.X = minX
|
||||
p.Velocity.X = 0
|
||||
}
|
||||
|
||||
if p.Position.X > maxX {
|
||||
p.Position.X = maxX
|
||||
p.Velocity.X = 0
|
||||
}
|
||||
|
||||
if p.Position.Y < minY {
|
||||
p.Position.Y = minY
|
||||
p.Velocity.Y = 0
|
||||
}
|
||||
|
||||
if p.Position.Y > maxY {
|
||||
p.Position.Y = maxY
|
||||
p.Velocity.Y = 0
|
||||
}
|
||||
|
||||
}
|
||||
} else {
|
||||
//find fluid cell in which the particle resides
|
||||
//compute distance of this cell from 'center' cell of our circle
|
||||
//if distance exceeds radius, clamp particle velocity and position
|
||||
//to edge cell
|
||||
|
||||
//find origin of circle, conveniently these dimensions are also our radius
|
||||
originX := float32(f.Field.Nx) * f.Field.H / 2
|
||||
originY := float32(f.Field.Ny) * f.Field.H / 2
|
||||
radius := originX - f.Field.H - f.particleRadius
|
||||
|
||||
//cxi := utils.Clamp(int(math.Floor(float64(originX*f.Field.InvH))), 0, f.Field.Nx-1)
|
||||
//cyi := utils.Clamp(int(math.Floor(float64(originY*f.Field.InvH))), 0, f.Field.Ny-1)
|
||||
//originCellIdx := cxi*f.Field.Ny + cyi
|
||||
|
||||
//find fluid cell for particle
|
||||
for i := range f.Particles {
|
||||
p := &f.Particles[i]
|
||||
|
||||
//xi := utils.Clamp(int(math.Floor(float64(p.Position.X*f.Field.InvH))), 0, f.Field.Nx-1)
|
||||
//yi := utils.Clamp(int(math.Floor(float64(p.Position.Y*f.Field.InvH))), 0, f.Field.Ny-1)
|
||||
//
|
||||
//dx := originX - px
|
||||
//dy := originY - py
|
||||
|
||||
dx := p.Position.X - originX
|
||||
dy := p.Position.Y - originY
|
||||
dist2 := dx*dx + dy*dy
|
||||
if dist2 < radius*radius || dist2 == 0 {
|
||||
continue
|
||||
}
|
||||
|
||||
dist := float32(math.Sqrt(float64(dist2)))
|
||||
|
||||
newx := originX + dx*radius/dist
|
||||
newy := originY + dy*radius/dist
|
||||
|
||||
p.Position.X = newx
|
||||
p.Position.Y = newy
|
||||
p.Velocity.X = 0
|
||||
}
|
||||
|
||||
if p.Position.X > maxX {
|
||||
p.Position.X = maxX
|
||||
p.Velocity.X = 0
|
||||
}
|
||||
|
||||
if p.Position.Y < minY {
|
||||
p.Position.Y = minY
|
||||
p.Velocity.Y = 0
|
||||
}
|
||||
|
||||
if p.Position.Y > maxY {
|
||||
p.Position.Y = maxY
|
||||
p.Velocity.Y = 0
|
||||
}
|
||||
|
||||
}
|
||||
@@ -348,8 +371,8 @@ func (f *Fluid) TransferVelocities(togrid bool) {
|
||||
for i := range f.Particles {
|
||||
x := f.Particles[i].Position.X
|
||||
y := f.Particles[i].Position.Y
|
||||
xi := clamp(int(math.Floor(float64(x*h1))), 0, f.Field.Nx-1)
|
||||
yi := clamp(int(math.Floor(float64(y*h1))), 0, f.Field.Ny-1)
|
||||
xi := utils.Clamp(int(math.Floor(float64(x*h1))), 0, f.Field.Nx-1)
|
||||
yi := utils.Clamp(int(math.Floor(float64(y*h1))), 0, f.Field.Ny-1)
|
||||
cellnumber := xi*n + yi
|
||||
if f.Field.CellType[cellnumber] == CellTypeAir {
|
||||
f.Field.CellType[cellnumber] = CellTypeFluid
|
||||
@@ -380,8 +403,8 @@ func (f *Fluid) TransferVelocities(togrid bool) {
|
||||
x := p.Position.X
|
||||
y := p.Position.Y
|
||||
|
||||
x = clamp32(x, f.Field.H, float32((f.Field.Nx-1))*f.Field.H)
|
||||
y = clamp32(y, f.Field.H, float32((f.Field.Ny-1))*f.Field.H)
|
||||
x = utils.Clamp32(x, f.Field.H, float32((f.Field.Nx-1))*f.Field.H)
|
||||
y = utils.Clamp32(y, f.Field.H, float32((f.Field.Ny-1))*f.Field.H)
|
||||
|
||||
//find cell components neighbouring x,y
|
||||
x0 := min(int(math.Floor(float64((x-dx)*h1))), f.Field.Nx-2)
|
||||
@@ -513,15 +536,15 @@ func (f *Fluid) UpdateParticleDensity() {
|
||||
h1 := f.Field.InvH
|
||||
h2 := 0.5 * h
|
||||
|
||||
clear(f.particleDensity)
|
||||
clear(f.ParticleDensity)
|
||||
|
||||
for i := range f.Particles {
|
||||
p := &f.Particles[i]
|
||||
x := p.Position.X
|
||||
y := p.Position.Y
|
||||
|
||||
x = clamp32(x, h, float32(f.Field.Nx-1)*h)
|
||||
y = clamp32(y, h, float32(f.Field.Ny-1)*h)
|
||||
x = utils.Clamp32(x, h, float32(f.Field.Nx-1)*h)
|
||||
y = utils.Clamp32(y, h, float32(f.Field.Ny-1)*h)
|
||||
|
||||
x0 := int(math.Floor(float64((x - h2) * h1)))
|
||||
x1 := min(x0+1, f.Field.Nx-2)
|
||||
@@ -536,32 +559,32 @@ func (f *Fluid) UpdateParticleDensity() {
|
||||
sy := 1.0 - ty
|
||||
|
||||
if x0 < f.Field.Nx && y0 < f.Field.Ny {
|
||||
f.particleDensity[x0*n+y0] += sx * sy
|
||||
f.ParticleDensity[x0*n+y0] += sx * sy
|
||||
}
|
||||
if x1 < f.Field.Nx && y0 < f.Field.Ny {
|
||||
f.particleDensity[x1*n+y0] += tx * sy
|
||||
f.ParticleDensity[x1*n+y0] += tx * sy
|
||||
}
|
||||
if x1 < f.Field.Nx && y1 < f.Field.Ny {
|
||||
f.particleDensity[x1*n+y1] += tx * ty
|
||||
f.ParticleDensity[x1*n+y1] += tx * ty
|
||||
}
|
||||
if x0 < f.Field.Nx && y1 < f.Field.Ny {
|
||||
f.particleDensity[x0*n+y1] += sx * ty
|
||||
f.ParticleDensity[x0*n+y1] += sx * ty
|
||||
}
|
||||
}
|
||||
|
||||
if f.particleRestDensity == 0.0 {
|
||||
if f.ParticleRestDensity == 0.0 {
|
||||
var sum float32 = 0.0
|
||||
numFluidCells := 0
|
||||
|
||||
for i := 0; i < f.Field.Nx*f.Field.Ny; i++ {
|
||||
if f.Field.CellType[i] == CellTypeFluid {
|
||||
sum += f.particleDensity[i]
|
||||
sum += f.ParticleDensity[i]
|
||||
numFluidCells++
|
||||
}
|
||||
}
|
||||
|
||||
if numFluidCells > 0 {
|
||||
f.particleRestDensity = sum / float32(numFluidCells)
|
||||
f.ParticleRestDensity = sum / float32(numFluidCells)
|
||||
}
|
||||
}
|
||||
}
|
||||
@@ -602,9 +625,9 @@ func (f *Fluid) SolveIncompressibility() {
|
||||
//divergence -> we want this to be zero
|
||||
div := f.Field.U[right] - f.Field.U[center] + f.Field.V[top] - f.Field.V[center]
|
||||
|
||||
if f.particleRestDensity > 0.0 {
|
||||
if f.ParticleRestDensity > 0.0 {
|
||||
var k float32 = 1.0
|
||||
var compression float32 = f.particleDensity[i*n+j] - f.particleRestDensity
|
||||
var compression float32 = f.ParticleDensity[i*n+j] - f.ParticleRestDensity
|
||||
if compression > 0.0 {
|
||||
div = div - k*compression
|
||||
}
|
||||
@@ -624,3 +647,7 @@ func (f *Fluid) SolveIncompressibility() {
|
||||
}
|
||||
|
||||
}
|
||||
|
||||
func (f *Fluid) ToggleShape() {
|
||||
f.Block = !f.Block
|
||||
}
|
||||
|
||||
Reference in New Issue
Block a user