# Shear strain correlations from displacement Fourier transform

Currently working with molecular dynamics simulations, I would like to compute shear strain correlations in my 2-dimensional system.

How I used to do things

Accumulated shear strain at position $\vec{r}$ between times $t$ and $t + \Delta t$ is defined as $$\varepsilon_{xy}(\vec{r}, t, t + \Delta t) = \frac{1}{2}\left(\frac{\partial}{\partial x} u_y(\vec{r}, t, t + \Delta t) + \frac{\partial}{\partial x} u_y(\vec{r}, t, t + \Delta t)\right)$$ with $\vec{u}(\vec{r}, t, t + \Delta t) = \begin{pmatrix} u_x(\vec{r}, t, t + \Delta t) \\ u_y(\vec{r}, t, t + \Delta t) \end{pmatrix}$ the displacement of the particle initially at position $\vec{r}$ at time $t$ between times $t$ and $t + \Delta t$. Hence the shear strain auto-correlation function $$C_{\varepsilon_{xy}\varepsilon_{xy}}(\Delta \vec{r}, \Delta t) = \frac{\int dt \int d^2\vec{r}~ \varepsilon_{xy}(\vec{r}, t, t + \Delta t) \varepsilon_{xy}(\vec{r} + \Delta \vec{r}, t, t + \Delta t) }{\int dt \int d^2\vec{r}~ \varepsilon_{xy}(\vec{r}, t, t + \Delta t)^2}$$ which I want to compute.

One can notice that $$\int d^2\vec{r}~ \varepsilon_{xy}(\vec{r}, t, t + \Delta t) \varepsilon_{xy}(\vec{r} + \Delta \vec{r}, t, t + \Delta t) = \mathcal{F}^{-1}\{\mathcal{F}\{\varepsilon_{xy}\}^* \times \mathcal{F}\{\varepsilon_{xy}\}\}(\Delta \vec{r}, t, t + \Delta t)$$ with $\mathcal{F}$ the Fourier transform operator. Computationally speaking, this identity is very useful to quickly evaluate correlations. Up to now, I have then followed the following method:

1. Coarse-grain shear strain at positions linearly distributed on a grid from particles' positions between times $t$ and $t + \Delta t$, following J. Chattoraj and A. Lemaître, Phys. Rev. Lett. 111, 066001 (2013) (available here) and Goldhirsch, I. & Goldenberg, C. Eur. Phys. J. E (2002) 9: 245 (available here).
2. Compute shear strain correlations using Fast Fourier Transform (FFT) then inverse FFT from the obtained grid.

This method works, but is unfortunately very slow despite my best efforts to enhance my code...

How I would like to do things

There is in B. Illing, S. Fritschi, D. Hajnal, C. Klix, P. Keim, and M. Fuchs, Phys. Rev. Lett. 117, 208002 (2016) (available here with supplemental material) a method to compute shear strain correlations from displacement Fourier transform.

For that they introduce — without much explanations — the transversal and longitudinal "collective mean-square displacement" in Fourier space, respectively $C^{\perp}(\vec{q}, \Delta t)$ and $C^{||}(\vec{q}, \Delta t)$, with $\vec{q} = \begin{pmatrix}q_x \\ q_y\end{pmatrix}$ the wave vector, and then claim that (see equation 10 in supplemental material) $$C_{\varepsilon_{xy}\varepsilon_{xy}}(\Delta \vec{r}, \Delta t) = \mathcal{F}^{—1}\left\{\left(C^{\perp}(\vec{q}, \Delta t) - C^{||}(\vec{q}, \Delta t)\right)\frac{-q_x^2q_y^2}{q^2} + C^{\perp}(\vec{q}, \Delta t) \frac{q_x^2 + q_y^2}{4}\right\}(\Delta \vec{r}, \Delta t)$$

What I don't understand

First of all, I have struggled to understand the significance of $C^{\perp}$ and $C^{||}$. Inspired by F. Leonforte, R. Boissière, A. Tanguy, J. P. Wittmer, and J.-L. Barrat, Phys. Rev. B 72, 224206 (2005) (available here), I used the following definitions \begin{aligned} C^{\perp}(\vec{q}, \Delta t) &= \frac{1}{q^2} \left< ||\vec{q}\wedge\mathcal{F}\{\vec{u}\}(\vec{q}, t, t + \Delta t)||^2 \right>\\ C^{\parallel}(\vec{q}, \Delta t) &= \frac{1}{q^2} \left< ||\vec{q}\cdot\mathcal{F}\{\vec{u}\}(\vec{q}, t, t + \Delta t)||^2 \right> \end{aligned} where $\left<\right>$ denotes an average over times $t$. Using these definitions works — almost — fine, and computing shear strain is now incredibly quicker.

