Yurim Lee
Welcome to my archive :D
Interests: Physics Based Modeling | Animation | Medical Imaging
Modeling Water Surface
images of the same water surface at different angles
Rendering Water Surface & Results
Once we obtain some set of coordinates of our fluid particles via this method, rendering comes into play.
Jeffrey Katzenberg, the producer of the movie Shrek, has said that the hardest shot in Shrek was the pouring of milk into a glass. [Hiltzik and Pham 2001] As such, simulating moving liquid in animations is a big task in the animation industry. Liquid surrounds us yet it is hard for the artists to draw moving liquid.
In 2018 summer, I got a chance to work on rendering physics based modeling of water waves created by Professor Gibou. Octrees were used in order to simulate small scale waves, and parallel implementation for free surface to capture speed.
Here I hope to share what I enjoyed learning throughout my summer.
water surface modeled by Professor Gibou
water surface rendered by me
Water has few optical properties we need to capture in order to make our eyes perceive these coordinates as a water surface. Rendering enables us to capture these properties.
Our rendered results using photon mapping and radiosity can be viewed above. Future work involves creating caustics patterns.
Comment
Navier Stokes Equations
First, in order for us to recreate the water waves, we need to understand how they behave and the physics behind fluid motion. Say we have some velocity field
for our fluid particle, which defines the velocity the particle passes through point at time . Let as write the velocity of the particle at time as and the velocity of the particle at time as
While I have explored few options for creating caustics patterns: i.e. (1) dividing the surface into triangles, mapping the triangles onto the bottom of the "water cup" using Snell's law, and setting each triangle's "brightness" proportional to the ratio of the area of the surface triangle to the area of the bottom triangle (2) recreating a very bright light source at a large distance - none of the methods turned out to be satisfying.
I find this problem very intriguing. This whole project was super fun. I would like to thank Professor Gibou for giving me this project.
Then
which gives us the acceleration of our particle
where each are velocities along x-axis, y-axis, and z-axis. Taking out the from above, we obtain
This is called the material derivative. We can rewrite this and see where the name comes from.
The second last term is local derivative which refers to the change in material with respect to time at certain location. The last term is convective derivative which refers to change in material with respect to time due to the movement of the particle.
So for instance, if we have
where is the temperature,
refers to the change in temperature from the flow field itself and refers to the change in temperature that arises from the particle moving from one point within the flow field to another.
Coming back to where
is the velocity field, gives us the change within the velocity field and gives us the change in velocity that results from the particle moving from one point within the velocity field to another.
We have now successfully obtained , the acceleration of our fluid particle.
Let us plug this into Newton's classic second law:
for our fluid particle. This is the right hand side of the Navier Stokes equation for our fluid.
For the left hand side of the Navier Stokes equation, we need to start by assuming water to be Newtonian fluid, or fluid in which the viscous stresses arising from its flow are linearly correlated to the local strain rate at every point.
"moving plate"
"fixed plate"
Simply, imagine a "moving plate" and a "fixed plate" and some fluid particles between them. A Newtonian fluid is a fluid that has a constant viscosity regardless of the velocity of the moving plate.
This leads to following equations of normal stress and shear stress:
The derivation will not be shown here as I only intend to go over the general idea but if you're interested, you can see the derivation on Hermann Schlichting's Boundary Layer Theory.
These equations allow us to write
Once again, plugging in the equations of normal stress and shear stress to the equation above will not be shown here for the sake of time and space but this gives us the left hand side of our Navier Stokes equation.
Equating the left hand side of the equation and the right hand side of the equation, we finally get
- our Navier Stokes equation for a Newtonian fluid. We can see how the first term comes from gravity, the second term, pressure, third term, normal stress, fourth term and fifth term, shear stress. The right hand side of the equation is simply Newton's second law.
Assuming our fluid is incompressible, we can simplify this equation further to
Then for inviscid fluid, we can simplify this even further as its viscosity is equal to zero.
Putting this in a vector format, we get
which is our final Navier Stokes equation for inviscid, incompressible fluid.
Navier Stokes Equations On Octrees
Now how do we solve our Navier Stokes equation?
We first constrain our Navier Stokes equation with:
conservation of mass (constant ) and conservation of momentum ( ). Absorbing the constant
into and generalizing as some external forces ,
we can rewrite our Navier Stokes equation as
Following [Losasso et. al 2004]s steps, we compute "an intermediate velocity " ignoring the pressure term. Then we compute the velocity update where the pressure is defined as the solution to the Poisson equation .
[Losasso et. al 2004] proposed a technique for discretizing the Poisson equation on an unrestricted octree data structure. The resulting linear system enables the use of fast solution methods. While we will not go into this technique in detail, the basic idea behind this technique follows:
An unrestricted octree data structure has five cells, one large cell and neighboring four smaller cells.
the octree data structure
The velocity components are defined on the cell faces, while the pressure is defined at the center of the cell. The density, temperature and level set function φ are stored at the nodes.
Note the paper uses u instead of v to denote velocity.
credit: Stanford PhysBAM
Discretization of the pressure gradient.
credit: Stanford PhysBAM
We solve for the pressure gradient by defining
at the cell face - See the figure above. Note this is the standard differencing with as the weighted average value. - We write
's contribution to row 1 and row 2 of the matrix representing the linear system of equations
in terms of the area of the cell face and and and , which yields a symmetric linear system that can be efficiently inverted with a PCG method [Saad. 2003].
Once we compute the pressure gradient at every cell face, we can compute the velocity update.
We can "scale" the Poisson equation by multiplying the volume of the large cell to the both handsides of the equation. We can then apply Green's theorem and compute the divergence of our intermediate velocity in terms of the areas of cell faces and the intermediate velocities defined on cell faces.
