I'm trying to improve the advection method in a 2D-windfield. The Navier-Stokes Equations (NSE) are currently used for the influence of pressure, viscosity,... I am just focusing on the convective acceleration term $(u\cdot\nabla)u$. So after rearranging the NSE I get $$\left(\frac{\partial \vec v}{\partial t}\right)_{\text{Adv}}=-(u\cdot\nabla)u.$$ For my case I can write this as: $$\begin{pmatrix} \frac{\partial u}{\partial t} \\ \frac{\partial v}{\partial t} \end{pmatrix}= \begin{pmatrix} -u(\frac{\partial u}{\partial x}+\frac{\partial v}{\partial y}) \\ -v(\frac{\partial u}{\partial x}+\frac{\partial v}{\partial y}) \end{pmatrix}.$$ My question is, can I now apply backward differencing $$\left(\frac{\partial \Phi}{\partial x}\right)_{i}=\frac{\Phi_i-\Phi_{i-1}}{\Delta x},\quad\left(\frac{\partial u}{\partial t}\right)_{i}=\frac{u_i^{t+1}-u_i^t}{\Delta t}~?$$ My final two equation then: $$u_{x,y}^{t+1}=u_{x,y}^{t}-\frac{u_{x,y}^{t}\Delta t}{\Delta x}(u_{x,y}^{t}-u_{x-1,y}^{t}+v_{x,y}^{t}-v_{x,y-1}^{t}) $$ $$v_{x,y}^{t+1}=v_{x,y}^{t}-\frac{v_{x,y}^{t}\Delta t}{\Delta x}(u_{x,y}^{t}-u_{x-1,y}^{t}+v_{x,y}^{t}-v_{x,y-1}^{t}). $$ Now I have an equation for the velocity $u$ in $x$-direction and $v$ in $y$-direction.

It's not working properly in the simulation. For example, if the direction of the wind should change over time, the direction stays the same but the velocity on the edges of the field (where the initial conditions of $u$ and $v$ are set) are getting bigger or smaller, depending on the "direction of the change". It seems like the information is not spreading properly over the whole wind field.

  • $\begingroup$ The expansion of your convective acceleration is wrong. Observe the brackets in your term $(\vec u\cdot \vec \nabla )\vec u$ with the vector $\vec u = (u,v)^T$ correctly. That will for sure give you wrong numerical results. $\endgroup$ Nov 22, 2023 at 15:59
  • $\begingroup$ @AtmosphericPrisonEscape What exactly do you say is wrong? First I take the scalar product of $ u $ and $ \nabla $. Then I multiply $\begin{pmatrix} -\frac{\partial u}{\partial x}+\frac{\partial v}{\partial y}\end{pmatrix}$ by the vector $\begin{pmatrix} u \\\ v \end{pmatrix}$ and get the equation above. $\endgroup$
    – CFD98
    Nov 22, 2023 at 17:19
  • $\begingroup$ The scalar $(\vec u \cdot \vec \nabla)$ is $u \partial/\partial_x + v \partial/\partial_y$, that is then multiplied into every element of $\vec u = (u,v)^T$. $\endgroup$ Nov 22, 2023 at 17:28

1 Answer 1