However I am unable to make the math and find the strain correlation expression from these definitions. Not having a solid mathematical proof also keeps from knowing if I forgot some factors or if I am completely mistaken.

If you know this proof or the correct definitions of the collective mean-squared displacements $C^{\perp}(\vec{q}, \Delta t)$ and $C^{||}(\vec{q}, \Delta t)$, or have seen either one elsewhere, this would help me a lot! Thank you!

I believe that your definitions of $C^\perp(\vec{q})$ and $C^\parallel(\vec{q})$ are correct, and hopefully I can shed some light. I think that the collective mean-square displacement tensor is defined $$\mathbb{C} = \langle \vec{u}(\vec{q})^* \, \vec{u}(\vec{q}) \rangle$$ i.e. as a dyadic product, a $2\times2$ matrix (in 2D). I'm omitting the time argument(s) throughout for clarity. Also we would normally use a hat or tilde to indicate Fourier transformed variables, but I'm omitting that as well. Now, for nonzero $\vec{q}$, this is not an isotropic tensor, even though the material (a glass) is taken to be isotropic. However, it is clear (by symmetry) that in a coordinate system based on unit vectors $(\vec{e}_\parallel,\vec{e}_\perp)$, defined so that $\vec{q}=q\vec{e}_\parallel$, and $\vec{e}_\perp$ is perpendicular to $\vec{q}$, the tensor will be diagonal. There will be a longitudinal component of $\vec{u}$, parallel to $\vec{q}$, and a transverse component, perpendicular to $\vec{q}$: $$\mathbb{C}' = \begin{pmatrix} \langle |u_\parallel(\vec{q})|^2 \rangle & 0 \\ 0 & \langle |u_\perp(\vec{q})|^2 \rangle \end{pmatrix} \equiv \begin{pmatrix} C^\parallel(\vec{q}) & 0 \\ 0 & C^\perp(\vec{q}) \end{pmatrix}$$ All the physics lies in those two functions $C^\parallel(\vec{q})$ and $C^\perp(\vec{q})$; there is no cross term. The definitions of these functions given here are the same (I believe) as the ones you took from the paper by the Barrat group.
The various coefficients in the complicated strain correlation expression, eqn (10) in the supplementary material for the Illing paper, are simply what's needed to rotate the matrix back from this diagonal form $\mathbb{C}'$ to the space-fixed form $\mathbb{C}$, and to use it to calculate the desired quantity related to the strain $C_{\varepsilon_{xy}\varepsilon_{xy}}$ rather than simply the displacement. The space-fixed $xy$ system is arbitrary, of course, but fixed; whereas you'll be considering a wide variety of $\vec{q}$ vectors. The cosine and sine of the rotation angle $\phi$ between the two coordinate systems are simply related to the components of the unit vector derived from $\vec{q}$. The conversion formula is $$\begin{pmatrix} u_x \\ u_y \end{pmatrix} = \begin{pmatrix} \cos\phi & -\sin\phi \\ \sin\phi & \cos\phi \end{pmatrix} \begin{pmatrix} u_\parallel \\ u_\perp \end{pmatrix} = \begin{pmatrix} q_x/q & -q_y/q \\ q_y/q & q_x/q \end{pmatrix} \begin{pmatrix} u_\parallel \\ u_\perp \end{pmatrix}$$ The strain is the symmetrized gradient of the displacement, so the appropriate term in Fourier space is (where $i=\sqrt{-1}$) \begin{align*} \varepsilon_{xy}(\vec{q}) &= \frac{1}{2}\left[ iq_x u_y(\vec{q}) + iq_y u_x(\vec{q})\right] \\ \langle |\varepsilon_{xy}|^2\rangle &= \frac{1}{4}\left[ q_x^2 \langle | u_y |^2\rangle + q_y^2 \langle | u_x |^2 \rangle + q_xq_y \langle u_x^* u_y\rangle + q_xq_y \langle u_y^* u_x\rangle \right] \end{align*} Substituting for $(u_x,u_y)$ using the above rotation formula is tedious but straightforward, and there is some simplification since all cross terms vanish: \begin{align*} C_{\varepsilon_{xy}\varepsilon_{xy}} = \langle |\varepsilon_{xy}|^2\rangle &= \frac{q_x^2q_y^2}{q^2} \langle | u_\parallel |^2 \rangle +\frac{q_x^4+q_y^4-2q_x^2q_y^2}{4q^2}\langle | u_\perp |^2\rangle \\ &= \frac{q_x^2q_y^2}{q^2} C^\parallel +\frac{q_x^4+q_y^4-2q_x^2q_y^2}{4q^2} C^\perp \end{align*} Bearing in mind that $q_x^2+q_y^2=q^2$, this is identical with the formula that you were concerned about.