Randomly stopped dynamics of $x(t)$: How can I find $\text{var}\{ x(t) \}$?

Question

Consider the simple dynamical equation

$$ \dot{x}(t) = u H(t-\tau),$$ where the timescale $\tau$ is an exponentially distributed random variable $\tau \sim \omega \exp\{\omega \tau\}$ and $H(t) = 1-\theta(t)$. That is, $$ H(t) = \begin{cases} 1 & \text{if} & t<0 \\ 0 & \text{if} & t> 0\end{cases}. $$ This describes a particle which moves with velocity $u$ for a random time $\tau$ at which point its motion is cut off.

I can integrate for its position conditional on the random variable $\tau$ (given $x(0)=0$): $$ x(t) = u \int_0^t dt' H(t-\tau),$$ which gives the mean value $$ \langle x(t) \rangle = u \int_0^t dt' \langle H(t'-\tau) \rangle$$ I believe I can evaluate the ensemble average over $H(t'-\tau)$ using the exceedance probability of $\tau$. Since $H(t-\tau)=1$ when $t<\tau$, $$ \langle H(t-\tau) \rangle = \text{Prob}(\tau > t') = \int_{t'}^\infty \omega \exp\{-\omega \tau\} = \exp\{-\omega t' \},$$ giving mean position $$ \langle x(t) \rangle = \frac{u}{\omega}(1-\exp\{-\omega t\}).$$ When the mean cut-off time $\langle \tau \rangle = 1/\omega \rightarrow \infty$, I obtain the expected result in absence of the cut-off: $$ \langle x(t) \rangle \approx ut + O(\omega),$$ so I believe my result is correct.

Now for my questions. First, have I done anything wrong? Second, given I'm on the right track, how can I compute the second moment of the position? That is, I'm seeking $\langle x(t)^2 \rangle$ which I believe should be given by:

$$ \langle x(t)^2 \rangle = u^2 \int_0^t dt_1 \int_0^t dt_2 \langle H(t_1-\tau) \rangle \langle H(t_2-\tau) \rangle \\ = u^2 \int_0^t dt_1 \int_0^t dt_2 \text{Prob}(\tau>t_1)\text{Prob}(\tau>t_2) \\ = \Big(u \int dt \text{Prob}(\tau > t) \Big)^2 \\ = \langle x(t) \rangle^2$$ but this is unexpected: how is the second moment just the square of the mean? This implies a zero variance... I expected a transition in the second moment when $t$ surpasses $1/\omega$, since this is the point when most particles in the ensemble will have been trapped.

Any guidance is appreciated !

Answer 1

Almost there. The two $\text{Prob}$'s are correlated: $$u^2 \int_0^t dt_1 \int_0^t dt_2 \langle H(t_1-\tau) H(t_2-\tau) \rangle = u^2 \int_0^t dt_1 \int_0^t dt_2 \text{Prob}(\tau>\max(t_1, t_2)) = \dots = \frac{2u^2}{\omega^2}\left(1-(1+\omega t)\exp\{-\omega t\}\right)$$

Which you can verify with a quick Python script if so inclined:

import numpy as np
omega = 2.5
t = 0.75
u = 4.
N = 1000000

tau = np.random.exponential(1./omega, N)        
x = u * ((tau > t) * t + (tau < t) * tau)

print (np.mean(x), u/omega * (1 - np.exp(-omega*t)))
print (np.mean(x**2), 2 * u**2/omega**2 * (1 - (1+omega*t) * np.exp(-omega*t)))

Answer 2

I revisted this problem to solve for the probability distribution of position. I'm posting this as an answer in case anyone is ever interested.

The equation of motion for the particle's displacement is $$ \dot{x} = u \Theta(\tau-t).$$ For a particular travel time, this integrates for $$ x(t) = u \big[ t - (t-\tau)\Theta(t-\tau) \big]. $$ The ensemble probability distribution of position can be defined as $p(x,t) = \langle \delta(x-x(t) \rangle, $ where the average is over all possible travel times $\tau$. Taking the Laplace transform of this definition and averaging over $\tau$ provides $$\tilde{p}(s,t) = \frac{\omega}{\omega + s u} + \frac{s u }{\omega + s u} e^{-(\omega + s u) t}.$$ Inverting this gives $$ p(x,t) = e^{-\omega t}\delta(x-ut) + \frac{\omega}{u}e^{-\omega x/u}\Theta(u t- x).$$

One can then calculate any statistical moments as $\int dx x^k p(x,t)$.