(whistle whistles)
- Hello and welcome to a Coding Challenge,
a calm, soothing, although somewhat turbulent,
Coding Challenge, called Fluid Simulation.
Now I have something to admit to you.
I don't really understand how any of this stuff works.
I did make this at one point, and so I'm hoping
this Coding Challenge to recreate exactly this,
which will be a basis on which hopefully a lot
of other interesting ideas will come about.
This idea came into mind when I recently saw
Smarter Every Day's video on laminar flow.
You know, I love laminar flow and all,
but Team Turbulence for life, muah.
Alright, so there's a wonderful 3Blue1Brown,
do I just reference the same other
YouTube channels every single time?
Yes, I do.
But there's also an excellent 3Blue1Brown video
on turbulence, which I would also recommend.
So, let me give you some background here.
So first of all, there is a GitHub issue thread
which started by deardanielxd, from 2016.
Lattice Boltzmann methods for fluid simulations.
So this is one method.
But what I want to highlight here is that this,
people have, oh, what I want to
highlight here are these three links.
So these seminal kind of canonical standard,
or the origins of doing fluid dynamics in computer graphics
in my research comes from this article by Jos Stam.
Real-Time Fluid Dynamics for Games.
I believe this was a Siggraph paper, 2003.
Somebody fact check me on that.
And it's built on top of this idea of these
Navier-Stokes Equations, which are partial differential
equations that describe fluid dynamics
and there was actually a $1 million challenge for proving
that this can or cannot be solved in three dimensions.
None of this is anything I'm capable of doing.
But this paper includes some of the formulas,
includes a lot of the code, and you can see one thing that's
sort of key concept here is a fluid simulation can be done
by thinking about fluid as kind of
particles that live in a grid.
And obviously it might be like an infinitely small grid
(laughs) in real life, but we can make that discrete.
Think about the grid of pixels and what the sort of density
or the velocity of the fluid is at every one of these spots
on the grid, that's ultimately what I'm going to do.
And there's some nice C code
and descriptions of some of these algorithms.
So, this article I believe serves as the basis
for Memo Atken's processing library,
called MSAFluid, and is also an open frameworks library,
which is a C++ engine.
And you can see here the way
that this ends up looking by sort of distorting
this vector field, ah, this is awesome.
Oh, I'm sure the YouTube compression
is totally ruining this.
But it's beautiful, check out that library.
You should also check out Gabriel Weymouth's,
I hope I pronounced that correctly, LilyPad project,
which was actually I believe used
in the 3Blue1Brown video on turbulence.
Gabriel writes here at the end a bunch of things
about Stan's approach, and what LilyPad does and a paper.
So this is a giant rabbit hole you could go down.
And I have spent some time poking around
in this rabbit hole in the last week.
The article that I found that I kind of enjoyed the most
in terms of style was Mike Ash's article called
Fluid Simulation for Dummies, which is actually a port,
not a port, but a version of Jos Stam's paper,
but actually turning it into 3-D and how to render that 3-D
with paralyzing computing power
is all here in his master's thesis.
That's another rabbit hole you could go down.
So what I would like to do is use this article,
and I try to make sure that I don't
repurpose somebody else's content without permission,
even if it's sort of on the web in an open-source way.
What I try to do, what I want to do,
I asked Mike Ash for permission on Twitter.
I think it's @MikeAsh, but I'll include a link in this
video's description, if I could go through this in a video.
So what I'm going to do is I'm going to
go through this in the video and mostly
just kind of like copy-paste this code,
which is written, I believe, C++ or C.
Something like that, some object-oriented C-flavor language.
And I'm going to copy-paste it into Processing,
which is a Java-based programming environment.
And kind of like adjust the code to work
in the way that I know and see if I can
get the result and play with the result.
So I don't feel obligated to understand or explain
all of the maths involved throughout all this.
And I would definitely encourage you to read this,
I would entirely stop and read this whole article
before you continue watching this.
And then, this is going to be in three parts at a minimum.
So this first part, I hope you just like get it working.
Just want to copy-paste the code, change the syntax around,
get it rendering, and play with it.
Number two is I want to kind of refactor the code,
this'll be in another video.
I'm actually going to refactor it in a more
kind of modern approach, I think modern's the wrong word,
I mean this was done over 10 years ago.
But I want to use like object-oriented programming,
vector, like the PVector class in processing,
I think there're some ways that I can redo the code
to make it a little bit more readable than this particular
style that uses a lot of esoteric variable meaning.
So that's number two, and then number three,
I want to then apply this logic to my,
my flow field example from the Nature of Code book.
So if I could take the fluid simulation,
turn it into a vector field,
I can just toss particles in that fly around.
I think some visual opportunities will come of that.
So that's three part.
Just get the thing working, that's what I'm about to attempt
to do right now, I'm sure it will go wrong. (laughs)
But I will try my best.
Number two, refactor the code to make it sort of fit
with how I think about coding and processing in P5GS today.
And then also, try to do some more stuff with it visually.
And I would say one of the things is, there's going to be some
performance issues, I'm going to keep things low-resolution.
But you'll see a lot of the implementations of this,
use shaders or webGL, all sorts of fancy tricks
that I'm not going to get into,
but if you know about that stuff and can build
on top of whatever I'm doing, then fantastic.
Alright, so I'm about to get started coding,
but before I do that, I've written in advance a bunch
of the concepts that are involved in this implementation
that I want to make sure that I don't forget to mention.
The first one that I think is really important
is this idea of an incompressible fluid.
An incompressible fluid is a fluid that density
must remain constant throughout, like water.
So for example if you have water in a balloon,
and you squeeze that balloon, the water's got to like
come out, it's not compressible, whereas air is
actually compressible, its density can change.
So this, apparently, from the little research
that I've done, simplifies a lot of the stuff.
So this fluid simulation is going to work
only for this idea of an incompressible fluid.
I should also mention that the goal here
is not necessarily to with scientific accuracy simulate
true fluid mechanics, but rather to create the illusion
and feeling of that through some remote connection
to that actual scientific accuracy.
And I'm sure whatever I do will be less accurate than what
people before, based on how I know the way I make things.
Alright, but that's an important thing.
So the first thing that we need to consider
is that the fluid is going to live inside of a box.
And I think, the way the math is sort of
tuned in these examples I believe is
so that this box should really be a square.
And sometimes, I don't know why,
like a power of two is maybe a nice thing.
So maybe I'll start with 512 by 512
or maybe 256 by 256 just to make it low resolution.
So you can think of this as a grid of pixels.
And there is going to be
inside of this grid a velocity vector.
A velocity vector that points in a given direction.
So if the velocity vector in every spot in the grid is zero,
it's like completely still water.
If I put a velocity vector moving to the right,
it's like the water, that would be
like Laminar Flow, by the way.
Because the velocity, everything is smooth and perfect
and all moving in exactly the same direction,
as opposed to turbulence where
everything's kind of going crazy.
So that's one thing that's going to exist in here.
So that is something I should say, there is this idea
of the velocity field, or the vector field.
And that's going to have Xs and Ys.
In Mike Ash's blog post, it's actually all done in 3-D
(laughs) with X, Ys, and Zs, and I'm going to take out
that third dimension while I'm doing it.
But you should add it back in and see what happens there.
Now, the other thing is there's going
to be this idea of dye.
Oh, I wrote this over here already,
'cause I wanted to explain that.
There's this idea of dye.
And we're going to talk about the density of the dye,
which is, in other words, this vector field,
we wouldn't be able to see anything moving,
we wouldn't be able to under see the flow through
the fluid without like putting something in it.
So you could imagine sprinkling a little dye in it
and having it diffuse,
maybe advect, diffuse, advect, project,
all around the fluid.
But what I want to make clear is when we start
talking about density in this code example,
and this Mike Ash makes very specific,
makes a very specific point about.
When we're talking about the density of this dye,
which is like an extra thing we're adding,
just so we can visualize it.
We're also able to visualize it just as a vector field,
which I will do at some point,
but the dye is what's going to give it more this
like smoky-like quality by visualizing the amount
of dye as it moves throughout the fluid.
So this is the basic idea.
So the first thing I need to do is get an array to store all
of the Xs and Ys of the vector field and the amount of dye,
for every single one of these spots.
And in the example, it's done with three separate arrays,
an array of Xs, an arrays of Ys, an arrays of densities,
and it kind of works like a cellular automata simulation,
where I need the previous state and the next state,
all the velocity of the previous and the next velocity,
all the density and the next density.
So I'm actually going to need two.
So this is why the code gets really confusing,
because I need x, Y, and density, then I need
what's named as x zero, Y zero, and density zero.
I forget, this might be called like S
or something in the code, but I need
the sort of previous of all of these.
So there's a whole bunch of arrays.
And this is what I want to later refactor it and see
if I can just use Pvector, an object that stores a PVector
and a density value in the previous, all that stuff.
Alright, let's go back to
Mike Ash's page.
And I'm going to start with this.
This is what I was talking about.
So this C++ structure, I'm going to take it into processing.
I'm going to add setup.
I'm going to add draw.
I'm going to say, size 250.
Actually, you know what I'm going to do,
I'm going to create a variable called N.
This'll appear in the code,
which is kind of like the width and height of the square.
So it'll be five, so N is going to be
in this case, what'd I say, 256.
And then actually in processing, if you use a variable for,
if you use a variable for the dimensions
that you want to put in size, you've actually got to put that
in the settings function so I can say this.
So this should get me like a 256 by 256 window.
And then, so I'm going to save that.
And I'm going to call this like FluidSimAsh.
And Stam, FluidSimAshStam. (laughs)
And then I'm going to create a class.
And I'm going to call this, make a new tab,
I'm going to say class fluid.
And in that class, I want to have all this stuff.
All of these.
Now the thing that's different here in processing is
this star, what is that star, what does that star mean?
And that's actually because it's some
C-flavored language in the code in Mike Stam's.
That is a pointer, meaning it's pointing to an area,
a memory address on the computer where all
of the density values will be stored.
But what I really want this to be is just an array.
So I'm going to change all these to an array.
Which is the same thing, in Java,
this is now a pointer to an array.
You'll see in a second.
And then, and by the way, I'm going to take out all the Zs.
I should mention, so, these variables,
this isn't a diffusion amount, like when she talks about,
which is like a variable to control how the velocity
and sort of the vectors and the dye
diffuses throughout the fluid.
This is viscosity.
Viscosity is like the thickness of the fluid,
so playing with that can change the behavior as well.
DT is the time step.
In all of my physics simulations, I've always done it,
just have the time step of one.
But I think you need a smaller time step
to be able to get the simulation to behave
somewhat accurately, so that'll come up later.
And now I have density,
and I think this is previous density.
Velocity X, velocity Y,
previous velocity X, previous velocity Y.
Alright, so, now this is creating it,
so I basically want to do exactly this.
So this would be in the constructor of the fluid.
And basically, N I already declared somewhere else.
So I don't think I need size.
I'm going to say cube size, and these should all be this dot.
So this is referring to the actual variables
in the object itself.
So, and I don't need the Z.
And, I don't need Z.
And, now, this is like calloc is like a memory allocation,
'cause you're allocating a certain amount
of memory for all these fluid values.
But I just want an array that is size N times N.
So this is going to be this, and it's N times N times N,
because his is in three dimensions.
So I know I could do like a find and replace,
but this is like so crazy that I want to do this.
This is like very meditative for me.
Okay, and this receives a DT, three arguments,
when you create the fluid.
You create it with a time step,
a diffusion, and a viscosity.
So in here, I'd be saying something like, fluid,
fluid, and in setup I would now say, fluid equals
a new fluid, maybe like a time step of .1, and I'm just
going to make a density and viscosity of zero for right now.
But those values would get filled in, presumably, okay.
(bell dings)
We're gettin' there, we're movin' through the article.
Here we go.
This is described, okay, you need to be able
to destroy the thing, free all the memory.
Ahh, Java garbage collection.
I don't have to worry about that,
it will get cleaned up for me.
Add density.
Okay, this makes sense.
Now remember, add density is not talking about,
is not referring to the fluid itself,
it's referring to the dye that's going into the fluid.
So this is basically like an add dye.
I sort of feel like I want to rename that function to add dye.
But I'm going to just call it add density.
So, we're going to take this function,
I think this could be a function of,
that's part of the object.
What does it need?
We don't need a reference to the object itself,
we have that, we need a location
and an amount, that makes sense.
So we need a location and an amount.
So that's like the amount of dye
we're adding at this XY spot.
Now, here's the thing.
You're notice something, everywhere in the code
there's this like IX function, index equals IX.
Well, this is a two-dimensional, this is two-dimensional.
But you notice all the arrays I made were one-dimensional.
This is kind of a pretty typical thing to do.
But I need a way of going from X, comma Y,
to the single index that's a lookup into this grid.
So, I could write a function, you know, I think it's done
with like a macro or define or whatever in this.
But I'm going to write a function,
I'm going to name it, I'm going to call it IX.
And it gets an X and a Y, and it would just return,
return X plus Y times N.
So this is, oh, it needs to return an integer.
It's non-void.
So this basically says for any given XY,
give me the one-dimensional index.
And this formula's the same thing that I use in all
my image processing and pixel processing, it's a way
to get a 2-D location in a one-dimensional array, okay.
So now, I would be saying that we can go back to this code,
basically, and we can do this, but no Z,
and then we're basically going to say, hey,
this dot, what is it called?
VX.
Oh wait, no, that's velocity, sorry.
This dot, density,
add index,
add some amount.
This is really simple.
This is like a really simple function, just add some density
to this spot, this amount of density to this spot.
Then we can also do this add velocity,
which is basically the same thing.
But just with an amount X and an amount Y,
so I think we could probably just copy this.
And we could say, add velocity
at XY with an amount X and an amount Y.
See, I would prefer to use like PVectors for all this stuff.
And then we get the index.
And we say, VX, VX plus amount X,
and VY,
by the way, this is going to get (laughs)
portings on this code is going to get a lot worse soon enough.
(laughs) Okay, now, alright.
Ah, look at this, so here, now we can take a moment,
here are the three main operations.
Diffuse, project, advect.
Let's do them one at a time.
So we can read this, I mean it's useful
to read Mike Ash's description, put a drop of soy sauce,
soy sauce gives a way of thinking about dye.
And you'll notice that doesn't
stay still, but it spreads out.
So this happens, even if the water and the sauce
are perfectly still, it's called diffusion.
And so, the dye obviously diffuses, that makes sense,
but the velocity also diffuses, if some of it moves,
it's causing everything around it to move as well.
So that's a function, this is really going
to be the function that does the solving
of the Navier Stokes equation.
Let's do diffuse.
So I'm going to put, and this is not going to be an object
inside the class, and I'll explain why in a second.
So the first thing that I'm going to do
is I'm just going to try to port this.
So any time there is the pointer,
I want this to be an array.
Dif is a diffusion amount, DT, the number of iterations,
and I don't think I need N as an argument,
'cause I'm just using that as a global variable.
I also kind of personally would like
to just keep iterations as a global variable.
I think that'll make things a little bit simpler.
Where did I put that?
I'm actually going to put these in the fluid tab,
'cause that's mostly where I'm working.
Final int iterations, like we can just put that as one,
but let's leave that as like 10 right now.
So the idea here is that this function
knows how to diffuse any arbitrary array of numbers, X,
based on its previous values, X zero,
based on a diffusion amount and a time step.
But you'll notice that what it does is it
immediately calls another function called linear solve.
So, a linear equation is like something like this.
I don't know, two X plus Y minus three Z equals 10.
Right, this is a linear equation.
We got to have multiple linear equations.
And algorithms for solving the sort of set of solutions
to this equations, what are all the X and Ys and Zs that
make this true, is known as a process of linear solution.
There're different techniques, I was reading the comments,
in Mike Ash's paper, there's like this
Gauss Seidel technique, some people did his code using this,
I forget which technique his code is particularly using.
But it's a way of basically solving those linear equations,
for fluid dynamics,
within the space of this grid.
So, that's kind of all I want to say about it.
But it's needed for every single,
for this diffusion algorithm.
And it's really just a thing that like passes
all the values all around over and over and over again
with lots of iterations, sort of spread out.
Like a cellular automata,
to have the velocities or densities of neighbors
affect the other neighbors and so on and so forth.
So, I'm now going to go and just grab
the linear solve code from here, I'm going to grab this,
and you can sort of see what it's doing.
So I'm going to grab this function, I'm going to bring it in here,
by the way, this is kind of how I work.
Like I kind of want to understand this more.
But I feel like I need to just port the code
and play around with it, and then maybe I could do
some more research about what the equations actually are
and how it's working, but sometimes this kind of like
messing with the code get in your hands and the code
could sometimes help you understand the math later.
So I'm going to look at this, I'm going to take out static,
I'm going to move this bracket over here.
I don't need iterations in N.
Those are going to be the same for everything.
These should be arrays.
And the kind of goofy thing here is, I don't need N.
So I need I and J for X and Y, but I can take out M.
And so that loses the bracket there.
And then I have to look at what's going on here.
And so I don't need M.
And then I also don't need to add these two components.
So let's see.
So the idea here is that I'm looking at,
and this, I've done this before,
I'm having like a crazy deja vu.
What this is now doing is saying,
the new value of a particular cell is based
on a function of itself and all of its neighbors.
And you can see that here.
It is equal to its old value
plus some combination of its current values.
Alright, so now we're going back to the paper.
And we're going to look at the next function.
Project.
Okay, so this project function is really
tied to the idea of this incompressible fluid.
The amount of fluid in each box has to stay constant.
So the amount of fluid going in has to be
exactly equal to the amount of fluid going out.
So this is kind of like a cleanup stage
to put the set thing back into equilibrium.
So we're going to be doing this
for all the different velocity arrays.
So let's grab this, oh, look at this craziness.
So let's grab the project function.
Paste and project.
These are arrays.
I don't need the Z.
These are also arrays.
And I don't need the iterations, and N, those are constants.
And then now, I also don't need K.
Only need I and J, 'cause I'm just doing two dimensions.
I don't need K.
I don't need these two.
And then, I don't have set boundaries yet, I'll add that.
Linear solve is the same thing but without this.
And I think, when I do set bounds,
I'll probably take out the N.
Then I don't need K.
Don't need K, don't need K.
I don't need K, I don't need Z.
And I don't need K.
I don't need K, I guess probably porting Jos Stam's
which was 2-D, might have made more sense. (laughs)
And I plus one J, I minus one J.
J plus one, J minus one.
Yes, yes, this all makes sense.
And then I don't need the Z here.
I think I lost a curly bracket, which would go here.
And definitely screwed up something
curly bracket-wise.
I think I have an extra curly bracket here.
This looks right.
Oh, I don't need the,
I don't need the Ks up here.
Missing left curly bracket.
And I have an extra curly bracket here.
Alright, we're good.
We don't know what set bounds is, but that's okay, alright.
Woo.
We now, don't worry, we're getting close.
So, now.
Advection.
The advect step is responsible
for actually moving things around.
To that end, it looks at each cell in turn.
In that cell it grabs the velocity,
follows that velocity back in time, whoa.
So let's grab this advect function,
and I could kind of scroll back up.
Let's look at what was, yeah.
So this is also important, every cell has a set
of velocities and these velocities make things move.
This is called advection.
As with diffusion, advection applies
both to the dye and to the velocity itself.
So what's really the difference
between diffusion and advection?
Well diffusion is just this idea of spreading out.
But advection is actually
the motion associated with the velocities.
They're obviously related and they both
happen together, but those are separate things.
So let's go grab the advection code.
Oh, what? (laughs)
Seriously?
Okay, wow, so this I really going to want to unpack
when I refactor this to understand what's going on.
But I'm just going to grab it right now.
This is crazy.
But bring it in here, oh boy, ooh.
Okay, right, everybody,
deep breath.
Deep breath.
This is an array, this is the current density.
This is the previous density, this is the current velocity.
Oh.
This is the current velocity.
I don't need Z.
Come on, scroll over.
So, what, are you kidding me?
So these are indexed in previous index values.
Don't need the Z.
We need a DT for X and a DT for Y,
we don't need the Z.
This has to do with like,
I think the Ss are density in this.
I'm not really sure.
We'll look through the code, but I'm going to take out the U.
I don't need the three, and I don't need the Z.
Oh my god, I don't need the K.
I don't need the K.
Alright, so I don't need the K, part of the loop.
And I don't need this last thing.
And I don't need this Z.
And I need the X and the Is, the Ys and the Js,
and I don't need the Zs and the Ks.
And I need the Ss and the Ts.
Don't need the Us, really, and I don't need the Ks.
And now I don't need the K,
and so what did I get rid of up here?
I'm using, this is what I mean.
Like I really want to refactor this and try to understand
what each of these variables are doing and rename them.
(laughs) This isn't refactoring,
this is just getting into work.
So I've got Is and Js, Xs and Ys, Is and Js, Ss and Ts.
And no Ks.
Oh, wait a second.
STU.
U was the thing, if I go back to Mike Ash,
U was the third dimension for the S and the Ts.
Okay.
So what I'm actually getting rid of
is this multiplication step.
So.
One last multiplication step.
Think I don't need this at all.
This is what I'm doing, yeah.
This is much simpler.
Oh, you silly third dimension.
Pretty sure this is right. (laughs)
What's the chances there,
going to definitely have to double check this.
Double check this.
And then I, one last thing here.
Do I have the right, okay, floor F
is just floor,
right, because those are integers.
There seems to be some extra flooring
that's unnecessary here.
'Cause these are floats, but I'm not sure.
I'm going to leave it as is.
Right, these would then have to be converted.
I think there's some extra unnecessary steps here
in how the numbers are being manipulated.
But I'm going to just leave everything in, and I will,
part two, we're going to get this to work.
Just fast forward like five or 10 minutes towards the end.
Okay.
Alright.
I'm pretty sure this is now the advect function.
I think I actually have everything, but this set boundary.
So this set boundary is actually something like,
really nice, and compared to all that,
much more manageable to understand.
But basically, we've got this kind of
unrealistic situation here, well it's not unrealistic,
but we're containing the fluid within a box.
And so we need some sort of mechanism for when,
what we do with the edge cells, and how we deal
with those velocity vectors or things moving.
And basically, we want to add some kind of like bouncing.
And then we also need to deal with the corner cells
even differently, in a different way, so we can read,
this is short for set bounds, and it's a way
to keep your fluid from leaking out of your box.
Not that it could really leak, writes Mike Ash,
it's just a simulation in memory,
but if there's no walls, obviously fluid's
not going to leak out of your computer screen.
That was a joke there, and I just like,
my body hurt so much, I couldn't even laugh, ha ha.
But the velocity, it needs to come nearer
when it hits the wall.
So if we look at what, and the reason
why this is done in this like really weird way
is the set bounds function is written in a way
that this B variable tells you like which,
which wall am I at?
Am I at the right wall, the left wall,
the top wall, the bottom wall?
And I think that I could rewrite this with vectors
in a way that it just kind of like knows.
And I don't need this extra variable B
that's passed through.
But once again, I am going to,
I am going to port it, so exactly as written.
X is any one of these arrays, we need to set,
we need to do the same things for all of the arrays
where there's velocity X, velocity Y,
or density, density of the dye.
And once again, I can actually eliminate,
this is for the third dimension,
so I can completely take this out.
Then I'm going to take out K.
And you can see like this, if B equals two,
you're dealing with the X stuff.
I mean it's actually to the Y, 'cause you can see
it's looking at the top row and the bottom row, sorry.
So, this is good.
And then this
is if B is one, think I might've messed this up.
(otherworldly music)
(bell dings)
Alright, strange edit point, but I'm jumping in here
from the future in the middle of this video
to explain something about this set bound function that I
kind of botched while I was actually doing the challenge
and I'm at the end now, but I'm going to fix it.
One thing that's really important is for me to point out,
what is written in the article by Mike Ash?
Every velocity in the layer next
to the outer layer is nearer.
So when you have some velocity going towards the wall
in the spot next to the wall,
the wall gets the velocity that perfectly counters it.
That's keeping it a closed box.
So that is the explanation for why, here, I'm saying,
hey, if I'm at I comma zero, any given X all along
the top row, if Y is zero, all the Xs along the top
counteract the velocity from one row down.
But only if I've sent in this B equals two, which is that
weird thing, and this is what I think would be easier
with Pvector, 'cause I could do this all at once.
But, the thing that I have in here, I have this like
leftover X room, 'cause in 3-D there's more papers.
So I actually could take this out,
and this is going to make this now look correct.
So this set boundary function is actually much simpler.
It is just reversing the velocities in the last,
in all the edge layers according to, edge columns.
Edge columns are rows.
According to the next two, the edge column or row.
And then for the corners,
it actually just does an average of the two neighbors.
And that keeps it this like walled box.
So, that's not, in the end, like it's kind of going to perform
the same way, I just had this like extra loops in there.
But this is like an important clarification.
So now, continue with the challenge wherever I left off.
But this will be the code that's
actually there when you look at this later.
If you can make it to the end of this video.
See you soon.
(otherworldly music)
We have set bounds now.
Alright.
Oh, alright, fine.
I don't think I need to pass in N.
N is just a constant everywhere.
So let's simplify that, so let's go through this.
Fluid is fine.
I had velocity diffusion.
Linear solve, this has to come back.
Without the N.
Project.
Linear solve, boundaries, no N.
No N.
(bell dings)
So, a couple things.
People in the chat are mentioning,
I have no idea what's going on.
I kind of don't either, maybe this isn't a video
that's worth watching, I'm not sure, I mean at this point,
this is just showing that I'm invested in trying to do.
And I will come back to it.
But this I think is an example of why I'm refactoring
and sort of thinking about code and variable meaning.
It's a good thing, and there's nothing wrong with what was
done before, this is probably perfectly appropriate
and quite common styles from like over 10 years ago.
But I think we could hopefully do better
and I could come back and refactor this.
But I realize that I totally go this wrong
and I have some code that I did the other day here.
There we go. (bell dings)
This is what I've been trying to do.
I want to add those two things together times S zero,
and those two things together times S one, okay.
Alright, no errors.
Woo, woo, no errors, no errors.
This is, by the way, what I really like to avoid doing.
I don't like to code for ever and ever and ever
without testing, so if this like actually works
when I start to render it, I will be
completely and totally god smacked. (laughs)
Okay.
So now, I need a time step function.
So I need a function, time step,
to step through every moment of time with the fluid.
Going back to Mike Ash's paper,
there is a step function, which I'm going to grab.
Course I just wrote it out here.
But what I want is void step.
And I have all these as variables, but let's,
let's change that to this dot.
Replace float space star.
I know you can't see this.
Replace float space star
with float bracket star, no, float bracket space.
I don't need, and I don't need the N,
or the number of iterations.
I don't need the Z.
So, iterations is four.
Okay, so let's look at that.
The process here is what?
Diffuse the velocities.
Right, diffuse the velocities.
Based on the time step and viscosity.
The Xs and the Ys, the one and the two
is controlling the set bounds function.
Then project, which is clean everything up
to make sure it's the same amount of fluid everywhere.
The velocities.
Then run advection, advect on the velocities, X and Y.
Then, clean all that up.
Then, diffuse the density and advect the density.
The density doesn't need the project step,
because the density doesn't remain consistent around it,
it's the dye that actually
is moving around and is inconsistent.
So we should really be in, we really got
this whole thing programmed now.
We should've just started with,
definitely should've started with memos library. (laughs)
So if I say background zero,
and I say fluid step,
let's run this.
Hey, no errors.
Oh, happy day, I mean, give me a,
I'm going to be shocked if I can actually render this.
But let's actually render this.
So first, let's do a thing where,
where as I drag the mouse, I'm going to add dye.
So I'm going to add dye in mouseX and mouseY,
some amount of dye.
Okay?
Oh, and I need to say fluid.addDensity.
Alright, no errors, that's good.
That function's called.
Now let's try to render.
Let's call renderD for rendering the density.
So how am I going to do that?
You know what I should actually really do?
I'm going to change this to 64.
And then I'm going to create a variable
called scale, and make that four.
And then,
and then I'm going to say size is N times scale.
N times scale, so I'm actually going to like,
lower it, and let's actually, let's make the scale 10.
So now we should have, and if I comment this off,
I should have 640 by 640 window, okay.
So now fluid renderD.
I'm going to go here, I'm going to add
a function in the fluid class.
I kind of want to take all these functions out
and put them in a separate tab.
I'm going to call this more fluid.
So just like the fluid class this year,
which has less stuff in it,
then I'm going to say, renderD.
So this should be pretty easy
in that I'm just going to say int I equals zero.
I is less than N.
I plus plus.
So I'm going to loop through.
J, then I'm going to use J.
Then X is going to be I times scale.
Y is J times scale.
I'm going to draw a rectangle at XY
scale, which is really a square, scale scale.
By the way.
Did you know that Processing
now has the function square in it?
I'm going to draw a square, get rid of that.
Oh, that's going to work, I've never used this before.
This is a momentous occasion four hours into this video.
Then I'm going to say, fill,
255, 100.
I just want to see, like just going to put arbitrary fill.
Great, so I can see all of those squares are there.
I'm going to say no stroke.
And then now what I want to do is the density
is this.density at the index I and J,
and I want to, I'm going to make that the alpha.
So D is D.
So the alpha of that square,
I mean it could just be the brightness,
but let's just make it the brightness, directly.
So it's black, oh, right out of bounds exception.
Okay, interesting.
Ah, I forgot.
So if I am adding fluid at mouseX and mouseY,
I have to divide that now by scale.
'Cause I scaled up the screen.
So that has to be divided by scale.
There we go, okay.
So you can see I'm adding the fluid.
Great.
Alright, so that's something.
Now there's, I guess the water is completely still.
Didn't I say it would diffuse anyway? (laughs)
Let's try adding some, let's try also adding some velocity.
So I need to just also add, you know what I'm going to do?
I have an idea.
When I do that, I'm also going to add velocity.
At that spot, and add velocity expects what?
Expects the X and Y and the amount X and Y, okay.
So the X and Y, and now, I'm going to do this.
Let amountX equal mouseX minus pmouseX.
PmouseX is the previous mouse possession.
So the velocity is going to kind of be
in the direction that I'm dragging.
Not let, float.
Java, not Javascript.
And an amount Y is mouseY.
PmouseX.
So let's add that amount.
Oh, Y, let's see what happens here.
Ooh, okay, I don't mind that so much.
We got an array out of bounds exception.
So what went out of bounds?
One thing I can do here, (laughs)
where's my index function?
X equals constrain X two zero N minus one.
Let's just put a constrain in here.
I don't know whether I have a bug
in my code or this is actually necessary.
But let's do that.
Oh, there we go.
Oh my god.
Ah.
(bell dings)
(whistle whistles)
(dramatic music)
I think what I want to do is add something that always like,
yeah, 'cause look at this, I'm going to always like fade out,
I'm going to add something to like fade.
I'm going to add a function called fadeD.
I should really render the velocity, where I'm going to say,
four and I equals zero,
I is less than this.density.
Length.
I plus plus.
Density index I minus equal one.
And then density index I equals constrain.
Density index I, actually no, I can just do this.
Equals constrain, here.
FloatD equals density index I.
And then constrain,
I want to set it equal to D minus one,
but stay within zero and 255.
And then let me add fluid render fadeD.
I also want to render the velocity field.
So just because if I add a lot,
I mean that's too much fade.
I just want it to like lose,
I just want it to like, the dye to like peter out.
Right now I think the dye is kind of like infinite.
It's additive.
So but this obviously maybe I just
want it to peter out more slowly.
Alright, this is a nice little swirl effect.
Now, here's what I want to do.
I would like to, I mean the other thing I could do
is make this 128 and change the scale to five.
Can it, yeah.
So that's nicer, it can still perform,
the frame rate is still fine.
So what I want to do is I actually want to now,
I'm going to just use Perlin noise really quickly.
I'm just going to add density to the middle.
This is what I should be doing.
So what happens if I add,
just continuously add density all the time
in Draw, and this has to be converted into an integer.
This is kind of silly what I'm doing, but it's fine.
Alright, so now, the dye is just sitting there
in the middle that I can kind of drag around.
But what I want to do now is,
I'm going to take this, and I'm going to use Perlin noise.
I'm going to say, noise of some time,
map noise of some time, which goes between zero and one
to between negative one and one.
And amount Y
is the same, but I'm going to look at like a different part
of the noise base, if you're not familiar with Perlin noise,
I would encourage you to check out my Perlin noise videos.
And then, I'm going to say this also
is going to be the same at this spot here.
So I'm going to say int center X equals this.
Int centerY equals that.
Then add some density at center X, center Y.
You could actually make this like
a random amount of density maybe.
And then add the velocity according to Perlin noise.
I am going to create a variable called T,
which is kind of like the offset through the noise space.
And I'm going to have T, oh, I should be using,
ah, I changed my mind, I'm going to use a Pvector.
And I'm going to say, angle equals noise of T
times two pi.
And then I'm going to say a PVector
is PVector from angle.
Create a vector
from angle,
that angle.
And then, then the time will go up.
And this is V.X, V.Y.
And let's try, I don't know what it should be,
but let's try like making the vector
a little bit stronger right now.
So let's see here.
There we go, yeah, this is what I was trying to do.
Maybe it is
too strong.
Yeah, this is kind of the idea.
I'm not seeing, I just wanted it to kind of like
shoot it out randomly.
Let's take out the fade, let me try,
someone's telling me I swapped VX and VY at line 96 and 97.
Oh, okay, thank you.
So that's an important fix. (laughs)
I mean visually we were getting some results anyway.
Let me go back to making it,
let me do this.
Let me add, alright, one thing I could try,
and I think you get the idea,
but one thing that I could try is I could try
kind of adding some density in a little grid pattern.
This might be sort of nice to do,
and then I could say, plus I.
So this is just adding a lot more density.
And let me turn that back to zero.
And,
yeah, that's kind of, this is more like what I wanted.
There we go, okay.
That was too, the velocity was too strong.
Okay.
So here we go, this is more what I was looking to try to do.
I've got sort of Perlin noise controlling this angle.
I'm just like dropping in dye constantly.
It'd be sort of interesting to like sprinkle dye everywhere.
There's so much you can do with color,
to drop particles in, so let's recap what I've done.
I wanted to have some type
of visualization of turbulent fluid.
I think I mostly have it.
And I have ported all the code from Mike Ash's page.
I'm actually going to release this code as is.
This was like super long.
I would really like to come back.
And I would like to come back and try, when this,
try refactoring this and try
to understand the pieces of it a bit more.
I would like to, and then also I would like to create
a particle simulation on top of this where I can drop in
particles that sort of fly around based on the vector field.
That is more like the reference on the GitHub thread.
So I'm kind of imagining something like this,
where I'm actually just going to draw particles moving around
according to the vector field and see what that looks like.
Then ultimately I think what I want to do
is try maybe moving to use one of these libraries,
which has more, a sophisticated solver,
I think this project by Gabriel Weymouth
I would like to take a look at.
So I'll probably do some followup videos.
I don't know if this was helpful or interesting.
It's useful for me to try to do it.
Took over an hour.
But this is the process, I think, but if anything,
this was useful to sort of see the process
of trying to figure out somebody else's code written
in another language and port it.
We could try porting this in JavaScript
and see if it performs even with just like
raw 2-D canvas, I'd be curious to see that.
So go ahead and do that.
Thank you to some comments from many people in the chat,
but specifically Kaye Wiegmann, who pointed out
that I do not have turbulence if the viscosity is zero.
So I now adjusted some things.
So I created a viscosity of, let me just zoom into this.
I created a viscosity, a pretty low viscosity,
but it's a nonzero viscosity.
I also sort of changed the way I'm adding the density,
and removed that fade, and you can see, this has,
you'll notice, this has more of the quality of turbulence.
And yes, I challenge you now, take this code and add,
it's filling up now, add rainbow colors.
So, what do you keep a separate RGB channel?
So the things that I did, just you see,
is like, that I'm going to publish to go with this,
is I actually have a pretty low viscosity.
Ooh, why did the time step get so big?
The other thing that I also did is I adjusted
the range of the angle to wrap around twice.
So two pi times two, because Perlin noise
is going to tend to give me values around .5,
so it was really sticking around pi.
So this now is my final version,
let me fix up some stuff here.
If you're at the end of this video, hashtag TeamTurbulence.
Okay?
Hashtag TeamTurbulence on Twitter.
Let's blow it up, of course, no one's going to
have made it to the end of this video.
But let's do this, color mode, HSB one.
So, oh, my range is between zero and 255.
So let's do a fill.
D 255, 255, 100, I think this is not what I want to do.
Yeah, there we go, except,
there we go. (laughs)
Hue, saturation, brightness, that's pretty interesting.
Team Turbulence all the way.
This is my turbulence song.
(driving electronic music)
Team Turbulence!
Ch ch ch ch ch.
(driving electronic music)