# Simulating atmospheric turbulence

This problem comes from simulation of atmospheric turbulence. It starts with the following transfer function in the frequency domain: \begin{equation} G_v(s)=\sigma k_1 \frac{1 + k_2 s}{\left(1+k_3 s\right)^2} \end{equation} This comes from the so-called Dryden turbulence model, which says that the speed of the turbulence perturbation on the spatial dimension $v$ is $G_v(s)$ applied to band-limited white noise. As usual in block-based system simulation, there are black-box assumptions which I cannot blindly accept and use, so I am trying to replicate the model with my own code.

Still in the frequency domain, I completed the equation with the meaning of the model, that is \begin{equation} y_v(s)=G_v(s)sW(s) \end{equation} where $y_v(s)$ is the turbulence speed and $sW(s)$ is the derivative of a Wiener process in the frequency domain, which should be equivalent to a white noise. Treating $W(t)$ as a generalized random process allows the existence of its $n$-th order derivatives and, by consequence, their equivalents in the frequency domain. Substituting the formula of $G_v(s)$ from above in the latter, we obtain \begin{equation} y_v(s)\left(1+k_3 s\right)^2=\sigma k_1 \left(1 + k_2 s\right)sW(s) \end{equation} Doing the inverse Laplace transformation to the time domain and abstracting the constants, we obtain a stochastic differential equation of the following form: \begin{equation} a_2 y_v^{\prime\prime}(t)+a_1 y_v^{\prime}(t)+a_0 y_v(t)+b_2W^{\prime\prime}(t)+b_1W^{\prime}(t)=0 \end{equation} where $W(t)$ is the Wiener process in the time domain and, by consequence, $W^{\prime}(t)$ should be white noise and $W^{\prime\prime}(t)$ should be the derivative of white noise.

Now, to simulate this equation numerically, it would be relatively easy if not for the $W^{\prime\prime}(t)$. With a We may interpret $W^{\prime}(t)$ either in the straightforward sense, \begin{equation} \frac{dW}{dt} \approx \frac{1}{h} \Big[ W(t+h) - W(t) \Big] \end{equation} and we can recur to Euler-Maruyama to solve this numerically. However, how to solve numerically for $W^{\prime\prime}(t)$?

Am I correct that you are asking how one could interpret a second-order derivative of a white noise?

You really need to consult a mathematician. I have an idea how you can define a second order derivative of white noise, but rely on it only if you are desperate, I don't know how correct this is.

First, recall how you view the first derivative through a finite difference: $$\frac{dW}{dt} = \frac{W(t+\Delta t) - W(t)}{\Delta t} = \frac{N_t(0,1) \sqrt{\Delta t}}{\Delta t} = N_t(0,1) \frac{1}{\sqrt{\Delta t}},$$ where $N_t(0,1)$ is a standardized Gaussian random variable with zero mean and variance of $1$. The subscript $t$ is just to remind that at each time $t$ you calculate a different random variable.

Now we can go by analogy and use a forward finite difference approximation for the second derivative: \begin{align} \frac{d^2W}{dt^2} &= \frac{W(t+2\Delta t) - 2W(t+\Delta t) + W(t)}{\Delta t^2}\\[1ex] &= \frac{\bigl[W(t+2\Delta t) - W(t+\Delta t)\bigr] - \bigl[W(t+\Delta t) - W(t)\bigr]}{\Delta t^2}\\[1ex] &= \frac{N_{t+\Delta t}(0,1)\sqrt{\Delta t} - N_{t}(0,1)\sqrt{\Delta t}}{\Delta t^2} \\[1ex] &= N_t(0,2) \frac{1}{\Delta t^{3/2}} \\[1ex] &= N_t(0,1) \frac{\sqrt 2}{\Delta t^{3/2}}. \end{align} Here I used that the sum of two standard Gaussian random variables is again Gaussian with variance $2$ (standard deviation $\sqrt 2$).

I am using forward differences due to the use of Ito's calculus.

Now you can take you differential equation and turn it into a finite difference scheme.

• Thanks, this solves my practical problem. I will try further theoretical verification, but this is what I need for now. – IRO Nov 25 '17 at 10:40
• One question that still bothers me is the correlation between the first order derivative and the second order derivative. Should they be correlated? How? – IRO Nov 25 '17 at 11:04