(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)