I advocate the view that the advection component is an operator, $$\left(\mathbf{u}\cdot\nabla\right)\mathbf{u}=\operatorname{adv}\mathbf{u};$$ hence you have, depending on what is most natural to you, $$\operatorname{adv}\mathbf{u}\equiv\operatorname{adv}\left(\begin{array}{c}u \\ v \end{array}\right)\equiv\left(\begin{array}{c}\operatorname{adv}u \\ \operatorname{adv}v \end{array}\right)\equiv \operatorname{adv}u_i\mathrm{e}_i.\tag{1}$$ Since $\operatorname{adv}=\mathbf{u}\cdot\nabla$, then Eq (1) is, $$\operatorname{adv}\mathbf{u}\equiv\left(u\partial_x+v\partial_y\right)\left(\begin{array}{c}u \\ v \end{array}\right)\equiv\left(\begin{array}{c}\left(u\partial_x+v\partial_y\right)u \\ \left(u\partial_x+v\partial_y\right)v \end{array}\right).\tag{2}$$ Then applying the backwards finite difference to the derivatives, $$\operatorname{adv}\mathbf{u}=\left(\begin{array}{c} u_{i,j}\left(u_{i,j}-u_{i-1,j}\right)/\Delta x +v_{i,j}\left(u_{i,j}-u_{i,j-1}\right)/\Delta y\\ u_{i,j}\left(v_{i,j}-v_{i-1,j}\right)/\Delta x +v_{i,j}\left(v_{i,j}-v_{i,j-1}\right)/\Delta y \end{array}\right)$$ So along with your time derivative component, you'd see something like, $$u_{i,j}^{n+1}=u_{i,j}^n-\frac{\Delta t}{\Delta x}u_{i,j}^n\left(u_{i,j}^n-u_{i-1,j}^n\right)-\frac{\Delta t}{\Delta y}v_{i,j}^n\left(u_{i,j}^n-u_{i,j-1}^n\right)\tag{3}$$ and similarly for the $v$ component of the velocity.

If you wanted to increase the accuracy by using a 2nd order scheme (still with backwards differencing), you'd have something like, $$u_{i,j}^{n+1}=u_{i,j}^n-\frac{\Delta t}{\Delta x}u_{i,j}^n\left(\frac{3}{2}u_{i,j}^n-2u_{i-1,j}^n+\frac{1}{2}u_{i-2,j}^n\right)-\frac{\Delta t}{\Delta y}v_{i,j}^n\left(\frac{3}{2}u_{i,j}^n-2u_{i,j-1}^n+\frac{1}{2}u_{i,j-2}^n\right)$$ It might also be useful to look over some of the predictor-corrector methods, such as the MacCormack method, wherein you compute an intermediate state, $u^p_{i,j}$, using Equation (3) and then use that state to evaluate the next step using, $$u^{n+1}_{i,j}=\frac{1}{2}\left(u^n_{i,j}+u^p_{i,j}\right)-\frac{\Delta t}{\Delta x}u_{i,j}^p\left(u_{i,j}^p-u_{i-1,j}^p\right)-\frac{\Delta t}{\Delta y}v_{i,j}^p\left(u_{i,j}^p-u_{i,j-1}^p\right)$$

  • $\begingroup$ Is it possible to get a more accurate solution by using more points on the grid to evaluate the velocity at the next time step? Which methods of finite differencing are preferable? $\endgroup$
    – CFD98
    Nov 23, 2023 at 13:01
  • $\begingroup$ Yes, there are higher-order stencils available (cf. this Wikipedia entry or ENO/WENO methods). The downside is that you need more "ghost cells" to apply boundaries. As to which is preferable, that's really a question for experimentation: if low-res schemes don't work, try different ones until you find one that does. $\endgroup$
    – Kyle Kanos
    Nov 23, 2023 at 13:15
  • $\begingroup$ But if I take for example the coefficients for accuracy 2 (Backward finite difference) and put 3 points in the brackets of your last equation with these coefficients, it doesn't give a proper solution anymore. And that does make sense to me because there is more "weight" on the "-1" point. So if the velocity just changes a bit from point "0" to "-1" its a bigger change than with just taking the two points without any "weighting". Am I missing something out? $\endgroup$
    – CFD98
    Nov 26, 2023 at 19:22
  • $\begingroup$ I've added some additional color on the higher-order schemes. You should still get a proper solution using the higher order scheme, unless there are bugs in the code. $\endgroup$
    – Kyle Kanos
    Nov 27, 2023 at 11:49
  • $\begingroup$ Also, are you sure that your time step satisfies the Courant condition? $\endgroup$
    – Kyle Kanos
    Nov 27, 2023 at 11:54

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.