I am running a 2D fluids simulation with a stochastic forcing $f$ in a doubly-periodic box, i.e. solving

$$ \frac{\partial \nabla^2 \psi}{\partial t} = J(\psi,\nabla^2 \psi) +f,$$ where $J$ is a Poisson bracket.

The forcing I've chosen is of the form

$$f= \sum_{k,l} c \sin(k x + \alpha_k) \sin(l y + \beta_l),$$ where the wavenumbers $k$ and $l$ are selected from a thin annulus centered at $k_f$ (a common prescription), and $\alpha_k$ and $\beta_l$ are phases which are in some sense stochastic (more on this later). (When $k=0$ or $l=0$ one needs a slightly different treatment but this is a minor detail).

One can show that $\int d^2 x \, \psi \frac{\partial \nabla^2 \psi}{\partial t} = - \frac{\partial E}{\partial t}$ where $E$ is the total kinetic energy. We can use this to compute the energy injection rate from the forcing:

$$\varepsilon = -\langle f \xi \rangle = \frac{N c^2}{4 k_f^2},$$ where $\nabla^2 \xi = f$ and $N$ is the number of wavevectors in the annulus. This motivates the choice $c = 2 k_f\sqrt{\frac{\varepsilon}{N}}$.

Now, we need to prescribe a temporal correlation for the forcing. The standard choice is white noise, i.e. $f$ is delta-correlated in time. In (for example) Appendix A of Srinivasan and Young (2012), the authors select the phases i.i.d. from a uniform distribution and aver that the forcing must be normalized by $1/\sqrt{\delta t}$ ($\delta t$ being the time-step of the integration algorithm) in order to ensure that it is delta-correlated. This raises two questions with which I'm struggling:

  1. How, precisely, does this lead to delta-correlated forcing? I'm having some trouble showing it analytically.
  2. What now of the energy injection rate? Isn't it altered by a factor of 1/$\delta t$? And aren't the dimensions now compromised?

Furthermore, as is pointed out in the same appendix, in a Runge-Kutta algorithm, the forcing must be kept reasonably smooth during the course of a time-step, so in that paper, they select the phases within a time-step by linear interpolation. I'm finding this tricky to implement with the library I'm using, so I had the idea to instead update the phases by their own random walk: $$\alpha_k(t+\delta t) = \alpha_k(t) + \sqrt{\delta t} \eta$$ with $\eta\sim {\cal N}(0,\sigma^2),$ and the same for the $\beta_l$. Then one can show that this leads to the correlation function

$$ \langle f(\mathbf{x},0)f(\mathbf{x},t) \rangle = \frac{Nc^2}{4} \exp(-|t|/\tau)$$ where $\tau= 1/\sigma^2$. This forcing is nice in that it's nice and smooth and you can control the correlation time...unless you want white noise. Thus, a third question:

  1. Can this forcing prescription be adjusted simply so that, in the limit $\tau \to 0$, the forcing is temporally white? Does, say, normalizing $f$ by $1/\sqrt{\tau}$ work?

Thanks in advance to anyone who can help me with this.

  • 1
    $\begingroup$ To answer (what I understood is) the first point, you must make units right. Since $<\eta(t)\eta(t')>\propto\delta(t-t')$ units of $\eta$ are $1/\sqrt{t}$. When you inject it to the stochastic equation $dx=adt + bdt\eta$ the last term is actually proportional to $\sqrt{dt}$ $\endgroup$
    – Alexander
    Jun 4, 2020 at 1:02
  • $\begingroup$ I've thought about this more and it's easier for me to conceptualize in terms of the central limit theorem, which shows that, if $\partial_t \zeta = \xi(t)$ for some random $\xi$, $\zeta(\tau)-\zeta(0) \sim \sqrt{\tau \delta t} {\cal N} (0, \langle \xi^2 \rangle)$, so we need the normalization to ensure that $\langle \zeta^2 \rangle \sim t$ like a good Brownian motion. I'm still trying to figure out how this affects the rate of energy injection... $\endgroup$ Jun 4, 2020 at 1:38
  • $\begingroup$ I think I now roughly understand the energy injection issue. The rate of energy injection is not, as I wrote, $-\langle \xi f \rangle,$ but $-\langle \psi f \rangle$. It stands to reason, at least in a hand-waving sort of way, that $\psi \sim \sqrt{t} \xi$ (where again $\nabla^2 \xi = f$), so the rate of energy injection in fact is only constant in time once the normalization is introduced. $\endgroup$ Jun 4, 2020 at 1:59
  • $\begingroup$ Reading this paper aip.scitation.org/doi/pdf/10.1063/1.870050 helped me resolve these conceptual problems and identify my error $\endgroup$ Jun 4, 2020 at 2:55
  • $\begingroup$ consider posting it as an answer, with short description of the solution $\endgroup$
    – Alexander
    Jun 4, 2020 at 5:13

1 Answer 1


The need for the normalization can be illustrated as follows. Consider only the contribution of the forcing: $$ \partial_t \zeta = f(t).$$

If we discretize time in steps $\delta t$, $f$ is a sequence of iid random variables $\{f_i\}$ with mean 0. Then we have (setting $\zeta(0)=0$) $$\zeta(t) = \delta t \sum_i f_i.$$

According to the central limit theorem, for $t\gg \delta t$, we have $$\sum_i f_i \sim \sqrt{N}{\cal N}(0,\langle f^2 \rangle)$$ where $N=t/\delta t$. Thus

$$\zeta(t) \sim \sqrt{t \delta t} {\cal N}(0,\langle f^2 \rangle)$$ and it is evident that $f$ must be normalized by $1/\sqrt{\delta t}$ in order for its evolution to be independent of the timestep.

My issue with the energy injection rate largely stemmed from an error in calculation. The rate is $$\varepsilon = - \langle \psi f\rangle,$$ not $-\langle \xi f \rangle$ (where again, $\nabla^2 \xi = f$). Following K. Alvelius PoF 11, 1880 (1999), we have (considering only the forcing's contribution to the dynamics) $$\psi = \int_0^t d\tau \, \xi(\tau), $$ so $$\varepsilon = - \int_0^t d\tau \, \langle \xi(\tau)f(t) \rangle.$$

My previous, erroneous calculation did not have the correct units, and clearly missed a factor with dimensions of time. In the case of my implementation of the forcing, the energy injection rate follows as

$$\varepsilon =\frac{Nc^2}{4k_f^2} \tau_c (1- \exp(-|t|/\tau_c)$$.

After many correlation times, the second term is negligible and choosing $$ c= 2 k_f \sqrt{\frac{\varepsilon}{N \tau_c}}$$ fixes a chosen rate $\varepsilon$. Note that factor $1/\tau_c$, previously missed. This normalization also ensures that in the limit $\tau_c \to 0$, the volume of the correlation function is correctly conserved and a delta function is obtained, as desired.


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.