<p>Suppose I have a fluid which, like water, is incompressible; that is, its density is constant. I can find out how this fluid is flowing; let's say I inject dye at various points and follow how the dye particles move. I can estimate a velocity at each point, a velocity field <span class="math-container">$\mathbf{w}$</span>. But water is incompressible; this means that the true velocity field <span class="math-container">$\mathbf{u}$</span> has to be mass-conserving. I want to "correct" my estimate <span class="math-container">$\mathbf{w}$</span> to find the closest mass-conserving velocity field.</p>

<p>The local accumulation or depletion of the fluid is computed with the <em>divergence</em> of the vector field, <span class="math-container">$\operatorname{div}\mathbf{w}$</span>. (Stam writes the divergence of <span class="math-container">$\mathbf{w}$</span> as <span class="math-container">$\nabla\cdot\mathbf{w}$</span>.) The divergence is a <em>derivative</em>, which is important here because that means it's <em>linear</em>; the divergence of the sum of two vector fields is the sum of their divergences:
<span class="math-container">$$ \operatorname{div}(\mathbf{u} + \mathbf{v}) = \operatorname{div}\mathbf{u} + \operatorname{div}\mathbf{v}. $$</span>
We want to find a mass-conserving field; that is, we want to find a field with <em>zero divergence</em>. That means that we want to construct a field <span class="math-container">$\mathbf{v}$</span> with the same divergence as <span class="math-container">$\mathbf{w}$</span>; then
<span class="math-container">$$ \operatorname{div}\mathbf{w} - \operatorname{div}\mathbf{v} = \operatorname{div}(\mathbf{w} - \mathbf{v}) = 0 $$</span>
so we can get our mass-conserving field <span class="math-container">$\mathbf{u} = \mathbf{w} - \mathbf{v}$</span>.</p>

<p>The only question is how we construct a field with a given divergence. For reasons I'll explain a few paragraphs down, we say that <span class="math-container">$\mathbf{v}$</span> is the gradient of some scalar field <span class="math-container">$q$</span>, <span class="math-container">$\mathbf{v} = \operatorname{grad}q$</span>. (Stam writes this as <span class="math-container">$\nabla q$</span>.) Then we have
<span class="math-container">$$ \operatorname{div}\mathbf{w} = \operatorname{div}\mathbf{v} =
\operatorname{div}\operatorname{grad}q.$$</span>
It's worth looking at what this says. The divergence of a vector field is a scalar field; at each point, it gives the local accumulation or depletion of the fluid. Both the divergence and the gradient are derivatives (which is why they can also both be written with <span class="math-container">$\nabla$</span>). The equation is thus saying that the scalar field (the divergence of <span class="math-container">$\mathbf{w}$</span>) is a particular second derivative of the scalar field <span class="math-container">$q$</span>. This derivative is called the <em>Laplacian</em> (Stam writes it as <span class="math-container">$\nabla^2 q$</span>). The differential equation relating scalar fields <span class="math-container">$q$</span> and <span class="math-container">$a$</span>
<span class="math-container">$$ \operatorname{div}\mathrm{grad}\ q = a $$</span>
is the <em>Poisson equation</em>. It's not so much that you use it to <em>estimate</em> <span class="math-container">$q$</span> as that it's the actual <em>definition</em> of <span class="math-container">$q$</span>.</p>

<p>So, to "correct" our estimate <span class="math-container">$\mathbf{w}$</span>, we can compute its divergence (a scalar field)
<span class="math-container">$$ f = \operatorname{div}\mathbf{w}, $$</span>
solve the Poisson equation to find the scalar field <span class="math-container">$q$</span>
<span class="math-container">$$ \operatorname{div}\operatorname{grad} q = f, $$</span>
then subtract <span class="math-container">$\operatorname{grad}q$</span> from <span class="math-container">$\mathbf{w}$</span> to find our mass-conserving velocity field
<span class="math-container">$$ \mathbf{u} = \mathbf{w} - \operatorname{grad}q. $$</span></p>

<hr>

<p>Above we just assumed that <span class="math-container">$\mathbf{v} = \operatorname{grad}q$</span> for some <span class="math-container">$q$</span>. What we actually want is to find a field <span class="math-container">$\mathbf{v}$</span> that is in some sense the "smallest" field with the given divergence, so we can find the "closest" mass-conserving field <span class="math-container">$\mathbf{u} = \mathbf{w} - \mathbf{v}$</span>; here I'll explain how the gradient of a scalar field satisfies this requirement.</p>

<p>We have a fixed amount of fluid moving through our space. 
If our flow is mass-conserving, then we have a conserved amount of fluid, without local accumulation or depletion. If we have any motion at all, then, the fluid must move in <em>cycles</em>. In the same way as we measure the local accumulation with the divergence operator, we can measure the local cyclic behaviour of the flow by using another vector operator, the <em>curl</em> (which can be written <span class="math-container">$\nabla\times\mathbf{v}$</span>). The curl at a point gives, broadly speaking, the axis about which the flow is circling near that point.</p>

