# NEW: dissipative EDMD
###### tags: `Dissipative self-assembly`

<style>.markdown-body { max-width: 1950px; }</style>
---

---

### Owners (the only one with the permission to edit the main test)

EF, AP, FS, GF

---

---

# I'm creating a new page, the other one has reached maximum lenght...

## LIST OF RELATED PAPERS:

| Paper name | SIMILARITY/WHAT THEY SAY | DIFFERENCES/WHAT THEY DONT SAY | NOTES|RELATED/INTERESTING CITATIONS|
| -------- | -------- | -------- |------|------|
| [Hydrodynamics of random-organizing hyperuniform fluids](https://www.pnas.org/doi/10.1073/pnas.1911596116) (theory + simulations) | 1. Competition between energy injection at collision and gamma <br />2. Manna universality class.<br /> 3. Hyperuniformity at transition and above with hydrodynamical theory for above. Hyperuniformity here because of $\gamma$. <br /> 4. Everything is simpler for them, but they have the mean free path argument | 1. They don't have a realistic granular collision (they have $\Delta E$ = cst) <br /> 2.They don't have synchronization of course | 1. Hyperuniformity in this system is caused by conservation of momentum at collision modeled in the fluctuating NS equation <br /> 2. They are aware that their model could apply to the quasi 2D setup: *"Thus, we anticipate that other driven-dissipative systems studied previously can also produce the fluidic HU: for example, the vibrated granular monolayer that allows the kinetic energy transferring from vertical to horizontal freedom degree*"|1. [Linear and nonlinear hydrodynamics of low-friction adsorbed systems](https://journals.aps.org/pra/abstract/10.1103/PhysRevA.26.1735) <br /> 2. [Enhanced hyperuniformity from random reorganization](https://www.pnas.org/syndication/doi/10.1073/pnas.1619260114)|
|[Hyperuniform states generated by a critical friction field](https://arxiv.org/pdf/1805.07408.pdf) (experiments + idealized simulations + theory)|1. Exactly same quasi 2D setup. <br /> 2. They find hyperuniformity at the liquid-solid transition (BUT not far above or below, see next column)|1. They never talk about phase transition/small density regime <br /> 2. They have a setup that can [lead to square crystals](https://arxiv.org/pdf/cond-mat/0312232.pdf) due to high height of plate<br />3. They don't find hyperuniformity far from the coiexistence region while we expect to do? For very large system size at least. - To them, the hyperuniformity is explained by: *"Prior to the transition, patches of the solid phase form, with length scales and mean lifetimes that diverge critically at the transition point. When reducing the wavenumber, density fluctuations encounter increasingly more patches that block their propagation, resulting in a static structure factor that tends to zero for small wavenumbers at the critical point, which is a signature of hyperuniformity."* While we expect it even in the liquid 4. Their theory of hyperuniformity is really based on the coupling between the density field and the local $q_4$. They model the dissipation with roof collision by a diffusive equation for the density field, instead of advection. And use a $D$ dependent on $q_4$.| 1. They have subdiffusion for the solid part of the system (which coexists with the fluid). I hope I will just see caging (plateau) in the fully cristallized system to prove that I have a real 2D crystal! <br /> 2. They do the following fit: $S(k, \phi) = S_0(\phi) +S_1(\phi)(ka)^\alpha$ to $S(k)$ with $k$ small to prove that when $\phi$ goes to the critical packing fraction of the solid liquid, $S_0/S_1$ goes to 0 with $k^\alpha$ which means hyperuniformity. In any case I find it a bit dubious. They don't do any scaling analysis. But their analysis seems careful...|1. [Liquid–solid-like transition in quasi-one-dimensional driven granular media](https://www.nature.com/articles/nphys884)(Nice theoretical treatment of liquid-solid phase transition/Coarse graining)|

Reference about synchronization:
[Dynamics of noncohesive confined granular media](https://arxiv.org/pdf/1601.04518.pdf)

[Phase behavior and perturbation response of quasi-2D driven granulates employing molecular dynamics simulations](https://opus4.kobv.de/opus4-fau/files/23690/thesis-ThomasSchindler-published.pdf)

[Sudden chain energy transfer events in vibrated granular media](https://arxiv.org/pdf/1101.5365.pdf)

______

I did a scan of (amp, height, freq) to see if I could obtain a high frequency/high amplitude/small height (I think, it's better for long range order to not have particles jumping high) asynchronized system.


Blue is asynchronization. The map is very coarse but we can see some async zones even in the high frequency limit. I think that high frequency doesn't prevent asynchronization but that the chaotic zone are smaller in a fractal way

(here is an updated 1000Hz map:

)

The small amp zone without synchronization corresponds to a zone where, if we start with no xy velocity, the motion of the particles is really dependent on the initial conditions and some of them might never touch the ceiling. I don't know if we want to place ourself here. Especially since we want a lot of collisions with the top plate for a strong drag.

______________

Here are some results of the MSD vs N for freq = 50 in the async region:


We see that at high N we might have a plateauing of the MSD (or at least we are losing the log(N) scaling of the MSD vs N). But we need to go at really high N... Perhaps to reach a $q_{min}\in$ the hyperuniform zone.

Nonetheless, it doesn't seem to be that dependent on $\gamma$ for the EDMD model, so...:

Here is with $\gamma = 0.01$:
We might see a logarithmic increase in N. And keep in mind that with lammps, $\phi_c\sim 0.1$ while for these values, in EDMD, $\phi_c$ would be close to $0.2$. So there still hope.

_______
### [27/09/23] WE MIGHT HAVE MERMIN WAGNER BREAKDOWN WITH LAMMPS:
The legend is wrong, it's g = 495 amp = 0.3 h = 1.25
Somehow the system starts diffusing in the region where all MSD plateau. So... I don't know what to think of this. I will let the simulations run longer...

Here is my two cents guess: We have a breakdown of mermin wagner, but this configuration is not a crystal. So we still see diffusion. Nonethess, we see a residue of the suppression of phonons (which are still somehow there, even wihtout a crystal.)
Now, I will run a simulation with parameters I know lead to a true crystal (because MSD is constant) with a very large system size. Above, we kind of already see the good trend with f = 53Hz

### [28/09/23] synchronization induced coexistence:

| f = 600hz Amp = 0.2377sig h = 0.13sig | f = 600hz Amp = 0.2377sig h = 0.13sig |
| ---- | ---- |
|  |  |

Colored by position:


Will the system die? Is it known?

### [04/10/23] WTFFFFF?

$f = 825$Hz
| WHEN PARTICLES DONT TOUCH THE CEILING AND ARE ASYNC BECAUSE OF THIS: | WHEN PARTICLES TOUCH THE CEILING AND ARE ASYNC DUE TO CHAOTIC ATTRACTOR: |
| ---- | ---- |
|  |  |

this is UNEXPECTED!

______


f = 53, A 0.053, h= 1.47
We can't reach high packing fraction without clustering..

Interestingly, this cliustering is not visible on the longitudinal velocity velocity structure factor:

| $\langle u_\parallel u_\parallel\rangle$ | $\langle\rho\rho\rangle$ |
| -------- | -------- |
|  |  |

______
Recall that: $$F_{wall} L_x = \lim_{t\to \infty}\frac{m}{t}\left[
\sum_{col} \delta x^{col} \delta{v}_x^{col} +
\sum_{kick}x^{kick} \frac{\gamma}{m}v_x + \sum_{all} \sum_i v_x^2\delta t_{last}
\right].$$
in the limit where we apply the damping by a time interval $dt_{noise}$ (kick). We see that the damping will be averaged out (we hope that this still holds in the continuous limit: $dt_{noise}\rightarrow 0$) so that the pressure with $\Delta$ + $\gamma$ can be found using only the delta pressure:

$$p = \dfrac{NT}{L^2}-\frac{1}{2L^2}\left\langle\sum_{i< j}\vec{r}_{ij}\cdot\vec{F}_{ij} \right\rangle$$

which gives:

$$p = \dfrac{NT}{L^2}-\frac{1+\alpha}{2tL^2}\sum_{coll-ij}\dfrac{m_im_j}{m_i + m_j}\vec{r}_{ij}\cdot\vec{v}_{ij}+\dfrac{\Delta}{tL^2}\sum_{coll-ij}\dfrac{m_im_j}{m_i + m_j}|r_{ij}|$$


With $\Delta$ model, the transition is shifted to lower density compared to eq. It seems to also be a first order one sine the loop disappear with size (to be sure we'd have to check is the area under and above the maxwell construction decays as $1/\sqrt{N}$ as expected from a loop emerging from an interface) (also keep in mind that the system has probability not perfectly thermalized..). The coexistence for the $\Delta$ model is also well seen in the snapshots if we color the particles by their local $q_6$ angle:

| $\phi = 0.679$ | $\phi = 0.687959$ | $\phi = 0.6989$ | $\phi = 0.71$ |
| -------- | -------- | --- | -------- |
|  |  |  |  |

We also observe the hexatic phase at intermediate $\phi$ with an exponentially decaying translational order parameter as in the disordered phase (while the high dnesirty case exhibit power law decay behavior):
Note that $\phi= 0.7075$ is clearly out of the coexistence zone estimated from the lever rule. Thus the system is clearly not in a coexisting phase.

By looking at the snapshots, both looks crystalline:

| $\phi = 0.7075$ | $\phi = 0.71$ |
| -------- | -------- |
|  |  |

However, the lowest density one is not a true crystal since it's in hexatic phase, thus orientational order is only quasi long range but we can't really expect to see it from the orientational correlation function since the decay is very small:

__________

Assuming that in the steady state, the rate at which different tiles in a quasi crystals can be described by a Markovian process. One can ask if the equilibrium quasi-crystal and the non equilibrium quasi crystals follow a detailed balance principle for the stochastic matrix or not.
To do this, one has to somehow track the tiles with time and deal with the fact that the number of tiles is not conserved with time. I do this by tracking the centroid of each tiles at each snapshot. At time $t$ I have a ensemble of $N$ centroid for different tiles. At time $t +dt$, I have an other of ensemble of $M$ centroid for different tiles. Then I use a linear sum assignment algorithm on the cost matrix defined as the distance between every centroid at time $t$ and $t + dt$ of size $N\times M$. This gives us a an array of min$(N, M)$ elements corresponding to the index of the min$(N, M)$ atoms at $t + dt$ that matches the positions of the centroid at $T
If the transition matrix is $M_{ij}$ with $i$ and $j$ corresponding to some tiles (for example small triangle with a given orientation, some other triangle with an other orientation, ambigous tile, etc) and $P_i$ is the probability of observing the tile $i$. Now can try to see if the relation $$M_{ij}P_i = M_{ji}P_j$$ is verified in different systems.
The first issue is of course the tracking. If a tile changes it will mostly becomes a defect and then an other one. So, I think it's preferable to look at long time scale for these kind of things but we can't track a tile change on long time scale...

______

### RESULTS:

Here are the transitions matrix from one tile to another:

For equilibrium qc8 and delta model qc8.
In order of appearance, the symbols correspond to the following tiles:
- equilateral triangle
- defect
- ambiguously oriented tiny square
- tiny squares (2)
- ambiguously oriented large square
- large squares (2)
- ambiguously oriented isocele triangle
- isocele triangle (8)

The matrix doesn't change much eq vs non eq. We see that defects often change to an ambiguous tiling and that isocele triangles often rotate to their neighboring possible tiling (maybe it's an issue of the tracking..)

____

If we look at the distribution of $M_{ij}P_i - M_{ji}P_j$ for every tiles (except i = j) we find the distribution of the deviation from detailed balance:

It's not gaussian (probably) but at least it's centered in 0. No clear difference between neq and eq can be seen here. Maybe Eq. is a bit more trimodal....

____

To ask wether we are indeed in a detailed balanced system or not, for these parameters, we can do a statistical analysis. That is, see how $M_{ij}P_i - M_{ji}P_j$ evolves when we increase the statistics (the number of frames on which we track):

This is the evolution of the variance of $M_{ij}P_i - M_{ji}P_j$ (thus how much the distribution above get closer to a dirac delta in 0). We see that for both eq and neq systems, the variance diminishes as 1/(number of snapshots) (green line), as expecgted from the central limit theorem, I think? Thus we can expect that in the limit of perfect statistics, we indeed observe detailed balance for both systems.

We observe also a decrease for the maximum (or minimum) value $M_{ij}P_i - M_{ji}P_j$ with a decrease as $1/\sqrt(N)$:

______
NOW FOR THE DEPRESSING RESULTS:

This is for a crystal crystallizing:
{%youtube Augc18KpOjc %}




### Continuous transition theoretically.

$\Delta$ model + drag. Collisions bring an energy $m\Delta^2$ into the system. And the two particles will lose a energy equal to: $mv_0^2e^{-2\gamma\tau_f}$ where $v_0$ is the exit velocity and $\tau_f = -\log(1 - \gamma l(\phi)/v_0)/\gamma$ is the time to reach a neigboring particle, on average, if the mean free path is $l(\phi)$. In the steady state we have:$$\Delta^2 = v_0^2(1 - e^{-2\gamma\tau_f})$$ which gives us $$v_0 = \dfrac{\Delta^2 + \Delta^2_c}{2\Delta_c}$$ where we defined $\Delta_c = \gamma l(\phi)$ which is the critical value of $\Delta$ for a given $\phi$ and $\gamma$ according to the mean free path argument. Note that we indeed find that the exit velocity is equal to $\Delta$ at the critical point. As expected. The temperature is defined as: $$ \begin{split}
T &= \dfrac{1}{2}m\dfrac{1}{\tau_f}\int_0^{\tau_f}|v^2(t)|\mathrm{d}t\\
&=\dfrac{v_0^2}{4}\dfrac{(1 - \gamma l(\phi)/v_0)^2 - 1}{\ln (1 - \gamma l(\phi)/v_0)}.\\
& = \frac{\Delta^2}{4 \log \left(1-\frac{2 \Delta_c^2}{\Delta^2+\Delta_c^2}\right)}
& = -\frac{\Delta^2}{4\log((\Delta^2 - \Delta_c^2)/(\Delta^2 + \Delta^2_c))}
\end{split}$$

Trying to fix the logarithmic power law by defining

\begin{split}
T(T') =& \dfrac{1}{2}m\int_{-\infty}^{\infty}d\mathbf{v_0}\dfrac{P(\mathbf{v_0}, T')}{\tau(\mathbf{v_0})}\int_0^{\tau(\mathbf{v_0})}dt\mathbf{v}^2(t)\\
\end{split}
with $P(v_0)\sim \dfrac{(v_0 - \Delta)}{T'}e^{-\dfrac{1}{2}(v_0-\Delta)^2}$ does not fix the issue.

_____
### Neq (orange) and equilibrium defects/ambiguous tiles area proportion
when steady state

when in the process of reaching the steady state (the only meaningfull things are the P_defects and P_miss)

______
Also concerning hyperuniformity due to the critical point. We also see it from the MSD, but it's hard to tell whether it's due to the decrease of temperature or not...

This is quasi2D EDMD. I keep $\phi_c$ fixed but change the closeness to criticality by changing the tangential coefficient of restitution of the walls!

__________

Nice image for turbulence:


This is interesting because it shows the thermalization of the system (note hoever that $k^{-3}$ is not a usual $k^{-3}$ cascade and is more likely due to the initial vortex. In any case, it would be interesting to see what happens with our granular system ($\Delta$ with or without $\gamma$). Altough, it's harder to correctly do the analysis since energy is not conserved.

_____

A bit less expected to me, the same kind of evolution is obtained with the $\Delta$ model even though it does not conserve energy!!! With the $\Delta$ + $\gamma$ model, the structure is quickly dissipated and we easily lose the $k^3$ spectrum...
_____________

### [28/11/23] Rheology
IMPORTANT: $\gamma$ is the same whatever the mass is!!!!!!!!!

Big intruder in bath. We probably want a passive intruder as Giuseppe said, otherwise it might be pushing particles a bit too much (due to conservation of momentum) but it looks fine even with a $\Delta$ for the intruder: 

Now, we want to study the autocorrelation velocity of the tracer.
First, let's try to get right the autocorrelation (which we never get right from the first time, this is like fourier transform :)). I'm generating an AR(1) process and comparing the theo prediction with different realizations -- is Andrea's code (which works) and solid lines is numpy code which is fast. The two seem to agree well.

Strangely however, increasing the number of realization is not helping much.

but increasing the sample size is helping, $T$ is the size of the sample (w t f?):

A nice thing to learn from that is that, to observe backscattering, we should be sure that the well is not changing with the duration of the observation:

_______
Puglisi:


Some other guy 3D:

| Column 1 | Column 2 | Column 3 |
| -------- | -------- | -------- |
|  |  |  |

________

As expected Eq is close to restituion coefficient + $\Delta$. We can't see any backscattering (I think)... Maybe we can?

We also see that at low density, the profile is very exponential:

We have backscattering close to the transition(TYPO, READ: $t^{-1}$):

This is really due to the ansorbing state/avalanches and not due to high packing fraction or small $\alpha$(TYPO, READ: $t^{-1}$):
_________

3d hard spheres (Velocity Autocorrelation Functions of Hard-Sphere Fluids: Long-Time Tails upon Undercooling):

_____

If the tracer has a mass close to the one of the other particle, it is "easy" to observe backscattering in the liquid:
If the tracer is massive, we can only see the backtracking in the solid.
_________

Evolution of quasi crystal at $\phi = 0.86$. Blue is eq and orange is delta model.
To bne compared with same simu but at $\phi = 0.85$


So $\Delta$ really looks like higher pressure. Also, here is the evolution of the energy.
_____________

### 15/12/23 very tired me, doing tiresome bookeeping analysis!

$\phi = 0.6$, $r_{intruder} = 5r_{bath}$, $m_{intruder} = 5^2m_{bath}$

Hard disk vs langevin: HUGE DIFFERENCE BACKSCATTERING AT LARGE COUPLING WITH BATH (wtf?)

Pure hard disk equilibrium and $\Delta$: almost no difference

Equilibrium and $\Delta$ + $\gamma$: huge difference (backscattering)
this is independent on $\alpha$ or $\phi$:

_____________________

Now to get backscattering with eq hard disks wihtout noise $\phi = 0.6$:

| Density intruder constant changing size (and mass consistently)| radius intruder constant (mass varying). But mass changing | Column 3 |
| -------- | -------- | -------- |
|  | same size as bath  | same mass as bath  |
| | intruder x5 size as bath | mass 10x compred to bath |

___________
More on langevin noise:

$\phi = 0.2$

_______

Puglisi model of velocity autocorreelation:

$$m\ddot{\mathbf{u}}_{k}(t) = - K\omega_k^2\mathbf{u}_{k}- \int_{-\infty}^{t} \mathrm{d} t' \Gamma_k(t -t')\dot{\mathbf u}_{k}(t') + \mathbf{\mathcal{E}}_k(t)$$
with:
$$ \Gamma(t) = 2\gamma_0\delta(t) + \dfrac{\gamma_1}{\tau}e^{-t/\tau}$$
and $$\langle\mathcal{E}_k^{\alpha}(t)\mathcal{E}_q^{\beta}(0)\rangle = 2\omega_k^2\delta_{q, -k}\delta^{\alpha, \beta}\left(\gamma_0T_0\delta(t)+\dfrac{\gamma_1}{2\tau}T_1e^{-|t|/\tau}\right)$$

Solved through Laplace transform or markovianization.

Some equilibrium results: $T_0 = T_1$

| case | picture |
| -------- | -------- |
| $\gamma_1 = 0$ |  |
| $\gamma_0 = 0$ |  |
| $\gamma_0\neq0$ and $\gamma_1\neq0$ (underdamped)|  |
|$\gamma_0\gg\gamma_1$ (overdamped)|
|$\gamma_0\gg\gg\gg\gg\gamma_1$||
_______

Note to future be: beware of derivative and fourier convention in mathematica :) with -iw and (FourierParameters->(1, 1)):

END OF THAT
_____

As Andrea said, to see the last result it is easier to look at the power spectrum density:

$$\langle v(w)v(-w)\rangle=\frac{2(\gamma_1 T_1+\gamma_0 T_0(1 + \tau^2\omega^2))}{(\gamma_0 + \gamma_1)^2 + (1-2\gamma_1\tau + \gamma_0^2\tau^2)w^2+\tau^2w^4}$$
_____
### $\tau \rightarrow 0$

In the limit $\tau = 0$ we recover the simple lorentzian:

$$\langle v(w)v(-w)\rangle=\frac{2(\gamma_0T_0 + \gamma_1T_1)}{(\gamma_0 + \gamma_1)^2+w^2}$$
And $v^2 =\dfrac{1}{2\pi}\int dw \langle v(w)v(-w)\rangle = \dfrac{\gamma_0T_0 + \gamma_1T_1}{\gamma_0 + \gamma_1}$

### $\gamma_1 \rightarrow 0$
If we put $\gamma_1 = 0$ before computing the PSD we obtain:
$$\langle v(w)v(-w)\rangle=\frac{2\gamma_0T_0}{\gamma_0^2+w^2}\Rightarrow\langle v(t)v(0)\rangle = T_0 e^{-\gamma_0 t}$$
Same thing if we take the limit after (the functions are the same):
$$\langle v(w)v(-w)\rangle=\frac{2\gamma_0T_0(1+\tau^2w^2)}{\gamma_0^2+(1 + \gamma_0^2\tau^2)w^2 + \tau^2w^4}$$

and $v^2(0) = T_0$ because the FDT is respected!

### $\gamma_0 \rightarrow 0$

The PSD is:
FDT is respected.

### $\gamma_0 \rightarrow \infty$

The PSD is:
$$\langle v(w)v(-w)\rangle=\frac{2\gamma_0 T_0(1 + \tau^2\omega^2)}{(\gamma_1+\gamma_0)^2 + \gamma_0^2\tau^2w^2+\tau^2w^4}$$

We do not recover a lorentzian! Hence, in the limit of very srong delta damping, we still see the effect of memory

FDT is off course, respected.

______

Back to full system case: the autocorrelation function is:
$$\dfrac{\langle v(t)v(0)\rangle}{\langle v(0)^2\rangle}=\frac{\left(\left(\sqrt{(\gamma_0 \tau-1)^2-4 \gamma_1 \tau}-\gamma_0 \tau+1\right) e^{\frac{t \sqrt{(\gamma_0 \tau-1)^2-4 \gamma_1 \tau}}{\tau}}+\sqrt{(\gamma_0 \tau-1)^2-4 \gamma_1 \tau}+\gamma_0 \tau-1\right) \exp \left(-\frac{t \left(\sqrt{(\gamma_0 \tau-1)^2-4 \gamma_1 \tau}+\gamma_0 \tau+1\right)}{2 \tau}\right)}{2 \sqrt{(\gamma_0 \tau-1)^2-4 \gamma_1 \tau}}$$

We see that it goes from underdamped to overdamped when: $(\gamma_0\tau-1)^2-4\gamma_1\tau =0$

Moreover, when $\gamma_0\rightarrow\infty$, a singular expansion gives:

$$\dfrac{\langle v(t)v(0)\rangle}{\langle v(0)^2\rangle}=\left(\frac{\gamma_1}{\gamma_o^2 \tau}+1\right) e^{\displaystyle{\dfrac{\gamma_1 t}{\gamma_o^2 \tau^2}+\dfrac{\gamma_1 t}{\gamma_o \tau}-\gamma_o t}}-\dfrac{\gamma_1 e^{-\frac{t}{\tau}}}{\gamma_o^2 \tau}$$

We see that there is a backscattering and that, as $\gamma_0$ approaches $\infty$ the smallest value of the autocorrelation goes to 0 (losing the backscattering) but also get closer to the y axis. Making it easier to see :) because at small time we have better statistics.
____
Out of equilibrium we obtain the same kind of results except that $\langle v(t)v(0)\rangle\neq\langle\delta v(t)\rangle$
The delta noise is of course not needed for the oscillations but it is needed to have a genuine non-equilibrium dynamics (with the damping form choosen above).

Note also that oscillation is the rule rather than the exception with (unphysical) noises.

Still equilibrium:

| case | Pic |
| -------- | -------- |
|noise = constant||
|power law correlated noise expo < 0.5 ||
| power law correlated noise exponent < 0.5 + dirac delta (kills the oscillations) |  |

_____________
Markovianization

A funny thing is that, any multiply exponentially correlated process can be converted to a markovian process as presented during the sieste by Andrea (because of Doob's theorem, it seems that the Orenstein Uhlenbeck process is the only exponentientially correlated continuous random process with zero mean. So the equivalent markovian system has to be a Ornstein-Uhlenbeck process. See: Stochastic Processes and Applications Grigorios A. Pavliotis
Diffusion Processes, the Fokker-Planck and Langevin Equations).

It is way easier to convert an equation with memorry to a markovian process than doing the laplace transform (even if it is equivalent. This is a nice way to obtain easily formula for hard integrals by the way :^).

More generally the following equation:


can be mapped to:

If the memory is exponential like:

or if the laplace transfrom of the memory has a continuous fraction form as in the original Mori's article;

Notably, power law correlated noise can be approximated by a continuous fraction or with an infinite sum of exponential (see: Fractional Kinetics in Kac–Zwanzig Heat Bath Model and https://arxiv.org/pdf/1804.00202.pdf) as andrea already noted. Some people look at these kind of models to study the (sub-supra)diffusive behavior of **equilibrium** system with power law correlation in time. One might wonder what happens if we assign a different temperature to each exponential making the power law.

_____________

### Entropons from memory?!

$$m\ddot{\mathbf{u}}_{k}(t) = - K\omega_k^2\mathbf{u}_{k}- \int_{-\infty}^{t} \mathrm{d} t' \Gamma_k(t -t')\dot{\mathbf u}_{k}(t') + \mathbf{\mathcal{E}}_k(t)$$
with:
$$ \Gamma(t) = \omega_k^ 2(2\gamma_0\delta(t) + \dfrac{\gamma_1}{\tau}e^{-t/\tau})$$
and $$\langle\mathcal{E}_k^{\alpha}(t)\mathcal{E}_q^{\beta}(0)\rangle = 2\omega_k^2\delta_{q, -k}\delta^{\alpha, \beta}\left(\gamma_0T_0\delta(t)+\dfrac{\gamma_1}{2\tau}T_1e^{-|t|/\tau}\right)$$

Equivalent to:

\begin{equation}
\begin{split}
\dot{v}_k(t) &= \dot{u}_k(t)\\
\dot{u}_k(t) &= -\gamma_0\bar\omega_k^2v_k(t)-K\bar\omega_k^2u_k(t)+\Omega_k(t)\\ &~~~+ \sqrt{2\gamma_0\omega_k^2T_0/m}\zeta(t) \\
\dot{\Omega}_k(t) &= -\frac{\gamma_1\bar\omega_k^2}{\tau}v_k(t) - \frac{1}{\tau}\Omega_k(t) + \sqrt{2\gamma_1\bar\omega_k^2T_1/m\tau^2}\eta(t),
\end{split}
\end{equation}
With better notation and forgetting the velocity equation which won't play any role:
\begin{equation}
\dot{\begin{bmatrix}
u_k\\
\Omega_k
\end{bmatrix}}=
\begin{bmatrix}
-\gamma_0\bar\omega_k^2v_k ~~~~~~\Omega_k\\
-\gamma_1\bar\omega_k^2v_k/\tau ~~~-\Omega_k/\tau\\
\end{bmatrix}
+
\begin{bmatrix}
\eta_k(t)\\
\zeta_k(t)
\end{bmatrix}
\end{equation}
Or equivalently:

$$\dot {\mathbf z}_k = \mathbf F[\mathbf z_k] + \mathbf \xi_k(t)$$
with $\langle\eta_k(t)\eta_{k'}(t)\rangle = 2\gamma_0\bar\omega_k^2T_0/m\delta(k + k')\delta(t - t')$ and $\langle\zeta_k(t)\zeta_{k'}(t')\rangle = 2\gamma_1\bar\omega_k^2T_1/\tau^2m\delta(k + k')\delta(t - t')$ and $\xi_{k,0} = \eta_k$ and $\xi_{k, 1} = \zeta_k$.
$\langle\xi_k(t)\xi_k^t(s)$ is diagonal and its (functional) inverse is called $L_k(t-s)$. We will forget the $^t$ transpose in what follow for simplocity.

Probability of a path going from -infinity to infinity:

$$P[u_{\{k\}}]\propto e^{-\dfrac{1}{2}\displaystyle{\int_{-\infty}^\infty dt\sum_k}\xi_k(t) L_k(t-s)\xi_{-k}(s)}$$

Since $L_k(t-s) = \tilde L_k\times\delta(t-s)$ with $\tilde L_k$ the appropriate diagonal matrix, we obtain:

$$
\begin{split}
P[u_{\{k\}}]&\propto e^{-\dfrac{1}{2}\displaystyle{\int_{-\infty}^\infty dw\sum_k}\xi_k(w) \tilde L_k\xi_{-k}(-w)}\\
&\propto e^{-\dfrac{1}{2}\displaystyle{\int_{-\infty}^\infty dw\sum_k}(-iw\mathbf z_k - \mathbf F[\mathbf z_k]) \tilde L_k (iw\mathbf z_{-k} - \mathbf F[\mathbf z_{-k}])}
\end{split}
$$

The entropy production is proportional to the log of the ratio of the probaiblity of a trajectory with its reversed one ($v\rightarrow-v$, $\Omega\rightarrow\Omega$ and $w\rightarrow-w$).

We find $\langle\sigma\rangle = m\dfrac{T_1-T_0}{T_0T_1}\langle \Omega_k(w) u_{-k}(-w)\rangle + c.c.$

Which is Andrea's result ! $m=1$ and $K = 1$ because we all hate the mass. The sum of average is:.
Here is a plot of the entrop poduction as a function of w and k.

The ratio of the entropy production with $<uu>$ is:

which does not take a nice form.

In the origianl entropon article, with two separates temperature acting on two different axis, it looks like we can separate well the displacement correlation function to obtain clear entropons:
Otherwise, when two temperatures act on a single degree of freedom, you're screwed, I think and cannot obtain nice entropons!

_________

#### Entropy production in an equilibrium like system?!
Let's do something simpler:
\begin{equation}
\begin{split}
\dot{v}(t) &= -\gamma_0^2v(t)+\Omega(t) + \sqrt{2\gamma_0T_0}\zeta(t) \\
\dot{\Omega}_k(t) &= -\frac{\gamma_1}{\tau}v(t) - \frac{1}{\tau}\Omega(t) + \sqrt{2\gamma_1T_1\tau^2}\eta(t),
\end{split}
\end{equation}

This describe a free particle with two gaussian thermostats $T_0$ delta correlated and $T_1$ exponentially correlated with decay rate $\tau$.

The detailed entropy production for this system is:

$\sigma(w) = 2\dfrac{(T_1-T_0)^2}{T_0T_1}\dfrac{\gamma_1\gamma_1}{(\gamma_0 + \gamma_1)^2 + (1-2\gamma_1\tau + \gamma_0^2\tau^2)w^2 + \tau^2w^4}$

And for $\tau = 0$, the total entropy production is:

$$\sigma = \dfrac{1}{2\pi}\int_{-\infty}^\infty \sigma(w) = \dfrac{(T_1-T_0)^2}{T_0T_1}\dfrac{\gamma_0\gamma_1}{\gamma_0 + \gamma_1}$$

Which is simply, the Fourier law of heat transfert! This is very neat, because the limit $\tau = 0$ corresponds to two delta correlated bath at different temperature. If we did the computation of the entropy production,, starting from $\tau = 0$ or equivalently, two delta correlated bath, we would have found an entropy production equal to 0 because the system would have juste felt a single bath at temperature $T_{eff} = (\gamma_0T_0 + \gamma_1T_1)/(\gamma_0 + \gamma_1)$ but here, with the limit of 0 memory the system/math, recalls that the two bathg are well separated and can exchange heat!!!!!!! Very cool! (This result can be also derived with the fokker planck where fluxes from different thermostats can be easily divided into different contribution, but to my knowledge this is very complicated to opbtain this result from a path integral entropy production approach using purely a langevin equation without taking some kind of limit and strarting directly from two delta correlated bath)

This also immediately shed some dark light upon the desesperate entropons! Because, in the limit of an equilibrium-**like** system where excitations should be purely phononic, a non zero entropy production can exist.
________

Now let's go back to the system with a lattice (there might be an issue with the m, so in the next formula, m = 1 :) ): $$\sigma(w, k)=\dfrac{1}{T_0T_1}\frac{4 \gamma_1 \gamma_0 w^2 \omega_k^4 (T_1-T_0)^2}{m^2 \left(\tau^2 w^2+1\right) \left(w^2-K \omega_k^2\right)^2-2 \gamma_1 m \tau w^2 \omega_k^2 \left(w^2-K \omega_k^2\right)+w^2 \omega_k^4 \left(\gamma_1^2+2 \gamma_1 \gamma_0+\gamma_0^2 \left(\tau^2 w^2+1\right)\right)}$$

When $\tau = 0$ we obtain:

$$\sigma(w, k)=\dfrac{1}{T_0T_1}\frac{4 \gamma_1 \gamma_0 w^2 \omega_k^4 (T_1-T_0)^2}{\left(w^2-K \omega_k^2\right)^2+w^2( (\gamma_0 + \gamma_1)\omega_k^2)^2 }$$

From which we also obtain the exact same entropy production as for the off lattice case:

$$\sigma(k) = \dfrac{1}{2\pi}\int_{-\infty}^\infty \sigma(k, w) = \dfrac{(T_1-T_0)^2}{T_0T_1}\dfrac{\omega_k^2\gamma_0\gamma_1}{\gamma_0 + \gamma_1}$$

Let's make a tab!

| Lattice model | $\langle u_k(w)u_{-k}(-w)\rangle$ | $\sigma(k, w)$ |
| -------- | -------- | -------- |
| Equilibrium at $\gamma$ and $T$ | $\dfrac{2\gamma T}{(w^2 - K\omega_k^2)^2 + \gamma^2w^2}$ | $0$ |
| Non equilibrium $\gamma_0, \gamma_1, T_0, T_1$ | $\dfrac{\gamma_0T_0 + \gamma_1T_1}{(w^2-K\omega^2_k)^2 + (\gamma_0+\gamma_1)^2w^2}$ | $0$ or $\dfrac{1}{T_0T_1}\frac{4 \gamma_1 \gamma_0 w^2 (T_1-T_0)^2}{\left(w^2-K \omega_k^2\right)^2+(\gamma_0 + \gamma_1)^2w^2 }$ |

There is not any entropon in this case. This is expected from the equilibrium-like nature of the equations.

A follow up would be to analyze the presence of entropons when $\tau$ is finite. In this case we would take $\langle u_k(w)u_{-k}(-w)\rangle^{eq} = \dfrac{\gamma_0T_0 + \gamma_1T_1}{(w^2-K\omega^2_k)^2 + (\gamma_0+\gamma_1)^2w^2}$. We would still have an entropy production proportional to $(T_{0}-T_{1})^2/T_1T_2$ which proves problematic if one want to relate this entropy production to $\langle u_k(w)u_{-k}(-w)\rangle^{neq}$.

________
Lets do this then :):

$$\sigma(w, k)=\dfrac{1}{T_0T_1}\frac{4 \gamma_1 \gamma_0 w^2 (T_1-T_0)^2}{m^2 \left(\tau^2 w^2+1\right) \left(w^2-K \omega_k^2\right)^2-2 \gamma_1 m \tau w^2 \left(w^2-K \omega_k^2\right)+w^2 \left(\gamma_1^2+2 \gamma_1 \gamma_0+\gamma_0^2 \left(\tau^2 w^2+1\right)\right)}$$

$\langle u_k(w)u_{-k}(-w)\rangle=$

$\langle u_k(w)u_{-k}(-w)\rangle^{neq}=\langle u_k(w)u_{-k}(-w)\rangle-\langle u_k(w)u_{-k}(-w)\rangle^{eq}=$
$=$

WWell this does not decompose well :)

For what its worse:

$\langle u_k(w)u_{-k}(-w)\rangle^{neq}/\sigma(k,w) =$ 

________
In any case the fact that some entropy (heat exchange between reservoir) is not taken into account in the dynamic displceentcorrelation function already kills the idea of entropy production.

____________
________
### [29/01/24]

This is false because with the separation above:
$\langle u_k(w)u_{-k}(-w)\rangle^{neq}=\langle u_k(w)u_{-k}(-w)\rangle-\langle u_k(w)u_{-k}(-w)\rangle^{eq}$, $\langle u_k(w)u_{-k}(-w)\rangle^{neq}$ is made to vanish with $\tau\rightarrow 0$, but not $T_0\rightarrow T_1$.
From $\langle u_k(w)u_{-k}(-w)\rangle=$
, it is hard (impossible?) to find a meaningful Eq and Neq contribution...
The Eq contribution is also hard to rationalize intellectually, should it be the SDF (structure displacement factor) of an equilibrium like system with an effective temperature? Or should it be simply the sum of two equilibrium SDF at $T_0$ and $T_1$. The fact that the entropy production is ambiguous makes the task even harder :(.

To recapitulate:

| Lattice model | $\langle u_k(w)u_{-k}(-w)\rangle$|$\langle u_ku_{-k}\rangle$ | $\sigma(k, w)$ |
| -------- | -------- | -------- |-------- |
| Equilibrium at $\gamma$ and $T$ | $\dfrac{2\gamma T}{(w^2 - K\omega_k^2)^2 + \gamma^2w^2}$ | $\dfrac{2}{K\omega_k^2} T$ | $0$ |
| Non equilibrium $\gamma_0, \gamma_1, T_0, T_1$ but equilibrium like | $\dfrac{2(\gamma_0T_0 + \gamma_1T_1)}{(w^2-K\omega^2_k)^2 + (\gamma_0+\gamma_1)^2w^2}$ | $\dfrac{2}{K\omega_k^2} T_k = \dfrac{2}{K\omega_k^2}\dfrac{\gamma_0T_0 + \gamma_1T_1}{\gamma_0 + \gamma_1}$ | $0 ~\textrm{or}~ \dfrac{1}{T_0T_1}\frac{4 \gamma_1 \gamma_0 w^2 (T_1-T_0)^2}{\left(w^2-K \omega_k^2\right)^2+(\gamma_0 + \gamma_1)^2w^2 }$ |
| Non equilibrium $\gamma_0, \gamma_1, T_0, T_1$ and not equilibrium like|$\dfrac{2(\gamma_1T_1+\gamma_0T_0(1+\tau^2w^2))}{\gamma_1^2w^2+2\gamma_1w^2(\gamma_0-\tau w^2+K\tau\omega_k^2)+(1+\tau^2w^2)(\gamma_0^2w^2 + (w^2 -K\omega_k^2)^2)}$| $\dfrac{2}{K\omega_k^2}\dfrac{\gamma_0T_0+\gamma_1T_1+\tau\gamma_0(\tau K \omega_k^2 T_0+\gamma_0T_0+\gamma_1T_1)}{\gamma_0 + \gamma_1+\tau\gamma_0(\tau K\omega_k^2 + \gamma_0 + \gamma_1)}$ |$0 ~\textrm{or}~\dfrac{1}{T_0T_1}\frac{4 \gamma_1 \gamma_0 w^2 (T_1-T_0)^2}{m^2 \left(\tau^2 w^2+1\right) \left(w^2-K \omega_k^2\right)^2-2 \gamma_1 m \tau w^2 \left(w^2-K \omega_k^2\right)+w^2 \left(\gamma_1^2+2 \gamma_1 \gamma_0+\gamma_0^2 \left(\tau^2 w^2+1\right)\right)}$|

It might look like we can define an effective temperature with $\tau\neq 0$ however the velocity velocity structure factor does not have the same temperature. Indeed,
For an equilibrium system equipartition dictates that: $\dfrac{K\omega_k^2\langle u_ku_{-k}\rangle}{2} = \dfrac{m\langle \dot u_k\dot u_{-k}\rangle}{2}=T$. However, for the case $\tau\neq 0$, we have:

$\dfrac{\omega_k^2K\langle u_ku_{-k}\rangle}{2}=\dfrac{\gamma_0T_0+\gamma_1T_1+\tau\gamma_0(\tau K \omega_k^2 T_0+\gamma_0T_0+\gamma_1T_1)}{\gamma_0 + \gamma_1+\tau\gamma_0(\tau K\omega_k^2 + \gamma_0 + \gamma_1)}$
and
$\dfrac{m\langle \dot u_k\dot u_{-k}\rangle}{2} = \dfrac{\gamma_0T_0 + \gamma_1T_1+\tau\gamma_0T_0(\tau K \omega_k^2+\gamma_0+\gamma_1)}{\gamma_0+\gamma_1+\tau\gamma_0(\tau K \omega_k^2 + \gamma_1 + \gamma_0)}$

These two are different, thus the system cannot be seen as an equilibrium system. At least not in the sense of a system behaving equilibrium like with an effective temperature. This is already visible in puglisi's article where the autocorrelation does not match the linear response.

_______________

### Separate Non equilibrium terms with equilibrium ones

We saw that it might be difficult to separate $C^{eq}_{uu}(k, w)$ into a eq and neq contribution based on entropy production. Let's even try to separate $C_{uu}$ into eq and neq without using the entropy production.

________________

Let's start with equilibrium models:

\begin{equation}
\left(-mw^2 - iw \gamma + K\omega_k^2\right) \tilde{u}_k = \sqrt{2\gamma T}\tilde{\xi}_k
\end{equation}

| Without memory equilibrium | $\langle u_k(w)u_{-k}(-w)\rangle$ | $\dfrac{K\omega_k^2}{2}\langle u_ku_{-k}\rangle$ |
| -------- | -------- | -------- |
| | $\dfrac{2\gamma T}{(w^2 - K\omega_k^2)^2 + \gamma^2w^2}$ | $T$ |

___

\begin{equation}
m\ddot{\tilde{{u}}}_{k}(t) = - K\omega_k^2\tilde{{u}}_{k}- \int_{-\infty}^{t} \mathrm d t' \Gamma_k(t -t')\dot{\tilde{u}}_{k}(t') + \tilde{{\mathcal{E}}}_k(t),
\end{equation}
with:
\begin{equation}
\Gamma(t) = \gamma/\tau e^{-t/\tau}
\end{equation}
and:
\begin{equation}
\langle\tilde{\mathcal{E}}_k^{\alpha}(t)\tilde{\mathcal{E}}_q^{\beta}(0)\rangle = \delta_{q, -k}\delta^{\alpha, \beta}T\gamma/\tau e^{-t/\tau}.
\end{equation}
| With memory equilibrium | $\langle u_k(w)u_{-k}(-w)\rangle$ | $\dfrac{K\omega_k^2}{2}\langle u_ku_{-k}\rangle$ |
| -------- | -------- | -------- |
| | $\dfrac{2\gamma T}{-2\gamma\tau w^2(w^2-K\omega_k^2) +(1+\tau^2w^2)(w^2-K\omega_k^2)^2 + \gamma^2w^2}$ | $T$ |

The limit $\tau\rightarrow 0$ corresponds to the equilibrium limit without memory, of course!
______
Note that the shape of the dynamic SDF (wioth $w$ and $k$) depends on the fine details of the driving while the integrated version does not because of equipartition of energy.
________________
------------------
Non equilibrium!
Let's start with:

\begin{equation}
\left(-mw^2 - iw (\gamma_0 + \gamma _{1})+ K\omega_k^2\right) \tilde{u}_k = \sqrt{2\gamma_{0} T_{0}}\tilde{\xi}_k+ \sqrt{2\gamma_1T_1} \tilde {\eta}_{k},
\end{equation}

| Without memory neq | $\langle u_k(w)u_{-k}(-w)\rangle$ | $\dfrac{K\omega_k^2}{2}\langle u_ku_{-k}\rangle$ |
| -------- | -------- | -------- |
| Full | $\dfrac{2(\gamma_0 T_0 +\gamma_1T_1)}{(w^2 - K\omega_k^2)^2 + (\gamma_0+\gamma_1)^2w^2}$ | $\dfrac{\gamma_0T_0 + \gamma_1T_1}{\gamma_0+\gamma_1}$ |
| FDT decomposition 1 | $\dfrac{2\gamma_0 T_0}{(w^2 - K\omega_k^2)^2 + \gamma_0^2w^2} +$ $\dfrac{2\gamma_1 T_1}{(w^2 - K\omega_k^2)^2 + \gamma_1^2w^2}+$awful | $T_0 + T_1 - \dfrac{\gamma_0T_1+\gamma_1T_0}{\gamma_0+\gamma_1}$ |
| FDT decomposition 2 | $\dfrac{\gamma_0 T_0}{(w^2 - K\omega_k^2)^2 + \gamma_0^2w^2} +$ $\dfrac{\gamma_1 T_1}{(w^2 - K\omega_k^2)^2 + \gamma_1^2w^2}+$awful | $\dfrac{T_0}{2} + \dfrac{T_1}{2} + \dfrac{\gamma_0-\gamma_1}{\gamma_0+\gamma_1}\dfrac{T_0-T_1}{2}$ |

And of course equivalently, with $\gamma_{com}$ and $\gamma\omega_k^2$

\begin{equation}
\left(-mw^2 - iw (\gamma + \gamma _{com}\omega_k^2)+ K\omega_k^2\right) \tilde{u}_k = \sqrt{2\gamma_{com} T_{com}}\tilde{\xi}_k+ \sqrt{2\gamma \omega_k^2T} \tilde {\eta}_{k},
\end{equation}

| Without memory neq | $\langle u_k(w)u_{-k}(-w)\rangle$ | $\dfrac{K\omega_k^2}{2}\langle u_ku_{-k}\rangle$ |
| -------- | -------- | -------- |
| Full | $\dfrac{2(\gamma\omega_k^2 T +\gamma_{com}T_{com})}{(w^2 - K\omega_k^2)^2 + (\gamma\omega_k^2+\gamma_{com})^2w^2}$ | $\dfrac{\gamma\omega_k^2T + \gamma_{com}T_{com}}{\gamma\omega_k^2+\gamma_{com}}$ |
| FDT decomposition 1 | $\dfrac{2\gamma\omega_k^2 T}{(w^2 - K\omega_k^2)^2 + (\gamma\omega_k^2)^2w^2} +$ $\dfrac{2\gamma_{com} T_{com}}{(w^2 - K\omega_k^2)^2 + \gamma_{com}^2w^2}+$awful | $T + T_{com} - \dfrac{\gamma\omega_k^2T_{com}+\gamma_{com}T}{\gamma\omega_k^2+\gamma_{com}}$ |
| FDT decomposition 2 | $\dfrac{\gamma\omega_k^2 T}{(w^2 - K\omega_k^2)^2 + (\gamma\omega_k^2)^2w^2} +$ $\dfrac{\gamma_{com} T_{com}}{(w^2 - K\omega_k^2)^2 + \gamma_{com}^2w^2}+$awful | $\dfrac{T}{2} + \dfrac{T_{com}}{2} + \dfrac{\gamma\omega_k^2-\gamma_{com}}{\gamma\omega_k^2+\gamma_{com}}\dfrac{T-T_{com}}{2}$ |

________
Finaly, strong neq (oh god!):

$$m\ddot{\mathbf{u}}_{k}(t) = - K\omega_k^2\mathbf{u}_{k}- \int_{-\infty}^{t} \mathrm{d} t' \Gamma_k(t -t')\dot{\mathbf u}_{k}(t') + \mathbf{\mathcal{E}}_k(t)$$
with:
$$ \Gamma(t) = 2\gamma_0\delta(t) + \dfrac{\gamma_1}{\tau}e^{-t/\tau}$$
and $$\langle\mathcal{E}_k^{\alpha}(t)\mathcal{E}_q^{\beta}(0)\rangle = 2\delta_{q, -k}\delta^{\alpha, \beta}\left(\gamma_0T_0\delta(t)+\dfrac{\gamma_1}{2\tau}T_1e^{-|t|/\tau}\right)$$

| Without half neq | $\langle u_k(w)u_{-k}(-w)\rangle$ | $\dfrac{K\omega_k^2}{2}\langle u_ku_{-k}\rangle$ |
| -------- | -------- | -------- |
| Full | $\dfrac{2(\gamma_1T_1+\gamma_0T_0(1+\tau^2w^2))}{\gamma_1^2w^2+2\gamma_1w^2(\gamma_0-\tau w^2+K\tau\omega_k^2)+(1+\tau^2w^2)(\gamma_0^2w^2 + (w^2 -K\omega_k^2)^2)}$| $\dfrac{2}{K\omega_k^2}\dfrac{\gamma_0T_0+\gamma_1T_1+\tau\gamma_0(\tau K \omega_k^2 T_0+\gamma_0T_0+\gamma_1T_1)}{\gamma_0 + \gamma_1+\tau\gamma_0(\tau K\omega_k^2 + \gamma_0 + \gamma_1)}$ |
| FDT decomposition 1 | $\dfrac{2\gamma_0 T_0}{(w^2 - K\omega_k^2)^2 + \gamma_0^2w^2} +$ $\dfrac{2\gamma T}{-2\gamma\tau w^2(w^2-K\omega_k^2) +(1+\tau^2w^2)(w^2-K\omega_k^2)^2 + \gamma^2w^2}+$awful | $T_0 + T_1$ + awful |
| FDT decomposition 2 | $\dfrac{\gamma_0 T_0}{(w^2 - K\omega_k^2)^2 + \gamma_0^2w^2} +$ $\dfrac{2\gamma T}{-2\gamma\tau w^2(w^2-K\omega_k^2) +(1+\tau^2w^2)(w^2-K\omega_k^2)^2 + \gamma^2w^2}+$awful | $\dfrac{T_0+T_1}{2}+\dfrac{\gamma_0-\gamma_1}{\gamma_0+\gamma_1}\dfrac{T_0-T_1}{2}+\tau^2\dfrac{\gamma_0\gamma_1K(T_0 - T_1)\omega_k^2}{(\gamma_0+\gamma_1)(\gamma_0 + \gamma_1 +\gamma_0\gamma_1\tau + \gamma_0\tau(\gamma_0+K\tau\omega_k^2))}$ |

## [31/01/2024] Generalized entropons

As we saw, non equilibrium excitations with multiple baths are hardly related to the entropy production. There is a "general" (a genralization of the Harada-Sasa relation) formula that takes into account these collective excitations that are clearly not equilibrium like but still cannot be realted to entropy production. Almost everything is in Puglisi's article" about cross correlation.

We start again with:

$$\dot{\mathbf{v}}_{k}(t) = \overbrace{- K\omega_k^2\mathbf{u}_{k} + f(\mathbf u_k,\mathbf k) -\int_{-\infty}^{\infty} \mathrm{d} t' \Gamma_k(t -t')\mathbf v_{k}(t')}^{D} + \mathbf{\mathcal{E}}_k(t)$$

where $f$ is some force and $\mathbf v = \dot {\mathbf u}$. For simplicity we drop the $k$ dependence, the vectorial dependence, $\langle\mathcal{E}(0)\mathcal{E}(t)\rangle=\nu(t)$ and $\Gamma(t)=\gamma(t)\Theta(t)$.

The probability of a trajectory is:

$$P[u(t)]\propto e^{-\displaystyle{\dfrac{1}{2}\int_{-\infty}^{\infty}\int_{-\infty}^{\infty}dtdt'[\dot v(t) - D(u(t))]\nu^{-1}(t - t')[\dot v(t')) - D(u(t'))]}}$$

Equivalently:

$$P[u(w)]\propto e^{-\displaystyle{\dfrac{1}{2}\int_{-\infty}^{\infty}dw[-iw v(w) - D(u(w))]\nu^{-1}(w)[iwv(-w) - D(u(-w))]}}$$

- $\nu(w)=\nu(-w)$ because of time invariance of the noise correlation in real space.
- $D(w)=\overbrace{-K\omega_k^2u(w) + f(w)}^{F} - \Gamma(w)v(w)$
- $\Gamma(w) = \int_{-\infty}^{\infty}\cos(wt)\Gamma(t)dt + iw\int_{-\infty}^{\infty}\dfrac{\sin(wt)}{w}\Gamma(t)dt=\phi(w) + iw\psi(w)$. Conveniently, $\phi$ and $\psi$ are even in $w$.

$$P[u(w)]\propto e^{-\displaystyle{\dfrac{1}{2}\int_{-\infty}^{\infty}dw[-iw v(w)(1-\psi(w))+\phi(w) v(w) - F(w)][iwv(-w)(1-\psi(w))+v(-w)\phi(-w) - F(-w)]\nu^{-1}(w)}}$$

- If $\mathcal{T}$ is the time reversal operator: $\mathcal{T}u(w) = u(w)$ while $\mathcal{T}v(w) = -v(w)$
$\tilde V = -iw v(w)(1-\psi(w))\Rightarrow \mathcal{T}{\tilde V}(w)=\tilde V(w)$
$\mathcal{T} F(w) = F(w)$
$\tilde v(w) = v(w)\phi(w) \Rightarrow \mathcal{T} \tilde v(w) = -\tilde v(w)$

$$\boxed{\langle\sigma(w)\rangle=\left\langle\log\dfrac{P[u(w)]}{P[\mathcal{T}u(w)]}\right\rangle=\nu^{-1}(w)\langle F(w)v(-w)+F(-w)v(w) \rangle\phi}$$

Which, [I think, is in accord with Cugliandolo's result.](https://sci-hub.se/10.1088/1742-5468/2005/09/p09013) (h is a force and $T^{-1}$ is a convolution)
The issue is that Puglisi says that his result (which is different from mine, is also in accord with Cugliandolo's :'(.)

IT might be easier nicer to define the entropy production with the imaginary part of the correlation

\begin{equation}
\boxed{
\langle\sigma(w)\rangle= -2\mathcal{Im}(\langle F(w)u(-w) \rangle)\phi\nu^{-1}w=2\mathcal{Im}(\langle F(-w)u(w) \rangle)\phi\nu^{-1}w}
\end{equation}

Now, let's compute the linear response function by adding a small constant force $h$:

$R(w)\equiv \dfrac{\delta \langle u(w)\rangle}{\delta h(w)}=\displaystyle{\int \mathcal{D}[x]\dfrac{\delta P[u[w]]}{\delta h(w)}u(w)}=\displaystyle{\int \mathcal{D}[x]P[u[w]]\nu^{-1}\left(iwu(w)v(-w)-u(w)F(-w) + u(w)v(-w)(\phi(-w) - iw \psi(-w))\right)}$

We finaly find:
\begin{equation}
\begin{split}
R(w)\nu &=iw\langle u(w)v(-w)\rangle(1-\psi(-w)) - \langle u(w)F(-w)\rangle + \langle u(w)v(-w)\rangle\phi(-w)\\
&=-w^2(1-\psi(-w))\langle u(w)u(-w)\rangle-\langle u(w)F(-w)\rangle +iw \langle u(w)u(-w)\rangle\phi(-w)\\
R(w)\nu&=-w^2(1-\psi(-w))C(w)-\langle u(w)F(-w)\rangle +iw \phi(-w) C(w)
\end{split}
\end{equation}

where we defined the correlation function $C(w) = \langle u(w)u(-w)\rangle$.

Notably, we have:
$$\mathcal{Im}(R(w))\nu=-\mathcal{Im}(\langle u(w)F(-w)\rangle) +w \phi(w) C(w)$$

Which gives us the generalized equation for the entropons:
Notably, we have:

$$ \dfrac{C(w)}{\nu}=\dfrac{\mathcal{Im}(R(w))}{w\mathcal{Re}(\Gamma)}+\dfrac{\sigma(w)}{2\phi^2w^2} $$

($\mathcal{Re}(\Gamma) = \phi$)

Now that it is written, it is trivial :). This is simply the FDT + the entropy production. Equivalently the Harada Sasa equation.

If we reintroduce the $k$, we finaly obtain what we were looking for:

$$\boxed{\dfrac{C_{uu}(\mathbf k, w)}{\nu(\mathbf k, w)}=\dfrac{\mathcal{Im}(R(\mathbf k,w))}{w\mathcal{Re}(\Gamma(\mathbf k,w))}+\dfrac{\sigma(\mathbf k, w)}{2\mathcal{Re}(\Gamma(\mathbf k,w))^2w^2}}$$

This also tells us how should $C_{uu}$ behaves when different modes at different temperatures are coupled.

This is of course equivalent to Caprini's equation in the case of delta correlated noise and a FDT respecting damping (i think Caprini's lacking a factor 2 in his main formula, but it is here everywhere else):

$$\dfrac{C_{uu}^{caprini}}{T}=\dfrac{\mathcal{Im}(R)}{w}+\dfrac{\sigma}{\gamma w^2}$$

is indeed (almost) equal to:

$$\dfrac{C_{uu}^{total}}{2\gamma T}=\dfrac{\mathcal{Im}(R)}{\gamma w}+\dfrac{\sigma}{2\gamma^2 w^2}$$

________
(using Cugliandolo's article) Let's check the total relation in the case of non additional forces or equivalently 0 entropy production (linear forces don't contribute to the Entropy production):

$-mw^2u(w)-iw\Gamma(w)u(w)+Ku(w)=\eta(w)$ and $G(w)=\dfrac{1}{-mw^2-iw\Gamma(w)+K}$

This is the usual relation I use to compute $C_{uu}$ because $C_{uu}=|G(w)|^2\nu(w)$.

for later convenience, the imaginary part of $G(w)$ is simply $\mathcal{Im}(G)=w\mathcal{Re}(\Gamma(w))|G(w)|^2$.

Good, now: the linear response is defined as $R(t-t')=\left.\dfrac{\delta\langle u(t)\rangle}{\delta h(t')}\right|_{h = 0}=\chi(t'-t)$ with $u(t)=\int_{-\infty}^{\infty}dt' \chi(t' - t)(\eta(t')+h(t')))$. Thus for linear systems: $R(t) \equiv \chi(t)=G(t)$. Notably, $\mathcal{Im}(G)=\mathcal{Im}(R)$.

Thus: $C_{uu}=\dfrac{\nu\mathcal{Im}(R)}{w\mathcal{Re}(\Gamma)}$ which is the above relation for the entropons (without entropons).

### Entropons with Auxiliary variables

Now, I'd like to compute $C_{uu}$ but with auxiliary variables:

$$\dot X= A X + B \eta + F$$

With $X_0(t) = u(t), X_{\{1, n-1\}}(t) = \Omega_{\{2, n\}}(t)$ and $X_n(t) = \dot u(t)$.

$A =\begin{bmatrix}
-\gamma & 1 & 1 & \dots &1 &-K \\
-\frac{\gamma_1}{\tau_1} & -1/\tau_1 & 0 & 0 &0& 0 \\
-\frac{\gamma_2}{\tau_2} & 0 & -1/\tau_2 & 0 & 0&0 \\
\vdots & 0 & 0 & \ddots &0 &0 \\
-\frac{\gamma_{n-1}}{\tau_{n-1}} & 0 & 0 & 0 &-1/\tau_{n-1} &0 \\
1 & 0 & 0 & 0 &0 &0
\end{bmatrix}, B = \textrm{diag}(\sqrt{2\gamma T}, \sqrt{2\gamma_1T_1/\tau_1^2}, \dots, \sqrt{2\gamma_{n-1}T_{n-1}/\tau_{n-1}^2}, 0) \textrm{ and } F_0=F, F_{j\neq0}=0$

This is a multivariate ornstein uhlenbeck process.

And I'm too lazy to typeset everything :^). In any case, what we will find is that the entropy production of this model has more terms than the generalized langevin equation even for an equivalent system. What this shows is that entropons are really really generated by the reduced entropy production (as Andrea called it) and not the total one.

For example here is respectively the entropy production with auxiliary variable and genberalized langevin equyation for two bath with memory and external force.
We see that the entropy production due to VF/T is the same for both model. The auxirliary variable model also include a energy transfer term (in red). But, both model have a "heat transfer-Force" dependence (blue) that is different.

________________________________________________________

We note that in the formula above, we took $F$ to be reversible: $\mathcal{T}F = F$. With $F = Ku+f$, we know that the linear force does not contribute to the entropy production. Let's focus only only on $f$. What happens if we take a general form for $f=f_r + f_i$ such that $\mathcal{T}f = f_r - f_i$.

In this case the entropy production reads:

$$
\begin{split}\langle\sigma(w)\rangle\nu(w)=&\langle f_r(w)v(-w)+f_r(-w)v(w) \rangle\phi\\
&- \langle f_i(w)f_r(-w) + f_i(-w)f_r(w)\rangle\\
&-i\langle f_i(w)v(-w) - f_i(-w)v(w)\rangle w(\psi(w)-1)
\end{split}$$

Or equivalently:

$$
\begin{split}\dfrac{\langle\sigma(w)\rangle\nu(w)}{2}=&\mathcal{Im}\langle f_r(-w)u(w)\rangle w\phi\\
&+ \mathcal{Re}\langle f_i(w)f_r(-w)\rangle\\
&+\mathcal{Re}\langle f_i(w)u(-w)\rangle w^2(\psi(w)-1)
\end{split}$$

and the response to correlation relation reads the same as before:

$R(w)\nu=w^2(1-\psi(w))C(w)-\langle u(w)(Ku(-w) + f_r(-w) + f_i(-w))
\rangle +iw \phi(w) C(w)$

We see that there is no way to make it work in the presence of reversible and irrerversible forces. With reversible forces alone. We can simply take the real part of $R(w)$ to make it work:

$$\mathcal{Re}(R(w))\nu = -C(w)(w^2(1-\psi(w))-K)-\mathcal{Re}\langle u(w)f_i(-w)\rangle$$

Such that $$C(w)\nu^{-1}=\dfrac{-\mathcal{Re}(R(w))}{w(w-\mathcal{Im}(\Gamma(w)))-K}-\dfrac{\sigma(w)}{2w(\mathcal{Im}(\Gamma(w))-w)(K + w(\mathcal{Im}(\Gamma(w))-w))}$$

Note that the first term has the form of a fluctuation dissipation theorem. This probably comes from a Kramers-Kronig like relation (note for later me, [this is a good source)](https://lampx.tugraz.at/~hadley/ss2/lectures17/apr7.pdf)

Note also that it works :^):


If we have both forces, but they are uncorrelated $\langle f_if_r\rangle$ we can get closer to the entropons formula but we will still have contributions like $\mathcal{Im}\langle f_i u\rangle$ and $\mathcal{Re}\langle f_r u\rangle$ to the correlation function that is not in the entropy production.

More generally, it tells us that the Harada Sasa formula should not work with a mix of reversible and irreversible forces.

____________
### Can we divide the FDT effective temperature into a non equilibrium contribution and an equilibrium one?

Without entropy production (due to active forces), we have:
$$C_{uu}(\mathbf k, w)=\dfrac{\nu(\mathbf k, w)}{\mathcal{Re}(\Gamma(\mathbf k,w))}\dfrac{\mathcal{Im}(R(\mathbf k,w))}{w}$$

The first term: $\dfrac{\nu(\mathbf k, w)}{\mathcal{Re}(\Gamma(\mathbf k,w))}$ is the effective temperature $2\tilde T_{\textrm{eff}}(k, w)$ of the system.

In the general case, for the system described above, we cannot define an entropy production due to heat flow because we do not have two (or more) well defined thermostats that can exchange heat. But in the case of multiple reservoirs effectively creating an effective temperature, it is interesting to ask if the effective temperature can be in part linked to the entroppy production due to heat flows.

For example, two delta correlated bath at $\gamma_0, \gamma_1, T_{0}$ and $T_1$ effectively act as a single global bath at $T_{eff}=\dfrac{\gamma_0T_0+\gamma_1T_1}{\gamma_0+\gamma_1}$. The heat flows are $\dot{Q}_0=-\dot{Q}_1=\dfrac{\gamma_0+\gamma_1}{\gamma_0\gamma_1}(T_0-T_1)$ and thus the entropy production rate is $\dot s=\dot Q_0/T_0 + \dot Q_1/T_1=\dfrac{\gamma_0+\gamma_1}{\gamma_0\gamma_1}\dfrac{(T_0-T_1)^2}{T_0T_1}$. The effective temperature can be decomposed as:

$$T_{eff}=T_0 + \dfrac{T_1}{2} + \dfrac{\gamma_0-\gamma_1}{\gamma_0+\gamma_1}\dfrac{T_0-T_1}{2}=\dfrac{T_0}{2} + \dfrac{T_1}{2}+\dfrac{\gamma_0\gamma_1(\gamma_0-\gamma_1)}{(\gamma_0+\gamma_1)^2}\dfrac{\dot Q}{2}$$

Unatural and in any case, this is not a good decomposition because the neq contribution should not vanish in the limit $\gamma_0\rightarrow\gamma_1$ (it could... but let's say no!). An other decomposition would be:

$$T_{eff}=T_0 + \dfrac{\gamma_1}{\gamma_0+\gamma_1}\dfrac{T_1-T_0}{2}=T_0 + \dfrac{\gamma_0\gamma_1^2}{(\gamma_0+\gamma_1)^2}\dfrac{\dot Q}{2}$$

no really better :^.

___________
Let's see what happens in $w$ space:

$$\sigma(w, k)=\dfrac{1}{T_0T_1}\frac{4 \gamma_1 \gamma_0 w^2 (T_1-T_0)^2}{\left(w^2-K \omega_k^2\right)^2+w^2( \gamma_0 + \gamma_1)^2 }$$

$$\dot Q(w, k)=\frac{4 \gamma_1 \gamma_0 w^2 (T_1-T_0)}{\left(w^2-K \omega_k^2\right)^2+w^2( \gamma_0 + \gamma_1)^2 }$$

But rememeber, $T_{eff}$ is still the same in $w$ or $k$ space or whatever. IT IS T_EFF (it can depends on $k$ if the baths are coupled differently to $k$ wave vector or on $w$ of the baths have correlation):

$$T_{eff}=\dfrac{\gamma_0T_0+\gamma_1T_1}{\gamma_0+\gamma_1}$$

We see that there is no way to obtain $\dot Q$ proportional to $T_{eff}$. Maybe, we can include it in a smarter way:

Here, $C_{uu}(k, w)=2T_{eff}\mathcal{Im}(R) = 2(T_0+\Delta T^{neq})\mathcal{Im}(R)=2T\mathcal{Im}(R)+2\Delta T\mathcal{Im(R)}$

We have $$2\Delta T\mathcal{Im(R)}=2\dfrac{\gamma_1}{\gamma_1 + \gamma_0}(T_0-T_1)\dfrac{(\gamma_0+\gamma_1)w^2}{\left(w^2-K \omega_k^2\right)^2+w^2( \gamma_0 + \gamma_1)^2 }=\frac{\dot Q}{2\gamma_0w^2}$$

__________

Delta + Exponential:

$C_{uu}=\dfrac{2(\gamma_1T_1+\gamma_0T_0(1+\tau^2w^2))}{\gamma_1^2w^2+2\gamma_1w^2(\gamma_0-\tau w^2+K\tau\omega_k^2)+(1+\tau^2w^2)(\gamma_0^2w^2 + (w^2 -K\omega_k^2)^2)}$

Heat Transfer:


It will also work.

_________

Let's trind to find the relation analytically.

$\Gamma(t)=(2\gamma\delta(t)+\sum_{m=1}^{n}\dfrac{\gamma_m}{\tau_m}e^{-t/\tau_m})\Theta(t)$
$\nu(t)=2\gamma T\delta(t)+\sum_{m=1}^{n}T_m\dfrac{\gamma_m}{\tau_m}e^{-t/\tau_m}$

this is euivalent to:
$$\dot X= A X + B \eta + F$$

With $X_0(t) = u(t), X_{\{1, n-1\}}(t) = \Omega_{\{2, n\}}(t)$ and $X_n(t) = \dot u(t)$.

$A =\begin{bmatrix}
-\gamma & 1 & 1 & \dots &1 &-K \\
-\frac{\gamma_1}{\tau_1} & -1/\tau_1 & 0 & 0 &0& 0 \\
-\frac{\gamma_2}{\tau_2} & 0 & -1/\tau_2 & 0 & 0&0 \\
\vdots & 0 & 0 & \ddots &0 &0 \\
-\frac{\gamma_{n-1}}{\tau_{n-1}} & 0 & 0 & 0 &-1/\tau_{n-1} &0 \\
1 & 0 & 0 & 0 &0 &0
\end{bmatrix}, B = \textrm{diag}(\sqrt{2\gamma T}, \sqrt{2\gamma_1T_1/\tau_1^2}, \dots, \sqrt{2\gamma_{n-1}T_{n-1}/\tau_{n-1}^2}, 0) \textrm{ and } F_0=F, F_{j\neq0}=0$

This is a multivariate ornstein uhlenbeck process with Onsger-Machlup functional:

$$\log(P(X)) =-\dfrac{1}{2}\int_{-\infty}^{\infty}dw (-iw X(w) - AX(w) - F(w))^t(BB^t)^{-1}(iw X(-w) - AX(-w) - F(-w))$$

The linear response function is given by
\begin{equation}
\begin{split}
R(w) &= \dfrac{\delta \langle u(w)\rangle}{\delta F(w)}=\dfrac{1}{2\gamma T}\left[iw\langle u(w)v(-w)\rangle-\langle u(w)\left(-\gamma v(-w)-Ku(w)+\sum_{k=1}^{n-1}\Omega_k(-w)+F(-w)\rangle\right)\right]\\
&=\dfrac{1}{2\gamma T}(-w^2C_{uu}+i\gamma wC_{uu} + KC_{uu}-\sum_{k = 1}^{n-1}\langle u(w)\Omega_k(-w)\rangle -\langle u(w) F(-w)\rangle)
\end{split}
\end{equation}

From which we obtain the relation:

$$2\gamma T\mathcal{Im}(R(w))=\gamma w C_{uu}(w)-\sum_{k = 1}^{n-1}\mathcal{Im}\langle u(w)\Omega_k(-w)\rangle -\mathcal{Im}\langle u(w) F(-w)\rangle$$

Note the similarity with the formula written by puglisi

------

The entropy production for the model is:
\begin{split}\sigma^{OU}(w)&=\dfrac{\langle v(w)F(-w)\rangle}{2T}+\frac 1 2 \displaystyle{\sum_{k=1}^{n-1}\left(\dfrac{1}{T}-\dfrac{1}{T_k}\right)}\langle v\Omega_m\rangle+c.c\\
&=\dfrac{\langle v(w)F(-w)\rangle}{2T}+\frac 1 2 \displaystyle{\sum_{k=1}^{n-1}\left(\dfrac{1}{T}-\dfrac{1}{T_k}\right)}\langle v\rangle\langle\Omega_m\rangle+\frac 1 2 \displaystyle{\sum_{k=1}^{n-1}\left(\dfrac{1}{T}-\dfrac{1}{T_k}\right)}\langle\delta v(w) \delta \Omega_k(-w)\rangle+c.c
\end{split}
We define the heat production due to reservoir k as: $\dot{Q}^k=\langle \delta v(w) \delta \Omega_k(-w)\rangle+\langle \delta v(-w) \delta \Omega_k(w)\rangle=2w\mathcal{Im}(\langle \delta u(w) \delta \Omega_k(-w)\rangle)$

And the active entropy production as:
$\sigma^a(w)=\overbrace{\dfrac{\langle v(w)F(-w)\rangle}{2T}}^{\sigma ^{ideal}}+\dfrac 1 2 \displaystyle{\sum_{k=1}^{n-1}\left(\dfrac{1}{T}-\dfrac{1}{T_k}\right)}\overbrace{\langle v\rangle\langle\Omega_k\rangle}^{\dot Q_{a}^k} + c.c$

Such that the total entropy production is:

$$\sigma^{OU} =\sigma^a +\frac 1 2 \displaystyle{\sum_{k=1}^{n-1}\left(\dfrac{1}{T}-\dfrac{1}{T_k}\right)}\dot Q^k $$

I don't like the fact that there is no energy exchange between the $k$ thermostats. In the limit $\tau_k\rightarrow 0$ they are no different from the reservoir at $T$.

Anyway, back at the equation.
In the simple case where the external force is equal to 0. We observe calorons (shout out to Andrea, on top of that, the name would please Murray Gell-Mann because the root and the suffixe are both latin !!! ):

$$\boxed{C_{uu}=\dfrac{2T}{w}\mathcal{Im}(R)+\overbrace{\dfrac{1}{2w^2\gamma}\sum_{k}\dot Q^k}^{\textrm{calorons}}}$$

We recall that the system is equivalent to one with a generalized langevin equation:

$$C_{uu}=\dfrac{2T}{w}\mathcal{Im}(R)+\overbrace{\dfrac{1}{2w^2\gamma}\sum_{k}\dot Q^k}^{\textrm{calorons}}=\dfrac{\nu(w)}{w\mathcal{Re}(\Gamma(w))}\mathcal{Im}(R)$$

_______

When forces are present:

$$C_{uu}=\dfrac{2T}{w}\mathcal{Im(R)}+\dfrac{1}{2w^2\gamma}\sum_k(\dot Q^k + \dot Q_a^k)+ \dfrac{T\sigma^{ideal}}{w^2\gamma}$$

The first embarassing thing is that we have to include $\dot Q_a^k$ as a calorons even though it is more like an active source of entropy production. Indeed, the entropy production associated to this term is: $\tilde\sigma=\sum_k (\dfrac 1 {T_k} - \dfrac 1 T)\dot Q^k_a$ while the terms present in the reponse function are: $\sim \sum_k \dot Q^k_a$. We surely could divide the calorons by the appropriate terms to make the entropy production appear. This would give:

$$C_{uu}=\dfrac{2T}{w}\mathcal{Im(R)}+\dfrac{T}{w^2\gamma}\sum_k\overbrace{\dfrac{T_k}{T - T_k}(\sigma^{Q^k}}^{calorons} + \sigma^{ Q_a^k})+ \dfrac{T\sigma^{ideal}}{w^2\gamma}$$

But this is pretty ugly :) And of top of that thge ratio of temperature diverges at equilibrium (eveyrthing is fine of course because the entropy producxtion goes to 0)

We might want to compare the two equivalent Harada-Sasa relations we found. Where here:

$\Gamma(t)=(2\gamma\delta(t)+\sum_{m=1}^{n}\dfrac{\gamma_m}{\tau_m}e^{-t/\tau_m})\Theta(t)$
$\nu(t)=2\gamma T\delta(t)+\sum_{m=1}^{n}T_m\dfrac{\gamma_m}{\tau_m}e^{-t/\tau_m}$

\begin{split}
\dfrac{C_{uu}(\mathbf k, w)}{\nu(\mathbf k, w)}&=\dfrac{\mathcal{Im}(R(\mathbf k,w))}{w\mathcal{Re}(\Gamma(\mathbf k,w))}+\dfrac{\sigma^{gene}(\mathbf k, w)}{2\mathcal{Re}(\Gamma(\mathbf k,w))^2w^2}\\
\dfrac{C_{uu}(\mathbf k, w)}{2\gamma T}&=\dfrac{\mathcal{Im}(R(\mathbf k,w))}{w\gamma}+\dfrac{1}{4w^2\gamma^2}\sum_k\dfrac{\dot Q^k + \dot Q_a^k}{T}+ \dfrac{\sigma^{ideal}}{2w^2\gamma^2}
\end{split}

I would say that, on a general basis, we cannot assert that:

$$\dfrac{2T}{w}\mathcal{Im}(R)+\overbrace{\dfrac{1}{2w^2\gamma}\sum_{k}\dot Q^k}^{\textrm{calorons}}=\dfrac{\nu(w)}{w\mathcal{Re}(\Gamma(w))}\mathcal{Im}(R)$$

because the entropy production of the generalized langevin equation is not always equal to the entropy production of the equivalent markjovian system without taking into account the "static/heat transfer between reservoir" term.

______________________

### Conceptual issues with the Markovianization.

We already saw that the entropy production is dependent on the level of description and thus sensible to the markovianization or inversely the coarse graining of a systemm even if every static and dynamical measurable quantities are the same.

Let's take a simple system (that Puglisi likes a lot):

$$\dot v = -\int_{-\infty}^{\infty}\Gamma(t'-t)v(t')dt'+F+\eta(t)$$

where

$$\Gamma(t) = \gamma\delta(t)+\gamma_1/\tau_1e^{-t/\tau_1}\Theta(t)$$
$$\nu(t) = 2\gamma T\delta(t)+\gamma_1/\tau_1T_1e^{-|t|/\tau_1}$$

_________

The equivalent markovian system is

\begin{split}
\dot v &= -\gamma v+\Omega + F+ \xi\\
\dot \Omega &=-\gamma_1/\tau_1 v-\dfrac{\Omega }{\tau_1}+\zeta
\end{split}

With $\langle\eta\eta\rangle=2\gamma T$ and $\langle\zeta\zeta \rangle=2\gamma_1T_1/\tau_1^2$

____________

The derivation of the entropy production of the generalized system is fucked up due to the constant force $F$. Let's go to real space:

$\sigma^{generalized}(t)=\int_{-\infty}^{t}dt'\tilde T^{-1}(t-t')(\dot r(t)F + \dot r(t')F)=F\int_{-\infty}^{t}dt'\tilde T^{-1}(t-t')(v(t) + v(t'))$

With $\tilde T^{-1}=\mathcal{TF}(2\mathcal{Re}(\Gamma(w))/\nu(w))$

We find $\tilde T^{-1}(t)=e^{-|t|/\tau}\dfrac{\tau\gamma_1(T-T_1)}{2\gamma T^2\tau_1^2 }+\dfrac{\delta(t)}{T}$, with $\tau^2=\dfrac{\tau_1^2}{1+\gamma_1T_1/\gamma T}$

In the large time limit, $v=\tilde v=F/(\gamma + \gamma_1)$

The entropy production is exactly computed as:

$$\sigma^{generalized}=F\tilde v\left(\dfrac{1}{T_f}+\dfrac{\tau^2\gamma_1(T-T_1)}{2\gamma T^2\tau_1^2 }\right)=F\tilde v/T_{eff}$$ with $T_{eff}=\dfrac{\gamma T + \gamma_1 T_1}{\gamma + \gamma_1}$. Note that this does not depend on $\tau_1$. In any case, in the trivial limit $\tau_1\rightarrow0$ (equivalent to the limit of a single bath at $T_{eff}$, it is interesting to see that we obtain the entropy production $\tilde vF/T_{eff}$

________

Let's do the computatio with the markovian system. In this case, the entropy production reads ($\langle \Omega\rangle=-\gamma_1\tilde v$):

$\sigma^{OU}=\overbrace{\tilde v F\left(\dfrac{1}{T}-\dfrac{\gamma_1}{\gamma_1 + \gamma}\left(\dfrac{1}{T}-\dfrac{1}{T_1}\right)\right)}^{\sigma^a}+\dfrac{\gamma\gamma_1}{(\gamma_1 + \gamma)(1+\tau_1\gamma_1)}\dfrac{(T_1-T)^2}{TT_1}$

This gives: $\sigma^a=\tilde vF/\tilde T_{eff}$. With $T_{eff}=T T_1/T_{eff}$, which is problematic because in the limit $\tau_1\rightarrow 0$ we should find $\tilde T_{eff}=T_{eff}$.

The limit is singular in the original equation which surely causes some issues.... But it is really strange that the effective temperature does not depend on $\tau_1$.

__________

### ENTROPONS AND CALORONS OF ACTIVE ORNSTEIN-UHLENBECK BROWNIAN PARTICLES

\begin{split}
\ddot{\mathbf{u}}_{k}(t) &= - K\omega_k^2\mathbf{u}_{k}- \int_{-\infty}^{t} \mathrm{d} t' \Gamma_k(t -t')\dot{\mathbf u}_{k}(t') + f_k(t)+ \mathbf{\mathcal{E}}_k(t)\\
\dot f_k &= -f_k/\tau_a + v_0\sqrt{2/\tau_a}\eta_k(t)
\end{split}
with:
$$ \Gamma(t) = 2\gamma\delta(t) + \dfrac{\gamma_1}{\tau}e^{-t/\tau}$$
and $$\langle\mathcal{E}_k^{\alpha}(t)\mathcal{E}_q^{\beta}(0)\rangle = 2\delta_{q, -k}\delta^{\alpha, \beta}\left(\gamma T\delta(t)+\dfrac{\gamma_1}{2\tau}T_1e^{-|t|/\tau}\right)$$

this is equivalent to:
$$\dot X_k= A X_k + B \eta $$

With $X_0(t) = v_k(t), X_1(t) = \Omega_k(t)$, $X_2(t)=u_k(t)$ and $X_3(t) = f_a(t)$.

$A =\begin{bmatrix}
-\gamma & 1 &-Kw_k^2 &1 \\
-\frac{\gamma_1}{\tau_1} & -1/\tau_1 & 0 & 0 \\
1 & 0 & 0 & 0 \\
0 & 0 & 0 & -1/\tau_a
\end{bmatrix}, B = \textrm{diag}(\sqrt{2\gamma T}, \sqrt{2\gamma_1T_1/\tau^2},0 , v_0\sqrt{2/\tau_a})$

The entropy production is (assuming that $f_a$ is even under time reversal):
$$2\sigma = \dfrac{\langle f(w)v(-w) \rangle}{T}+\dfrac{T - T_1}{T T_1}\langle \Omega(w)v(-w)\rangle + c.c. $$

We recall that:
$$C_{uu}=\dfrac{2T}{w}\mathcal{Im(R)}+\dfrac{\mathcal{Im}\langle u(-w) \Omega(w) \rangle}{w\gamma}+ \dfrac{\mathcal{Im}\langle u(-w) f(w) \rangle}{w\gamma}$$

______
We compute:

\begin{split}
\dfrac{\mathcal{Im}\langle u(-w) \Omega(w) \rangle}{w}&=-\dfrac{4\gamma_1 \tau_av_o^2}{\gamma_1^2w^2+2\gamma_1w^2(\gamma-\tau w^2+K\tau\omega_k^2)+(1+\tau^2w^2)(\gamma^2w^2 + (w^2 -K\omega_k^2)^2)}\dfrac{1}{1+\tau_a^2w^2}\\
&+\dfrac{4\gamma\gamma_1 (T_1-T)}{\gamma_1^2w^2+2\gamma_1w^2(\gamma-\tau w^2+K\tau\omega_k^2)+(1+\tau^2w^2)(\gamma^2w^2 + (w^2 -K\omega_k^2)^2)}
\end{split}

\begin{split}
\dfrac{\mathcal{Im}\langle u(-w) f(w) \rangle}{w}&=\dfrac{4 \gamma_1\tau_av_o^2}{\gamma_1^2w^2+2\gamma_1w^2(\gamma-\tau w^2+K\tau\omega_k^2)+(1+\tau^2w^2)(\gamma^2w^2 + (w^2 -K\omega_k^2)^2)}\dfrac{1}{1+\tau_a^2w^2}\\
&+\dfrac{4 \gamma(1+\tau^2w^2)\tau_av_o^2}{\gamma_1^2w^2+2\gamma_1w^2(\gamma-\tau w^2+K\tau\omega_k^2)+(1+\tau^2w^2)(\gamma^2w^2 + (w^2 -K\omega_k^2)^2)}\dfrac{1}{1+\tau_a^2w^2}
\end{split}

First, note that $\langle u, w, v, f\rangle=0$. Compared to what we did before, we cannot say that $\mathcal{Im}\langle u(-w) \Omega(w) \rangle$ are calorons between $T$ and $T_1$. Because it does not vanish with $T = T_1$ and is highly depends on active forces. Note however that the entropy production does vanish when $T_1 = T$:
$$2\sigma = \dfrac{\langle f(w)v(-w) \rangle}{T}+\dfrac{T - T_1}{T T_1}\langle \Omega(w)v(-w)\rangle + c.c. $$

However, we still need the term to make the formula about $C_{uu}$ work. This is very interesting and tells us that the entropy production is indeed not the right object. Because, the term $\langle \Omega(w)v(-w)\rangle$ does not contribute to the entropy production when $T = T_1$ but it is still needed in the formula for $C_{uu}$. The term iin $\langle \Omega(w)v(-w)\rangle$ not linked to heat transfer between $T$ and $T_1$ is linked to the heat tranfered from the active particle to the environnement.

________

Let's see $C_{uu}$ explicitly (I corrected the fucked up factor 2. I don't know where the issue is, fuck it):

\begin{split}
C_{uu} &= \dfrac{2T}{w}\mathcal{Im}(R)+\dfrac{2\gamma_1 (T_1-T)}{\gamma_1^2w^2+2\gamma_1w^2(\gamma-\tau w^2+K\tau\omega_k^2)+(1+\tau^2w^2)(\gamma^2w^2 + (w^2 -K\omega_k^2)^2)}\\
&+\dfrac{2 (1+\tau^2w^2)\tau_av_o^2}{\gamma_1^2w^2+2\gamma_1w^2(\gamma-\tau w^2+K\tau\omega_k^2)+(1+\tau^2w^2)(\gamma^2w^2 + (w^2 -K\omega_k^2)^2)}\dfrac{1}{1+\tau_a^2w^2}
\end{split}

The first term is the equilibrium term, the second is the caloron term, the last one is the active term. The response function is given by:

$$\mathcal{Im}(R)/w=\dfrac{\gamma_1+\gamma(1+\tau^2w^2)}{\gamma_1^2w^2+2\gamma_1w^2(\gamma-\tau w^2+K\tau\omega_k^2)+(1+\tau^2w^2)(\gamma^2w^2 + (w^2 -K\omega_k^2)^2)}$$

We can rewrite the two first term of $C_{uu}$ as:
$$2T\mathcal{Im}(R)+blabla=\underbrace{2\dfrac{\gamma_1T_1 + \gamma T ( 1 + \tau^2w^2)}{\gamma_1 + \gamma ( 1 + \tau^2w^2)}}_{T_{eff}}\mathcal{Im}(R)$$

As expected the effective temperature appears!

_______________________________

Now, let's do the same computation, the generalized langevin way.

We recall the equation of motion:

\begin{split}
\ddot{\mathbf{u}}_{k}(t) &= - K\omega_k^2\mathbf{u}_{k}- \int_{-\infty}^{t} \mathrm{d} t' \Gamma_k(t -t')\dot{\mathbf u}_{k}(t') + f_k(t)+ \mathbf{\mathcal{E}}_k(t)\\
\dot f_k &= -f_k/\tau_a + v_0\sqrt{2/\tau_a}\eta_k(t)
\end{split}
with:
$$ \Gamma(t) = 2\gamma\delta(t) + \dfrac{\gamma_1}{\tau}e^{-t/\tau}$$
and $$\langle\mathcal{E}_k^{\alpha}(t)\mathcal{E}_q^{\beta}(0)\rangle = 2\delta_{q, -k}\delta^{\alpha, \beta}\left(\gamma T\delta(t)+\dfrac{\gamma_1}{2\tau}T_1e^{-|t|/\tau}\right)=\nu$$

This equation would give us:

$$C_{uu}=\underbrace{\dfrac{\nu(\mathbf k, w)}{\mathcal{Re}(\Gamma(\mathbf k,w))}}_{T_{eff}}\dfrac{\mathcal{Im}(R(\mathbf k,w))}{w}+...\langle u_k(-w)f_k(w)\rangle$$

So far so good, this is the expression found by Caprini! Which is of course equal to what was found above. Except that above, we didn't have $T_{eff}$ but $T = T_{eff}-\Delta T$ in front of the $Im(R)$ and the term $\Delta T Im(R)$ appeared as an entropy/caloron term.

### What happen if we coarse grained this even more?

What happen if we integrate $f_k$? In this case, we obtain:

\begin{split}
\ddot{\mathbf{u}}_{k}(t) &= - K\omega_k^2\mathbf{u}_{k}- \int_{-\infty}^{t} \mathrm{d} t' \Gamma_k(t -t')\dot{\mathbf u}_{k}(t') + \mathbf{\mathcal{E}}_k(t)
\end{split}
with:
$$ \Gamma(t) = 2\gamma\delta(t) + \dfrac{\gamma_1}{\tau}e^{-t/\tau}$$
and $$\langle\mathcal{E}_k^{\alpha}(t)\mathcal{E}_q^{\beta}(0)\rangle = 2\delta_{q, -k}\delta^{\alpha, \beta}\left(\gamma T\delta(t)+\dfrac{\gamma_1}{2\tau}T_1e^{-|t|/\tau}+\dfrac{v_o^2}{2\tau_a}e^{-|t|/\tau_a}\right)=\nu$$

Off course, the entropy production is 0 :). And $C_{uu}$ reads:

$$C_{uu}=\underbrace{\dfrac{\nu(\mathbf k, w)}{\mathcal{Re}(\Gamma(\mathbf k,w))}}_{T_{eff}}\dfrac{\mathcal{Im}(R(\mathbf k,w))}{w}$$

Everything is comprised into an effective temperature:

\begin{split}T_{eff}&=\dfrac{(1+\tau_a^2w^2)(2\gamma T (1+\tau^2w^2)+2\gamma_1T_1)+v_o^2 (1 + \tau^2w^2)}{(1+\tau_a^2w^2)(\gamma_1 + 2\gamma(1 + \tau^2w^2))}\\
&=\dfrac{2\gamma T (1+\tau^2w^2)+2\gamma_1T_1}{\gamma_1 + 2\gamma(1 + \tau^2w^2)}+\dfrac{v_o^2 (1 + \tau^2w^2)}{(1+\tau_a^2w^2)(\gamma_1 + 2\gamma(1 + \tau^2w^2))}
\end{split}

The first term corresponds to the old effective temperature we found and the second term to the entropy production (up a factor $\tau_a$, fuck, it).

### Questions left to answer up to now:

Up to know, we saw that the entropons, for linear systems can really be seen as an effective temperature. Here are some questions:

* The formula linking $C_{uu}$ and $\mathcal{Im}(R)$ does not work for odd force w/ to time reversal symetries. The issue is deeper and affect the Harada Sasa relation. The relation misses odd forces. We can fix it if the only forces present in the system are odd ones. But, if a mix of odd and even forces are present, we cannot express the entropy production of the system as the departure from the FDT. I think this has never really been said. It is interesting.

* How to coarse grain a system with odd forces? For example with the active ornstein uhlenbeck process, $\dot v = F + u$ and $\dot u= -u + \eta$ with $v$ the velocity and $u$ the random force. If $u$ is even under time reversal, we can integrate $u$ into an exponentially correlated noise and put it inside $\dot v$ then compute the entropy production from it. But, if $u$ is odd under time reversal and we iontegrate it, we kind of lose its parity, messing up the residual entropy production. Which is already the case with the markovianization by the way. What are the consequence of such a thing.

* If we go back to the example of a generalized langevin equation that can then be markovianized. We obtain two different entropy production for these models. Even though, every observables will be the same. Interestingly, the entropy production of the markovian system is not made of the entropy production of the generalized langevin equation plus an additional simple term. Which is really unpleasant. I'd like to have a clear idea about how this kind of thing does not affect the Harada Sasa relation (I have some example but not a clear picture in mind: for example, we cannot divide the entropy of the markovian system as entropons and active entropy production. Everything is mixed..)... **Partial answer**: [this paper](https://arxiv.org/pdf/1809.02391.pdf) and partially [this one](https://arxiv.org/pdf/1904.00475.pdf). Since the entropy production changes with the level of coarse graining and the parity of the variable we keep track off. One must then be careful and evaluate the physical knowledge we have about the system at hand. See also my mathematica notebook about active ornstein uhlenbeck particles.

* WHAT HAPPENS FOR NON LINEAR SYSTEMS????????? Shame on me, but what happens to the linear response function? The equivalence we saw above between effective temperature and entropy production only hold because the linear response function does not care about this linear forces. But in the case of non linear forces, I would expecte the linear response function to change, alongside the entropy production. The Harada-Sasa relation still holds (I think? Notrhing in my derivation assumed linearity of the external forces) but I'd like to see numerically wxhat happens. Unforntunaltely this is a bit too complicated and unobtainable I think. And What about toda's chain, non linear-integrable systems?

____________________

### Graphical idea of the non equilibrium excitations
(notebook: entropons graphics)

Note that, the neq part is taken to be always positive (we chose a positive heat tranzsfer), we also chose parameters that exhibit clear neq effects.
Moreover,

#### Double delta + k dependent one

damping: $\gamma_0 + \gamma w_k^2$
noise: $\gamma_0 T_0+\gamma \omega_k^2 T$

| parameters | $C_{uu}(k, w)$ as a function of $w$ with $k$ varying | $C_{uu}(k, w)$ as a function of $k$ with $w$ varying |
| -------- | -------- | -------- |
| $\gamma_0$ and $\gamma_1$ small |  | | |
|$\gamma_0$ and $\gamma_1$ large|||

As expected, the $w$ dependence does not change whatever $k$ is, but the opposite is not true and the neq part vanishes at small $k$.

#### Delta + exponeitially correlated one

damping: $\gamma_0 + \gamma/\tau e^{-t/\tau}$
noise: $\gamma_0 T_0+\gamma/\tau e^{-t/\tau} T$

| parameters | $C_{uu}(k, w)$ as a function of $w$ with $k$ varying | $C_{uu}(k, w)$ as a function of $k$ with $w$ varying |
| -------- | -------- | -------- |
| $\gamma_0$ and $\gamma_1$ small and $\tau$ large |  | | |
|$\gamma_0$ and $\gamma_1$ large|||

As expected in this case the dependence on $k$ is the same for the neq and Eq part but the frequency part is different. Note however, that integrating over $w$ to obtain the structure factor or the static correlation function gives us $C_{uu}(k, t = 0)$. This function is really different from its equilibrium counter part because it will not be flat.

### Fixing the frequency of collision.

Using the boltzman equation with the molecular chaos assumption does not gives us a correct continuous transition. We can interpret this result undeer different facet:

- It gives the enskog frequency of collision which must be modifies in the presence of damping

- We forgot correlation between velocities which are ultra important close to the transition (even for a mean field theory?!) because we observed very strong backscattering, telling us that we were far away from molecular chaos.

- The temperature is not the most well defined and usable quantity around the transition because particles will have a high energy after colliding and a small one before colliding. Therefore, the "global" temperature mightbe a "too washed out, too coarse grained" quantity to deal with. Therefore, it might be easier to deal with exiting velocity for example! Something we already did in the long paper.

Taking into account these elements, we try to obtain an equation for the exiting bvelocity and then average it over some well defined quantity! We will also use Ran Ni model in 1D, because it it simpler to deal with :).

With Ran Ni model, the particle after a collision raise the energy of the system by $\Delta^2$ and then lose an energy: $v_o^2(1 - e^{-2\gamma\tau_f})$ with $v_0$ the modulus of the exiting velocity and $\tau_f$ the travel time before the next collision: $\tau_f=-\log(1-\gamma l(\phi)/v_0)/\gamma$. In the steady state:
$$\Delta^2 = v_o^2(1 - e^{-2\gamma\tau_f})\Rightarrow v_0 = \dfrac{\Delta^2 + \gamma^2l^2}{2\gamma l}$$

The order parameter is (is this the precollisinal velocity? It works from the definition of $\tau_f$ and $v(\tau_f) = $ precollisional but it does not seem to work from the collision rule... Probably because we are mixing kinetic theory and determinsitic mechanics) $\bar v_0 = v_0-\Delta=\gamma\dfrac{(l(\phi) - l(\phi_c))^2}{2l(\phi)}$ with $l(\phi_c) = \Delta/\gamma$ from which we get $\bar v_0(\phi)\sim \dfrac{\gamma l'(\phi_c)^2}{l(\phi_c)}(\phi - \phi_c)^2$
If we define the temperature from the collision average, it will be $T\sim v_0^2$ and $T$ will have an exponent 4.

Unforunately, the real system does not reach a zero exiting velocity at the transition 

the man free path argument is not exact.

NEXT :)

### POPULATION DYNAMICS FOR THE CONTINUOUS TRANSITION

$$\dfrac{dT}{dt} = (2P_A - P_A^2)\dfrac{\omega(T)}{2}\langle \Delta E\rangle_{coll} -2\gamma T $$
$$\dfrac{dP_A}{dt} = \omega(T)P_A(1-P_A) -we^{-\omega \tau_a}P_A$$

Where $\tau_a$ is the time needed to stop (more or less). Let's take it equal to $a/\gamma$ with $a$ around 5.

We will assume that $\langle\Delta E\rangle_{coll} = \Delta^2$ is independent on $T$. This is akin to Ran ni model! We also assumed that the frequency of collision between dead and alvie was the same than between alive and dead... This is off course false! Let's assume that at first.

The steady state is given by $P_A = 1-e^{-\omega \tau_a}$

Using $\omega(\phi, T) = \omega_0(\phi)\sqrt{T}$ We find the second steady state equation:

$$0=-4\gamma T + \Delta^2(1- e^{-2\sqrt{T}\tau_a\omega_0})\omega_0\sqrt{T}$$

$T=0$ is always a solution. The second solution is given by mathematica to be:

$$T = \left(\dfrac{\Delta^2 \tau_a \omega_0^2+2\gamma W_0(-\Delta^2e^{-\Delta^2\tau_a w_o^2/2\gamma}\tau_a\omega_0^2/2\gamma)}{4\gamma\tau_a\omega_0}\right)^2$$

We define $\Gamma = \Delta^2\tau_a\omega_0^2/2\gamma$

$$T = \left(\dfrac{\Gamma+ W_0(-\Gamma e^{-\Gamma})}{2\tau_a\omega_0}\right)^2$$

Where $W_0(z)$ is the principal lambert function; one of the solution of $z = xe^x\Rightarrow x = W_0(z)$. It is defined when $z>1/e$ and moreover $W_0(1/e) = -1$. Hence, the transition is continuous and happens when $\Gamma = 1$. Equivalently, when $\Delta^2\tau_a/2\gamma=1/\omega_0^2$. From the mean free path argument, we can find $\tau_a$. Indeed, the mean free path argument tells us that the transiton happens at $l(\phi_c) = \Delta/\gamma$. Moreover, $\omega(\phi, T)=\langle |v|\rangle/l(\phi)\Rightarrow\omega_0 = \sqrt{\pi/2}/l(\phi)$. Thus to be consistent with the mean free path we must have:

$$l(\phi_c) = \Delta/\gamma = \Delta \sqrt{\dfrac{\pi}{4}\dfrac{\tau_a}{\gamma}}\Rightarrow \tau_a = \dfrac{1}{\gamma}\dfrac{4}{\pi}$$

With the full model (total copllision rule), here is what we find (numerically):

Wwith Ran Ni model where we have the full analytical solution:

________

Let's find the critical exponent of $T$ close to the transition ($\Gamma = 1$):

$\Gamma + W_0(-\Gamma e^{-\Gamma})\sim 2(\Gamma -1)\Rightarrow T \sim (\phi-\phi_c)^2$

2 times larger than the manna universality class one :(.

___________

By talor expanding everything, we can obtain the approximate steady state:

$$(-2\gamma + \Delta^2\tau_a \omega_0^2)T - AT^{3/2}+BT^2 -C T^{5/2} + \dots =0 $$

With $A$ and $B$ positive.

While the CDP is:

and the hydrodynamical model of Mari (not expected to be close, just note:)


Our model seems equivalent to the one proposed by Mari but, for sure, ours don't behave as his...
_________

### Casimir effect in fluid:

Let's at first look at the effect of having a wall on the pressure for different fluids:


Strange...

Maybe, we will see better things with a setup for depletion forces.

_______________

### Kuroda model

The chiral model is the following:

It is interesting because it obtains hyperunfirpmity through time correlated noise and not through spatial correlated noise.
This was found in: https://arxiv.org/pdf/2302.13666 https://arxiv.org/pdf/2309.03155 https://arxiv.org/pdf/2304.14235

Basically, if the noise correlation is $\langle \chi(k, w)\chi(-k, -w)\rangle=D(k, w)$. If $D(w\rightarrow 0)=0$ or $D(k\rightarrow 0) = 0$ the system ends up hyperuniform is the damping is a global simple stupid damping. Hikeda call that temporal and spatial hyperuniformity. Large lenght or time fluctiuations are suppressed in the noise while they are not suppressed in the damping. The mix of the two destroy fluctuations at large lenght scales.

In the case of kuroda, the noise is $\propto \delta(w \pm \Omega)$ which explains the hyperuniformity as a suppression of large (and small) time fluctuations.

What about an effective temperature view? In his case, we can define an effective temperature through the FDT:

$$C_{uu}(k, w) = T_{eff}(k,w)\Im(R(k,w))/2w$$

or through a possible equipartition of potential energy:
$$\dfrac{1}{2}Kw_k^2C_{uu}(k) = \tilde T_{eff}(k)$$

Note that his system is overdamped, hence we cannot look at the kinetic energy (I think?)? In our case, both $T_{eff}$ and $\tilde T_{eff}$ matched because it did not depend on $w$. But in kuroda case, it does. Which already poses some concern to the validity of giving an effective temperature to the system.

#### Potential energy

The potential energy is:

$$\dfrac{1}{2}K\omega_k^2C_{uu}(k)=\dfrac{v_0^2}{2}\dfrac{D+Kq^2}{\Omega^2 + (D + Kq^2)^2}$$

Which of course, goes to 0 at $k=0$ when $D=0$, but it has to do so since the term in the right is more or less the structure factor (because by construction $S = k^2 C_{uu}$) and we know that it should be because of the requirement of the hyperuniformity.

#### FDT

The linear response is:

$$R=\dfrac{1}{-iw+Kq^2}$$

the noise variance is:

$$\nu(k, w) = \dfrac{v_0^2}{2}\left(\dfrac{D}{D^2 + (w+\Omega)^2}+\dfrac{D}{D^2 + (w-\Omega)^2}\right)$$

The dynamic displacement factor is:

$$C_{uu}(k, w) = \dfrac{1}{w^2 + K^2 q^4}\nu$$

It must respect the FDT:

$$C_{uu}(k, w) = T_{eff}\Im(R(k,w))\Rightarrow T_{eff}(k,w)\propto \nu(w)$$

SO, we see that: $\tilde T_{eff}(k, w)=\nu(w)$ independent on $k$. We note that of course $\int dw T_{eff}(k, w)\neq\tilde T_{eff}(k)$.

In any case, here, every mode is forced at the same temperature. But it is the temporal correlation that leads to a suppression of excitations on large lenght scales.

<details>
<summary> **This does not work!**</summary>
<br>
### Beyond the equilibrium frequency of collision for the absorbing state

We always say that the issue with our theoretical description of the continuous absorbing transition is the fact that the frequency of collision is the equilibrium one. Perhaps, the problem can be seen differently, one can start from the boltzmann equation (homoegenous profile):

$$\dfrac{\partial f(v, t)}{\partial t} - \gamma v \cdot\nabla_v f=J(v|f, f)$$
with ($g = v_1-v_2$):

$$\dfrac{J(v_2|f, f)}{g(\sigma^+)} = \int dv_1\int d\hat\sigma\left[\Theta(-\hat \sigma \cdot g - 2\Delta)(-\hat \sigma \cdot g - 2\Delta)\dfrac{f(v_2'')f(v_1'')}{\alpha^2}- \Theta(\hat \sigma \cdot g)(\hat \sigma \cdot g)f(v_2)f(v_1)\right]$$

where $v''$ are the restitutive velocities (the velocities from which, if we start from, we obtain $v_1$ and $v_2$ after the collision):

$$v_1''=v_1 - \dfrac{1+\alpha^{-1}}{2}(\hat\sigma\cdot g)\hat\sigma - \Delta\alpha^{-1}\hat\sigma$$
$$v_2''=v_2 + \dfrac{1+\alpha^{-1}}{2}(\hat\sigma\cdot g)\hat\sigma + \Delta\alpha^{-1}\hat\sigma$$

Under this form, the boltzmann equation simply takes the form of a master equation. Interestingly, without $\gamma$ we see that $g(\sigma^+)$, the enksog correction simply renormalizes the time, but with drag, it plays a major role.

The transition should also happens when $\alpha = 1$.

\begin{split}
-\gamma\int dv_2 v_2 \partial_{v_2} f(v_2) &= g(\sigma^+)\int dv_2\int dv_1\int d{\hat\sigma}\left[\Theta(-\hat \sigma \cdot g - 2\Delta)(-\hat \sigma \cdot g - 2\Delta)f(v_2'')f(v_1'')- \Theta(\hat \sigma \cdot g)(\hat \sigma \cdot g)f(v_2)f(v_1)\right]\\
\end{split}

Moreover, we will assume that $f$ can be written as $f(v)= \dfrac{1}{\sqrt{2\pi T}}e^{-v^2/(2T(\Delta, \gamma, g(\sigma^+)))}$ we know that some gausiannity is lost around the critical point, but we hope for the best :) I think though, that it will be strictly equivalent to the computation by brito.

We perform the usual change of variable to go into center of mass $\vec V$ and relative velocity $\vec g$ coordinates.

\begin{split}
\gamma \langle \vec v^2\rangle/T =& \dfrac{g(\sigma^+)}{(2\pi T)^d}\int d\vec g \int d\vec V\int d{\hat\sigma}\Theta(-\hat \sigma \cdot\vec g - 2\Delta)(-\hat \sigma \cdot\vec g - 2\Delta)e^{-\dfrac{2\vec V^2 + \vec g^2/2}{2T}}e^{-\dfrac{2\vec g^2 + 2\Delta^2-2(\vec g \cdot \hat \sigma)^2 + 2\Delta \vec g \cdot \hat \sigma}{2T}}\\ &- \Theta(\hat \sigma \cdot \vec g)(\hat \sigma \cdot \vec g)e^{-\dfrac{2\vec V^2 + \vec g^2/2}{2T}}
\end{split}

Well :) I already see the Bessel and Struve function coming, so hell no, I'm not doing that. Let's go to 1D :) :) :)

Here is 1D:

\begin{split}
\gamma =& \dfrac{g(\sigma^+)}{2\pi T}\int d g \int d V\sum_{\hat\sigma = \pm 1} \Theta(-\hat \sigma g - 2\Delta)(-\hat \sigma g - 2\Delta)e^{-\dfrac{2 V^2 + g^2/2}{2T}}e^{-\dfrac{2\Delta^2 + 2\Delta g \hat \sigma}{2T}}\\ &- \Theta(\hat \sigma g)(\hat \sigma g)e^{-\dfrac{2 V^2 + g^2/2}{2T}}\\
&= \dfrac{g(\sigma^+)e^{-\Delta^2/T}}{2\sqrt{\pi T}}\int d g \sum_{\hat\sigma = \pm 1} \Theta(-\hat \sigma g - 2\Delta)(-\hat \sigma g - 2\Delta)e^{-\dfrac{g^2/2}{2T}}e^{-\dfrac{\Delta g \hat \sigma}{T}}- \Theta(\hat \sigma g)\hat \sigma ge^{-\dfrac{g^2/2}{2T}}\\
\end{split}

_________________
</details>

_____

<details>
<summary> **Don't read, this is dumb** </summary>
<br>

### Parity of auxiliary variable: VERY STUPID DONT READ

One idea, that I guess is wrong, would be that the parity of $\Omega$ is not clear. When we perform the computation of the entropy production from the markovianized system, we assume that $\mathcal{T}\Omega(t)=\Omega(t)$. But is it the case? It is defined as:

\begin{split}
\Omega(t) &= -\int_{-\infty}^{t}dt'e^{-\dfrac{t - t'}{\tau_1}}\left(\gamma_1v(t')/\tau_1 - \sqrt{2\gamma_1 T_1/\tau_1^2} \eta(t')\right)\\
&= -\int_{-\infty}^{t}dt'e^{-\dfrac{t - t'}{\tau_1}}\Theta(t-t')\left(\gamma_1v(t')/\tau_1 - \sqrt{2\gamma_1 T_1/\tau_1^2} \eta(t')\right)\\
&= -\int_{-\infty}^{\infty}dt' K(t-t')\left(\gamma_1v(t')/\tau_1 - \sqrt{2\gamma_1 T_1/\tau_1^2} \eta(t')\right)
\end{split}