<p>Now, consider the gradient of a scalar field, <span class="math-container">$\operatorname{grad} q$</span>. (For instance, <span class="math-container">$q(x,y) = \sin x \cos y$</span>.) The potential can be drawn as a height field:
<a href="https://i.sstatic.net/sVzgE.png" rel="noreferrer"><img src="https://i.sstatic.net/sVzgE.png" alt="plot of z = sin(x)*cos(y)"></a>
Then the gradient of this field puts at each point a vector pointing in the steepest uphill direction:
<a href="https://i.sstatic.net/Okc02.png" rel="noreferrer"><img src="https://i.sstatic.net/Okc02.png" alt="vector field grad(sin(x) * cos(y))"></a>
But since you're always going "uphill" as quickly as you can, I think it's clear that there's <em>no circulation</em> in this vector field! In fact, for <em>any</em> scalar potential <span class="math-container">$q$</span>,
<span class="math-container">$$ \operatorname{curl}\operatorname{grad}q = 0 ! $$</span></p>

<p>There's also a second similar identity: <span class="math-container">$\operatorname{div}\operatorname{curl}\mathbf{v} = 0$</span>; unfortunately, I can't think of a good intuition for that one. (But look at <a href="https://math.stackexchange.com/a/27010/152310">this Math.SX answer</a> for a more mathematical explanation.) We can use these two identities to explicitly write the Hodge decomposition: <span class="math-container">$\mathbf{w} = \mathbf{v} + \mathbf{u}$</span>, where <span class="math-container">$\mathbf{v}$</span> has no curl, and <span class="math-container">$\mathbf{u}$</span> has no divergence, or
<span class="math-container">$$ \mathbf{w} = \operatorname{grad}q + \operatorname{curl}\mathbf{\nu},¹ $$</span>
for some scalar field <span class="math-container">$q$</span> and vector field <span class="math-container">$\mathbf{\nu}$</span>.</p>

<hr>

<p>¹The Hodge decomposition is actually more complicated if your space isn't simple, i.e. finite and ball-like; on more complicated spaces, like say the surface of a torus, we have curl-free fields which aren't the gradient of anything, and divergence-free fields which aren't the curl of anything. For most fluid modelling in CG, though, this isn't usually an issue.</p>

CC BY-SA 4.0
2018-10-21T13:25:45.037Z
Thanks for the answer! your explanation of origin of the formula used in Stam's paper, really helps to get a deeper understanding of the math.
Mn9
CC BY-SA 4.0
2018-10-19T12:22:48.717Z
2018-10-19T12:22:48.717Z
gilgamec
<p>I have just started implementing a fluid simulator based on the navier stokes equation.
I am following Jos Stam's paper to do so.
In this paper, it says that the velocity field must be mass conserving. Which I understand as, the velocity that flows into a cell should be the same the velocity that flows out.
To get this mass conserving field, the hodge decomposition is used, and the text says:</p>

<pre><code>"every vector field is the sum of a mass conserving field and a gradient field." 
Velocity field = mass conserving field + gradient field
</code></pre>

<p>1) Is there a intuitive way to understand why this statement above is true? Maybe I have to read more on Hodge decomposition to understand this.</p>

<p>2) The text also says that the gradient field can be computed by solving the poisson equation, Is there a simple explantion of why the poisson equation can be used to approximate this field? I'm wondering what else could be used instead. </p>

<p>Could you point me to some resources that help understand my two questions intuitively, before delving into all the math?</p>

<fluid-sim>
Hodge Decomposition in fluid simulations
Mn9
# Hodge Decomposition in fluid simulations

(Asked by: username_0 on 2018-10-19)

I have just started implementing a fluid simulator based on the navier stokes equation.
I am following Jos Stam's paper to do so.
In this paper, it says that the velocity field must be mass conserving. Which I understand as, the velocity that flows into a cell should be the same the velocity that flows out.
To get this mass conserving field, the hodge decomposition is used, and the text says:

```
"every vector field is the sum of a mass conserving field and a gradient field." 
Velocity field = mass conserving field + gradient field

```

1) Is there a intuitive way to understand why this statement above is true? Maybe I have to read more on Hodge decomposition to understand this.

2) The text also says that the gradient field can be computed by solving the poisson equation, Is there a simple explantion of why the poisson equation can be used to approximate this field? I'm wondering what else could be used instead.

Could you point me to some resources that help understand my two questions intuitively, before delving into all the math?

**Tags:** fluid-sim


---

## 1 answers

### Answer by username_1 on 2018-10-19. Score: 7 | Accepted answer

Suppose I have a fluid which, like water, is incompressible; that is, its density is constant. I can find out how this fluid is flowing; let's say I inject dye at various points and follow how the dye particles move. I can estimate a velocity at each point, a velocity field $\mathbf{w}$. But water is incompressible; this means that the true velocity field $\mathbf{u}$ has to be mass-conserving. I want to "correct" my estimate $\mathbf{w}$ to find the closest mass-conserving velocity field.

The local accumulation or depletion of the fluid is computed with the *divergence* of the vector field, $\operatorname{div}\mathbf{w}$. (Stam writes the divergence of $\mathbf{w}$ as $\nabla\cdot\mathbf{w}$.) The divergence is a *derivative*, which is important here because that means it's *linear*; the divergence of the sum of two vector fields is the sum of their divergences:
$$ \operatorname{div}(\mathbf{u} + \mathbf{v}) = \operatorname{div}\mathbf{u} + \operatorname{div}\mathbf{v}. $$
We want to find a mass-conserving field; that is, we want to find a field with *zero divergence*. That means that we want to construct a field $\mathbf{v}$ with the same divergence as $\mathbf{w}$; then
$$ \operatorname{div}\mathbf{w} - \operatorname{div}\mathbf{v} = \operatorname{div}(\mathbf{w} - \mathbf{v}) = 0 $$
so we can get our mass-conserving field $\mathbf{u} = \mathbf{w} - \mathbf{v}$.

The only question is how we construct a field with a given divergence. For reasons I'll explain a few paragraphs down, we say that $\mathbf{v}$ is the gradient of some scalar field $q$, $\mathbf{v} = \operatorname{grad}q$. (Stam writes this as $\nabla q$.) Then we have
$$ \operatorname{div}\mathbf{w} = \operatorname{div}\mathbf{v} =
\operatorname{div}\operatorname{grad}q.$$
It's worth looking at what this says. The divergence of a vector field is a scalar field; at each point, it gives the local accumulation or depletion of the fluid. Both the divergence and the gradient are derivatives (which is why they can also both be written with $\nabla$). The equation is thus saying that the scalar field (the divergence of $\mathbf{w}$) is a particular second derivative of the scalar field $q$. This derivative is called the *Laplacian* (Stam writes it as $\nabla^2 q$). The differential equation relating scalar fields $q$ and $a$
$$ \operatorname{div}\mathrm{grad}\ q = a $$
is the *Poisson equation*. It's not so much that you use it to *estimate* $q$ as that it's the actual *definition* of $q$.

So, to "correct" our estimate $\mathbf{w}$, we can compute its divergence (a scalar field)
$$ f = \operatorname{div}\mathbf{w}, $$
solve the Poisson equation to find the scalar field $q$
$$ \operatorname{div}\operatorname{grad} q = f, $$
then subtract $\operatorname{grad}q$ from $\mathbf{w}$ to find our mass-conserving velocity field
$$ \mathbf{u} = \mathbf{w} - \operatorname{grad}q. $$

---

Above we just assumed that $\mathbf{v} = \operatorname{grad}q$ for some $q$. What we actually want is to find a field $\mathbf{v}$ that is in some sense the "smallest" field with the given divergence, so we can find the "closest" mass-conserving field $\mathbf{u} = \mathbf{w} - \mathbf{v}$; here I'll explain how the gradient of a scalar field satisfies this requirement.

We have a fixed amount of fluid moving through our space.
If our flow is mass-conserving, then we have a conserved amount of fluid, without local accumulation or depletion. If we have any motion at all, then, the fluid must move in *cycles*. In the same way as we measure the local accumulation with the divergence operator, we can measure the local cyclic behaviour of the flow by using another vector operator, the *curl* (which can be written $\nabla\times\mathbf{v}$). The curl at a point gives, broadly speaking, the axis about which the flow is circling near that point.

Now, consider the gradient of a scalar field, $\operatorname{grad} q$. (For instance, $q(x,y) = \sin x \cos y$.) The potential can be drawn as a height field:
[![plot of z = sin(x)*cos(y)](https://i.sstatic.net/sVzgE.png)](https://i.sstatic.net/sVzgE.png)
Then the gradient of this field puts at each point a vector pointing in the steepest uphill direction:
[![vector field grad(sin(x) * cos(y))](https://i.sstatic.net/Okc02.png)](https://i.sstatic.net/Okc02.png)
But since you're always going "uphill" as quickly as you can, I think it's clear that there's *no circulation* in this vector field! In fact, for *any* scalar potential $q$,
$$ \operatorname{curl}\operatorname{grad}q = 0 ! $$

There's also a second similar identity: $\operatorname{div}\operatorname{curl}\mathbf{v} = 0$; unfortunately, I can't think of a good intuition for that one. (But look at [this Math.SX answer](https://math.stackexchange.com/a/27010/152310) for a more mathematical explanation.) We can use these two identities to explicitly write the Hodge decomposition: $\mathbf{w} = \mathbf{v} + \mathbf{u}$, where $\mathbf{v}$ has no curl, and $\mathbf{u}$ has no divergence, or
$$ \mathbf{w} = \operatorname{grad}q + \operatorname{curl}\mathbf{\nu},¹ $$
for some scalar field $q$ and vector field $\mathbf{\nu}$.

---

¹The Hodge decomposition is actually more complicated if your space isn't simple, i.e. finite and ball-like; on more complicated spaces, like say the surface of a torus, we have curl-free fields which aren't the gradient of anything, and divergence-free fields which aren't the curl of anything. For most fluid modelling in CG, though, this isn't usually an issue.

> **Comments:**
> **username_0**: Thanks for the answer! your explanation of origin of the formula used in Stam's paper, really helps to get a deeper understanding of the math.  

---