with $K(t)=e^{-t/\tau}\Theta(t)$.

We fourier transform $\Omega(t)$ and obtain:

\begin{split}
\Omega(w)&=-\overbrace{\dfrac{\tau_1}{1-i\tau_1w}}^{K(w)}\left(\alpha v(w)+\beta\eta(w)\right)\\
&=-\dfrac{1+i\tau_1w}{1+(\tau_1w)^2}(\alpha v(w) + \beta \eta(w))
&=-\frac{1}{1+(\tau_1w)^2}\left[\underbrace{\alpha v(w)}_{\text{odd}}+i\tau_1 \beta w\eta(w)+\underbrace{i\alpha wv(w)}_{\textrm{even}}+\beta\eta(w)\right]
\end{split}

I do not really know what to say about the oddness or evenness of the random variable with time reversal (one can even see that it is problematic with a simple langevin equatrion $\ddot x = -\dot x - x +\eta$) and I don't know if taking the average messes up the idea of the time reversal symetry. But it looks like that it is made of an even and an odd part...

If we can decompose $\Omega(w)$ as follow: $\Omega(w)=\Omega^r(w) + \Omega^i(w)$, then the entropy production reads:

$2\sigma =\dfrac{T - T_1}{TT_1}\langle\Omega^r(w)v(-w)\rangle+\dfrac{\langle\Omega^i(w)\Omega^r(-w)\rangle\left(\gamma T + \gamma_1T_1 + \gamma \tau_1^2Tw^2\right)-iw\gamma_1\langle\Omega^i(w)v(-w)\rangle(\gamma\tau_1T+T_1)}{\gamma\gamma_1TT_1}+ c.c.$

It's really painful to check by hand the validity of this formula :=) (More so we don't ven know the right Odd/Even decomposition of $\Omega$)

In any case, I check, it does not work with the partiy I guessed for the noise.

___________
</details>

### Metastable fluid:
oscillation with period $\gamma$. otherwise the enveloppe is $S(w) \sim w^{-2}$ like a brownian noise.

_______________________

### Effect of noise of the system.

#### Discontinuous transition

Metastability seems very very hard to obtain because the region is narrow...

$\phi = 0.23$ and no noise (last active state)
{%youtube tn6f0QlQ9Yg %}

$\phi = 0.22$ and WITH noise $T = 0.001$
{%youtube sipgNdZf4Ig %}

$\phi = 0.22$ and noise $T = 0.0001$
{%youtube N4u_CjoPGO0 %}

$\phi = 0.19$ and $T = 0.001$
{%youtube Dh6L1OO1DE8 %}

### Continuous transition noise

Without noise critical point:

{%youtube 0tGadH6NVAA %}

With small noise critical point:
{%youtube 1lsBxfz8Upc %}

THERE IS HOPE FOR THE STRUCTURE FACTOR OF THE SYSTEM?