\documentclass[9pt]{beamer} \newcommand{\blue}[1]{\textcolor{blue}{#1}} % \newcommand{\red}[1]{\color{red}} \newif\iflong \longfalse \newcommand{\placeholder}{\mathord{\color{black!33}\bullet}}% \newcommand{\bu}{$\bullet \ $} \newcommand{\bi}{\begin{itemize}} \newcommand{\ei}{\end{itemize}} \renewcommand{\leq}{\leqslant} \renewcommand{\le}{\leqslant} \renewcommand{\geq}{\geqslant} \newcommand{\dt}{{\Delta t}} \newcommand\centerequation[1]{\par\smallskip\par \centerline{$\displaystyle #1$}\par \smallskip\par} \newcommand{\D}{\,\mathrm{d}} \newcommand{\cX}{\mathcal{X}} \newcommand{\E}{\expect} \newcommand{\wcL}{\widetilde{\mathcal{L}}} \newcommand{\Li}{\mathcal{K}} \newcommand{\I}{\mathrm{Id}} \newcommand{\dps}{\displaystyle} \newcommand{\red}{\color{red}} \input{header} \input{macros} \newcommand{\highlight}[2]{% \colorbox{#1!20}{$\displaystyle#2$}} \newcommand{\hiat}[4]{% \only<#1>{\highlight{#3}{#4}}% \only<#2>{\highlight{white}{#4}}% } \graphicspath{{figures/}} \AtEveryCitekey{\clearfield{pages}} \AtEveryCitekey{\clearfield{eprint}} \AtEveryCitekey{\clearfield{volume}} \AtEveryCitekey{\clearfield{number}} \AtEveryCitekey{\clearfield{month}} \addbibresource{main.bib} \title{Nonequilibrium systems and computation of transport coefficients\\[.3cm] \small \textcolor{yellow}{SINEQ Summer school}% } \author{% Urbain Vaes \texorpdfstring{\\\texttt{urbain.vaes@inria.fr}}{} } \institute{% MATHERIALS -- Inria Paris \textcolor{blue}{\&} CERMICS -- École des Ponts ParisTech } \date{\today} \begin{document} \begin{frame}[plain] \begin{figure}[ht] \centering % \includegraphics[height=1.5cm]{figures/logo_matherials.png} % \hspace{.5cm} \includegraphics[height=1.2cm]{figures/logo_inria.png} \hspace{.5cm} \includegraphics[height=1.5cm]{figures/logo_ponts.png} \hspace{.5cm} \includegraphics[height=1.5cm]{figures/logo_ERC.jpg} \hspace{.5cm} \includegraphics[height=1.5cm]{figures/logo_EMC2.png} \end{figure} \titlepage \end{frame} \begin{frame} {Some references} \begin{itemize} \itemsep.2cm \item \fullcite{MR3509213} \item \fullcite{pavliotis2011applied} \item Lecture notes by Gabriel Stoltz on computational statistical physics: \url{http://cermics.enpc.fr/~stoltz/Cours/intro_phys_stat.pdf} \end{itemize} \end{frame} \begin{frame} {Outline} \tableofcontents \end{frame} \section{Introduction} \begin{frame} {Transport coefficients} At the \alert{macroscopic} level, transport coefficients relate an external forcing to an average response expressed through some steady-state flux. \textbf{Examples:} \begin{itemize} \item The \emph{mobility} relates an external force to a velocity; \item The \emph{heat conductivity} relates a temperature difference to a heat flux; \item The \emph{shear viscosity} relates a shear velocity to a shear stress; \end{itemize} \vspace{.3cm} They can be estimated from molecular simulation at the \blue{microscopic level}. \begin{itemize} \item They are defined from \emph{nonequilibrium} dynamics; \item There are three main classes of methods to calculate them. \end{itemize} \vspace{.3cm} \textbf{Challenges we do not address:} \begin{itemize} \item Choose thermodynamical ensemble; \item Prescribe microscopic dynamics; \end{itemize} \end{frame} \begin{frame} \begin{center} \Large \color{blue} Part I: Definition and examples of nonequilibrium systems \end{center} \centering \begin{minipage}{.8\textwidth} \begin{itemize} \item Equilibrium vs nonequilibrium dynamics \item Existence of an invariant measure for nonequilibrium dynamics \item Convergence to the invariant measure \end{itemize} \end{minipage} \end{frame} \section{Equilibrium and nonequilibrium dynamics} \begin{frame} {Equilibrium and nonequilibrium dynamics} Consider a general diffusion process of the form \[ \d x_t = b(x_t) \, \d t + \sigma(x_t) \, \d W_t, \] and assume that it admits an invariant distribution $\mu$. \vspace{.2cm} \begin{definition} [Time-reversibility] A stationary ($x_0 \sim \mu$) stochastic process $(x_t)$ is time-reversible if its law is invariant under time reversal: the law of the \emph{forward paths} $(x_s)_{0 \leq s \leq t}$ coincides with the law of the \emph{backward paths} $(x_{t-s})_{0 \leq s \leq t}$. \end{definition} \vspace{.2cm} \begin{theorem} A stationary diffusion processes $x_t$ in $\real^d$ with generator $\mathcal L$ and invariant measure~$\mu$ is reversible if and only if $\mathcal L$ is self-adjoint in~$L^2(\mu)$. \end{theorem} \end{frame} \begin{frame} {Example of nonequilibrium dynamics} \begin{block}{Overdamped Langevin dynamics perturbed by a constant force term} \begin{equation} \label{eq:Langevin_F} \tag{NO} \d q_t = - \grad V(q_t) \, \d t + \alert{\eta F} + \sqrt{2} \, \d W_t \end{equation} \end{block} \begin{block}{Langevin dynamics perturbed by a constant force term} \begin{equation} \label{eq:Langevin_F} \tag{NL} \left\{ \begin{aligned} \d q_t & = M^{-1} p_t \D t, \\* \d p_t & = \bigl( -\nabla V(q_t) + {\red \eta F} \bigr) \D t - \gamma M^{-1} p_t \D t + \sqrt{\frac{2\gamma}{\beta}} \D W_t, \end{aligned} \right. \end{equation} \end{block} where \begin{itemize} \item $F \in \real^d$ with $\abs{F} = 1$ is a given direction \item $\eta \in \real$ is the strength of the external forcing. \end{itemize} Is there an invariant probability measure? \end{frame} \begin{frame} {Worked example in dimension one} Consider the perturbed overdamped Langevin dynamics with~$q_t \in \torus$ \[ \d q_t = - V'(q_t) \, \d t + {\red \eta} \, \d t + \sqrt{2} \, \d W_t, \] The associated Fokker--Planck equation reads \[ \frac{\d}{\d q}\left( \left(\frac{\d V}{\d q} - \eta\right) \rho_{\eta} + \frac{\d \rho_{\eta}}{\d q} \right) = 0. \] \begin{minipage}[t]{.45\textwidth} \vspace{.5cm} The solution is unique and given by \[ \rho_{\eta}(q) \propto \e^{-V(q)} \int_{\torus} \e^{V(q+y) - \eta y} \, \d y. \] \textbf{Example:} $\rho_{\eta}$ with $V(q) = \frac{1}{2} (1 - \cos q)$. \end{minipage} \begin{minipage}[t]{.5\textwidth} \end{minipage} \begin{minipage}[t]{.45\textwidth} \begin{figure}[ht] \centering \includegraphics[width=\linewidth]{figures/invariant_perturbed_ol.pdf} \end{figure} \end{minipage} \end{frame} \begin{frame} {Nonequilibrium overdamped Langevin dynamics} In general, can we prove existence of and convergence to an invariant measure for \[ \d q_t = - \grad V(q_t) \, \d t + \alert{\eta F} + \sqrt{2} \, \d W_t \, ? \] \begin{itemize} \item If the state space is compact (e.g. $\torus^d$), then we can apply Doeblin's theorem. \item If not, we can apply its generization, Harris' theorem. \end{itemize} \medskip Fix $t > 0$ and denote by $p\colon \mathcal E \times \mathcal B(\mathcal E)$ the Markov transition kernel \[ p(x, A) := \proba \left[ q_t \in A \, \middle| \, q_0 = x \right]. \] For an observable $\phi \colon \mathcal E \to \real$ and a probability measure $\mu$, we let \[ (\mathcal P \phi)(x) := \int_{\mathcal E} \phi(y) \, p(x, \d y), \qquad (\mathcal P^{\dagger} \mu)(A) := \int_{A} p(x, A) \, \mu(\d x). \] Note that $\mathcal P$ and $\mathcal P^{\dagger}$ are formally $L^2$ adjoints: \[ \int_{\mathcal E} (\mathcal P \phi) \, \d \mu = \int_{\mathcal E} \phi \, \d (\mathcal P^\dagger \mu). \] \end{frame} \begin{frame} {Existence of an invariant measure (1/2)} Remember that the set of probability measures with TV distance $d(\placeholder, \placeholder)$ is complete. \begin{theorem} [Doeblin's theorem] If there exists $\alpha \in (0, 1)$ and a probability measure $\eta$ such that \[ \mathcal P^{\dagger} \mu \geq \alpha \eta, \] then there exists $\mu_{\infty}$ such that $\mathcal P^{\dagger} \mu_{\infty} = \mu_{\infty}$. Furthermore $d(\mathcal P^{\dagger^n} \mu, \mu_{\infty}) \leq \alpha^n d(\mu, \mu_{\infty})$. \end{theorem} \emph{Sketch of proof.} Define the Markov transition kernel \[ \widetilde {p}(x, \placeholder) := \frac{1}{1-\alpha} p(x, \placeholder) - \frac{\alpha}{1 - \alpha} \eta(\placeholder), \] Let $\mathcal F$ denote the set of measurable functions $\phi \colon \mathcal E \to [-1, 1]$. We have \begin{align*} d(\mathcal P^\dagger \mu, \mathcal P^\dagger \nu) &= \sup_{\phi \in \mathcal F} \int_{\mathcal E} \phi(q) (\mathcal P^{\dagger} \mu - \mathcal P^{\dagger} \nu) (\d q) = \sup_{\phi \in \mathcal F} \int_{\mathcal E} \mathcal P \phi(q) \bigl(\mu - \nu\bigr) (\d q) \\ &= (1 - \alpha) \sup_{\phi \in \mathcal F} \int_{\mathcal E} \widetilde {\mathcal P} \phi(q) (\mu - \nu) (\d q) \leq (1 - \alpha) \, d(\mu, \nu). \end{align*} Conclude using Banach's fixed point theorem. \end{frame} \begin{frame} {Existence of an invariant measure (2/2)} \begin{itemize} \item Suppose that $\phi$ is uniformly bounded. Then \begin{align*} \left\lvert \mathcal P^n \phi(x) - \overline \phi \right\rvert &= \int_{\mathcal E} \mathcal P^n (\phi - \overline \phi) \, \d(\delta_x - \mu_{\infty}) = \int_{\mathcal E} (\phi - \overline \phi) \, (\mathcal P^{\dagger n} \delta_x - \mathcal P^{\dagger n} \mu_{\infty}) (\d q) \\ &\leq \norm{\phi - \overline \phi}_{L^{\infty}} (1-\alpha)^n d(\delta_x, \mu_{\infty}) \leq 2 \norm{\phi - \overline \phi}_{L^{\infty}} (1 - \alpha)^n. \end{align*} \item In molecular dynamics, this theorem can be employed for showing existence of and convergence to the invariant measure, provided that the \blue{state space is compact}. \item For \alert{noncompact state spaces}, an extension called \emph{Harris' theorem} \end{itemize} \end{frame} \begin{frame} {Existence of an invariant measure for noneq.\ dynamics} Consider the paradigmatic dynamics \begin{align*} \d q_t &= M^{-1} p_t \, \d t, \\ \d p_t &= - \bigl(\grad V(q_t) + \eta F\bigr) \, \d t - \gamma M^{-1} p_t \, \d t + \sqrt{\frac{2 \gamma}{\beta}} \, \d W_t, \end{align*} where $(q_t, p_t) = \torus^d \times \real^d$ and $F \in \real^d$ with $\abs{F} = 1$ is a given direction. \begin{figure}[ht] \centering \includegraphics[width=0.39\linewidth]{figures/intro_position.pdf} \includegraphics[width=0.39\linewidth]{figures/intro_velocity.pdf} \caption{% Marginals of the steady state solution of the Langevin dynamics with forcing } \end{figure} \end{frame} \begin{frame} {Existence of an invariant distribution} \begin{theorem} Fix~$\eta_* > 0$ and $n \geq 2$, and let $\mathcal K_n(q, p) := 1 + \abs{p}^n$. For any $\eta \in [- \eta_*, \eta_*]$, there exists a unique invariant probability measure, with a smooth density~$\psi_{\eta}(q, p)$ with respect to the Lebesgue measure. Furthermore there exists $C = C(n, \eta_*) > 0$ and $\lambda = \lambda(n, \eta_*) > 0$ such that \[ \forall \phi \in L^{\infty}_{\mathcal K_n}(\mathcal E), \qquad \left\lVert \e^{t \mathcal L_n} \phi - \int_{\mathcal E} \phi \, \psi_{\eta} \right\rVert_{L^{\infty}_{\mathcal K_n}} \leq C \e^{-\lambda t} \left\lVert \phi - \int_{\mathcal E} \phi \, \psi_{\eta} \right\rVert_{L^{\infty}_{\mathcal K_n}} \] \end{theorem} \textbf{Idea of the proof.} Show that \begin{align*} \mathcal L \mathcal K_n &\leq - c_1 \mathcal K_n(q, p) + c_2, \end{align*} for $c_1 > 0$ and $c_2 > 0$. Then apply the main theorem from~\footfullcite{MR2857021}. \end{frame} \begin{frame} {Linear response of nonequilibrium dynamics (1)} \bu The force $\eta F$ induces a non-zero velocity in the direction $F$ \medskip \bu Encoded by $\dps \expect_\eta(R) = \int_{\mathcal E} R \, \psi_\eta$ with $\dps R(q,p) = F^\t M^{-1}p$ \begin{block} {Definition of the mobility} \[ \rho_F = \lim_{\eta \to 0} \frac{\expect_\eta (R)-\expect_0 (R)}{\eta} = \lim_{\eta \to 0} \frac{\expect_\eta (R)}{\eta} \] \end{block} \medskip \bu It is {\red expected} that $\psi_\eta = f_\eta\psi_0$ with $\psi_0(q,p) = Z^{-1} \e^{-\beta H(q,p)}$ and \[ f_\eta = \mathbf{1} + \eta \mathfrak{f}_1 + \mathrm{O}(\eta^2) \] \medskip \bu In this case, $\dps \rho_F = \int_{\mathcal E} R \mathfrak{f}_1 \psi_0$ \bigskip \bu {\bf Questions:} Can the expansion for $f_\eta$ be made rigorous? What is $\mathfrak{f}_1$? \end{frame} \begin{frame} {Computation of transport coefficients} Three main classes of methods: \begin{itemize} \itemsep.2cm \item Non-equilibrium techniques \begin{itemize} \item Calculations from the steady state of a system out of equilibrium. \item Comprises bulk-driven and boundary-driven approaches. \end{itemize} \item Equilibrium techniques based on the Green--Kubo formula \[ \rho = \int_{0}^{\infty} \expect_{\mu} \bigl[\varphi(x_t) \phi(x_0)\bigr] \, \d t. \] We will derive this formula from linear response. \item Transient techniques: \end{itemize} \end{frame} \iffalse \begin{frame}\frametitle{Linear response of nonequilibrium dynamics (2)} \bu {\red Perturbative framework} where $\mathcal L_0$ considered on $L^2(\psi_0)$ is the reference \medskip \bu The invariance of $\psi_\eta$ can be written as \[ \int_{\mathcal E} (\mathcal L_\eta \varphi) \psi_\eta = 0 = \int_{\mathcal E} (\mathcal L_\eta \varphi) f_\eta \psi_0 \] \begin{block}{Fokker-Planck equation on $L^2(\psi_0)$} \centerequation{\mathcal L_\eta^* f_\eta = 0} \end{block} \bigskip \bu Formally, $\mathcal L_\eta^* f_\eta = (\mathcal L_0)^* \underbrace{\left(\I + \wcL \mathcal L_0^{-1}\right)^*f_\eta}_{=1 ?} = 0$ \medskip \bu To make the result precise, introduce $L_0^2(\psi_0) = \Pi_0 L^2(\psi_0)$ with \[ \Pi_0 f = f - \int_{\mathcal E} f \, \psi_0 \] \end{frame} \begin{frame}\frametitle{Linear response of nonequilibrium dynamics (2)} \begin{block}{Power expansion of the invariant measure} Spectral radius $r$ of the bounded operator $(\wcL \mathcal L_0^{-1})^* \in \mathcal{B}(L_0^2(\psi_0))$: \[ r = \lim_{n \to +\infty} \left\| \left[ \left(\wcL \mathcal L_0^{-1}\right)^* \right]^n \right\|^{1/n}. \] Then, for $|\eta| < r^{-1}$, the unique invariant measure can be written as $\psi_\eta = f_\eta\psi_0$, where $f_\eta \in L^2(\psi_0)$ can be expanded as \begin{equation} \label{eq:expansion_psi_xi_general} f_\eta = \left( 1+\eta (\wcL \mathcal L_0^{-1})^* \right)^{-1} \mathbf{1} = \biggl( 1 + \sum_{n=1}^{+\infty} (-\eta)^n [ (\wcL \mathcal L_0^{-1})^* ]^n \biggr) \mathbf{1}. \end{equation} \end{block} \medskip \bu Note that $\dps \int_{\mathcal E} \psi_\eta = 1$ \medskip \bu Linear response result: $\dps \rho_F = -\int_{\mathcal E} R \left[(\wcL \mathcal L_0^{-1})^* \mathbf{1}\right] \psi_0 $ \end{frame} \begin{frame}\frametitle{Elements of proof} \bu Since $\dps \frac{\gamma}{\beta} \| \nabla_p \varphi \|^2_{L^2(\psi_0)} = -\langle \mathcal L_0 \varphi,\varphi \rangle_{L^2(\psi_0)}$, it follows that \vspace{-0.2cm} \[ \| \wcL \varphi \|^2_{L^2(\psi_0)} \leq \| \nabla_p \varphi \|^2_{L^2(\psi_0)} \leq \frac{\beta}{\gamma} \| \mathcal L_0 \varphi \|_{L^2(\psi_0)} \| \varphi \|_{L^2(\psi_0)} \] \bu {\red $\mathcal L_0^{-1}$ is a well defined bounded operator on $L_0^2(\psi_0)$} (hypocoercivity + hypoelliptic regularization) \[ \| \wcL \mathcal L_0^{-1} \varphi \|^2_{L^2(\psi_0)}\leq \frac{\beta}{\gamma} \| \varphi \|_{L^2(\psi_0)} \| \mathcal L_0^{-1} \varphi \|_{L^2(\psi_0)}. \] \bu {\blue $\Pi_0 \wcL \mathcal L_0^{-1}$ is bounded on $L^2_0(\psi_0)$}, so $(\wcL \mathcal L_0^{-1})^* \Pi_0 = (\wcL \mathcal L_0^{-1})^*$ is also bounded on $L^2_0(\psi_0)$ \medskip \bu Invariance of $f_\eta$ by $\mathcal L_\eta^* = \mathcal L^* + \eta \wcL^*$ \vspace{-0.2cm} \[ \mathcal L_\eta^* f_\eta = \mathcal L_0^* \left(1 + \eta (\wcL \mathcal L_0^{-1})^* \right) f_\eta = \mathcal L_0^* \mathbf{1} = 0 \] \bu {\red Prove that $f_\eta \geq 0$} (use some ergodicity result to show that $\psi_\eta = f_\eta \psi_0$) \end{frame} \begin{frame}\frametitle{Reformulation as integrated correlation functions} \bu Conjugate response $S = \wcL^* \mathbf{1}$, equivalently $\dps \int_{\mathcal E} \left(\wcL \varphi\right) \psi_0 = \int_{\mathcal E} \varphi \, S\, \psi_0$ \medskip \begin{block}{Green--Kubo formula} For any $R \in L^2_0(\psi_0)$, \[ \lim_{\eta \to 0} \frac{\expect_\eta(R)}{\eta} = \int_0^{+\infty} \expect_0 \Big(R(q_t,p_t)S(q_0,p_0) \Big) d t, \] where $\expect_\eta$ is w.r.t. to $\psi_\eta(q,p)\,d q\, p$, while $\expect_0$ is taken over initial conditions $(q_0,p_0) \sim \psi_0$ and over all realizations of the equilibrium dynamics. \end{block} \medskip \bu For the dynamics \eqref{eq:Langevin_F}, it holds $S(q,p) = \beta R(q,p) = \beta F^T M^{-1} p$ so that \[ \rho_F = \lim_{\eta \to 0} \frac{\dps \expect_\eta (F\cdot M^{-1}p )}{\eta} = \beta \int_0^{+\infty} \expect_0 \Big( (F\cdot M^{-1}p_t) (F\cdot M^{-1}p_0) \Big) d t \] \end{frame} \begin{frame}\frametitle{Elements of proof} \bu Proof based on the following equality on $\mathcal{B}(L_0^2(\psi_0))$ \[ -\mathcal L_0^{-1} = \int_0^{+\infty} \mathrm{e}^{t \mathcal L_0} \, d t \] \bu Then, \begin{align*} \lim_{\eta \to 0} \frac{\expect_\eta(R)}{\eta} & = -\int_{\mathcal E} R \left[(\wcL \mathcal L_0^{-1})^* \mathbf{1}\right] \psi_0 = -\int_{\mathcal{E}} [\mathcal L_0^{-1}R ] [\wcL^* \mathbf{1} ] \, \psi_0 \notag \\* & = \int_0^{+\infty} \left( \int_{\mathcal{E}} \left(\mathrm{e}^{t \mathcal L_0} R\right) \, S \, \psi_0\right)dt \notag \\ & = \int_0^{+\infty} \expect \Big( R(q_t,p_t)S(q_0,p_0) \Big) d t \end{align*} \bu Note also that $S$ has average 0 w.r.t. invariant measure since \[ \int_\cX S \, d\pi = \int_\cX \wcL^* \mathbf{1} \, d\pi = \int_\cX \wcL\mathbf{1} \, d\pi = 0 \] \end{frame} \begin{frame}\frametitle{Generalization to other dynamics} \bu {\bf Possible assumptions to justify the linear response} \begin{itemize} \item existence of invariant measure with smooth density $\psi_\eta$ \item ergodicity $\dps \frac1t \int_0^t \varphi(x_s) \,d s \xrightarrow[t\to+\infty]{} \int_\cX \varphi \, \psi_\eta$ \item $\mathrm{Ker}(\mathcal L_0^*) = \mathbf{1}$ and $\mathcal L_0^*$ is invertible on~$L_0^2(\psi_0)$ \item the perturbation $\wcL$ is $\mathcal L_0$-bounded: there exist $a,b \geq 0$ such that \[ \| \wcL \varphi\|_{L^2(\psi_0)} \leq a \| \mathcal L_0 \varphi\|_{L^2(\psi_0)} + b \| \varphi\|_{L^2(\psi_0)} \] \end{itemize} \bigskip \bu {\bf When the perturbation is not sufficiently weak?} (thermal transport) \begin{itemize} \item compute $\dps \int_\cX [(\mathcal L_0+\eta\wcL)\varphi ] (1+\eta\mathfrak{f}_1)\psi_0 = \mathrm{O}(\eta^2)$ \item use a pseudo-inverse $Q_\eta = \Pi_0\mathcal L_0^{-1}\Pi_0 - \eta \Pi_0\mathcal L_0^{-1}\Pi_0\wcL\Pi_0\mathcal L_0^{-1}\Pi_0$ \item allows to prove that $\dps \int_\cX \varphi \, \psi_\eta = \int_\cX \varphi \, \psi_0 + \eta \int_\cX \varphi \, \mathfrak{f}_1 \, \psi_0 + \eta^2 r_{\varphi,\eta}$ \end{itemize} \end{frame} \begin{frame} \begin{center} \Huge{Other examples} \end{center} \end{frame} \begin{frame}\frametitle{Shear viscosity in fluids (1)} \bigskip 2D system to simplify notation: $\mathcal{D} = \left( L_x\mathbb{T} \times L_y\mathbb{T} \right)^N$ \begin{figure} \psfrag{x}{} \psfrag{z}{} \psfrag{F}{force} \center \includegraphics[height=7cm]{figures/osc_shear.eps} \end{figure} \end{frame} \begin{frame}\frametitle{Shear viscosity in fluids (2)} \bu Add a smooth {\blue nongradient force} in the $x$ direction, depending on~$y$ \begin{block}{Langevin dynamics under flow} \centerequation{\left \{ \begin{aligned} d q_{i,t} &= \frac{p_{i,t}}{m} \, dt,\\ d p_{xi,t} &= -\nabla_{q_{xi}} V(q_t) \, dt + {\red \eta F(q_{yi,t}) \, dt} - \gamma_x \frac{p_{xi,t}}{m} \, dt + \sqrt{\frac{2\gamma_x}{\beta}} \, dW^{xi}_t, \\ d p_{yi,t} &= -\nabla_{q_{yi}} V(q_t) \, dt - \gamma_y \frac{p_{yi,t}}{m} \, dt + \sqrt{\frac{2\gamma_y}{\beta}} \, dW^{yi}_t, \end{aligned} \right. } \end{block} \smallskip \bu {\red Existence/uniqueness of a smooth invariant} measure provided $\gamma_x,\gamma_y>0$ \smallskip \bu The perturbation $\dps \wcL = \sum_{i=1}^N \! F(q_{y,i}) \partial_{p_{x,i}}$ is $\mathcal{L}_0$-bounded \smallskip \bu {\blue Linear response}: $\dps \lim_{\eta \rightarrow 0} \frac{\left\langle \mathcal L_0 h \right\rangle_\eta}{\eta} = - \frac{\beta}{m} \! \left\langle \!h, \sum_{i=1}^N p_{xi} F(q_{yi}) \!\right\rangle_{L^2(\psi_0)} $ \medskip \end{frame} \begin{frame}\frametitle{Shear viscosity in fluids (3)} \bu Average {\red longitudinal velocity} $u_x(Y) = \dps \lim_{\varepsilon \to 0} \lim_{\eta \to 0} \frac{\left\langle U_x^\varepsilon(Y,\cdot)\right\rangle_\eta}{\eta}$ where \vspace{-0.3cm} \[ U_x^\varepsilon(Y,q,p) = \frac{L_y}{Nm}\sum_{i=1}^N p_{xi} \chi_{\varepsilon}\left(q_{yi}-Y\right) \] \vspace{-0.5cm} \bu Average {\red off-diagonal stress} $\dps \sigma_{xy}(Y) = \lim_{\varepsilon \to 0} \lim_{\eta \to 0} \frac{\left\langle ... \right\rangle_\eta}{\eta}$, where $... =$ \vspace{-0.4cm} \[ \hspace{-0.1cm} \frac{1}{L_x} \left( \sum_{i=1}^N \frac{p_{xi} p_{yi}}{m}\chi_{\varepsilon}\left(q_{yi}-Y\right) - \! \! \! \! \! \! \! \! \sum_{1 \leq i < j \leq N} \! \! \! \! \mathcal{V}'(|q_i-q_j|)\frac{ q_{xi}-q_{xj}}{|q_i-q_j|} \!\int_{q_{yj}}^{q_{yi}} \!\chi_{\varepsilon}(s-Y) \, ds \right) \] \bu {\blue Local conservation} of momentum\footnote{Irving and Kirkwood, {\it J. Chem. Phys.} {\bf 18} (1950)}: replace $h$ by $U_x^\varepsilon$ (with $\overline{\rho} = N/|\mathcal{D}|$) \[ \frac{d\sigma_{xy}(Y)}{dY} + \gamma_{x} \overline{\rho} u_x(Y) = \overline{\rho} F(Y) \] \end{frame} \begin{frame} \frametitle{Shear viscosity in fluids (4)} \bu {\blue Definition} $\sigma_{xy}(Y) := -\eta(Y)\dfrac{du_x(Y)}{dY}$, {\red closure} assumption $\eta(Y) = \eta > 0$ \begin{block}{Velocity profile in Langevin dynamics under flow} \centerequation{-\eta u_x''(Y) + \gamma_x \overline{\rho} u_x(Y) = \overline{\rho} F(Y)} \end{block} \bigskip \hspace{-0.5cm} \begin{minipage}{6cm} \psfrag{F}{{\scriptsize $F$}} \psfrag{U}{{\scriptsize $u$}} \psfrag{Y}{{\scriptsize $\ \ Y$}} \psfrag{v}{{\scriptsize value}} \includegraphics[width=6cm]{figures/ux5.eps} \end{minipage} \hspace{-0.5cm} \begin{minipage}{6cm} \psfrag{Y}{} \psfrag{v}{{\scriptsize value}} \psfrag{S}{{\scriptsize $\sigma_{xy}$}} \psfrag{D}{{\scriptsize $-\nu u'$}} \includegraphics[width=6cm]{figures/dux5.eps} \end{minipage} \end{frame} \begin{frame}\frametitle{Thermal transport in one-dimensional chains (1)} \bu Atoms at positions $q_0,\dots,q_N$ with $q_0 = 0$ fixed \medskip \bu Hamiltonian $\dps H(q,p) = \sum_{i=1}^N \frac{p_i^2}{2} + \sum_{i=1}^{N-1} v(q_{i+1} - q_i) + v(q_1)$ \begin{block}{Hamiltonian dynamics with Langevin thermostats at the boundaries} \centerequation{ \left\{ \begin{aligned} dq_i & = p_i \, dt \\ dp_i & = \Big( v'(q_{i+1}-q_i) - v'(q_i-q_{i-1}) \Big) dt,\qquad i\neq 1, N \\[-3pt] dp_1 & = \Big( v'(q_2-q_1) - v'(q_1) \Big) dt - \gamma p_1 \, dt + \sqrt{2\gamma (T{\red +\Delta T})} \, dW^1_t\\[-3pt] dp_N & = - v'(q_N-q_{N-1}) \, dt - \gamma p_N \, dt + \sqrt{2\gamma (T{\red -\Delta T})} \, dW^N_t\\[-5pt] \end{aligned} \right. } \end{block} \medskip \bu {\red Perturbation} $\wcL = \gamma( \partial_{p_1}^2 - \partial_{p_N}^2)$ (not $\mathcal L_0$-bounded...) \medskip \bu Proving the existence/uniqueness of the invariant measure already requires quite some work\footnote{P. Carmona, {\emph Stoch. Proc. Appl.} (2007)} \bigskip \end{frame} \begin{frame}\frametitle{Thermal transport in one-dimensional chains (2)} \bu Response function: {\blue Total energy current} \begin{block}{} \centerequation{J = \frac{1}{N-1}\sum_{i=1}^{N-1} j_{i+1,i}, \qquad j_{i+1,i} = -v'(q_{i+1}-q_i)\frac{p_i+p_{i+1}}{2}} \end{block} \smallskip \bu Motivation: Local conservation of the energy (in the bulk) \[ \frac{d\varepsilon_i}{dt} = j_{i-1,i} - j_{i,i+1}, \qquad \varepsilon_i = \frac{p_i^2}{2} + \frac12 \Big( v(q_{i+1}-q_{i}) + v(q_i-q_{i-1}) \Big) \] \bu Definition of the {\blue thermal conductivity}: linear response \[ \kappa_N = \lim_{\Delta T \to 0} \frac{\langle J \rangle_{\Delta T}}{\Delta T/N} = 2\beta^2 \frac{N}{N-1}\int_0^{+\infty} \sum_{i=1}^{N-1} \expect\Big(j_{2,1}(q_t,p_t)j_{i+1,i}(q_0,p_0)\Big)\, dt \] \medskip \bu {\blue Synthetic dynamics}: fixed temperatures of the thermostats but external forcings $\to$ {\red bulk driven dynamics} with $\wcL^* = -\wcL + c J$ \end{frame} \begin{frame} \begin{center} \Huge{Error estimates on} \\ \bigskip \Huge{the computation of} \\ \bigskip \Huge{transport coefficients} \end{center} \end{frame} \begin{frame}\frametitle{Reminder: Error estimates in Monte Carlo simulations} \bu General SDE $dx_t = b(x_t)\,dt + \sigma(x_t) \, dW_t$, invariant measure $\pi$ \bigskip \bu {\red Discretization} $x^{n} \simeq x_{n\dt}$, {\blue invariant measure $\pi_\dt$}. For instance, \[ x^{n+1} = x^n + \dt\,b(x^n)+\sqrt{\dt}\,\sigma(x^n)\,G^n, \qquad G^n \sim \mathcal{G}(0,{\rm Id}) \ \mathrm{i.i.d.} \] \medskip \bu {\blue Ergodicity} of the numerical scheme with invariant measure~$\pi_\dt$ \[ \frac{1}{N_{\rm iter}} \sum_{n=1}^{N_{\rm iter}} A(x^n) \xrightarrow[N_{\rm iter}\to+\infty]{} \int_\cX A(x) \, \pi_\dt(dx) \] \begin{block}{Error estimates for {\red finite} trajectory averages} \[ \widehat{A}_{N_{\rm iter}} = \frac{1}{N_{\rm iter}} \sum_{n=1}^{N_{\rm iter}} A(x^n) = \expect_\pi(A) + \underbrace{C\dt^\alpha}_{\rm bias} + \underbrace{\frac{\sigma_{A,\dt}}{\sqrt{N_{\rm iter}\dt}} \mathscr{G}}_\mathrm{statistical~error} \] \end{block} \smallskip \bu Bias $\expect_{\pi_\dt}(A)-\expect_\pi(A) \longrightarrow$ {\bf Focus today} \medskip \end{frame} \begin{frame}\frametitle{Weak type expansions} \bu Numerical scheme = {\red Markov chain} characterized by {\blue evolution operator} \[ P_\dt \varphi(x) = \expect\Big( \varphi\left(x^{n+1}\right)\Big| x^n = x\Big) \] where $(x^n)$ is an approximation of $(x_{n \dt})$ \bigskip \bu (Infinitely) Many possibilities! Numerical analysis allows to {\blue discriminate} \medskip \bu Standard notions of error: {\red fixed integration time $T < +\infty$} \begin{itemize} \item {\blue Strong error} $\dps \sup_{0 \leq n \leq T/\dt} \expect | X^n - X_{n\dt} | \leq C \dt^p$ \item {\blue Weak error}: $\dps \!\!\!\! \sup_{0 \leq n \leq T/\dt} \Big| \expect\left[\varphi\left(X^n\right)\right] - \expect\left[\varphi\left(X_{n\dt}\right)\right] \Big| \leq C \dt^p$ (for any $\varphi$) %\item ``mean error'' \emph{vs.} ``error of the mean'' \end{itemize} %\medskip %\bu Example: for Euler-Maruyama, weak order~1, strong order $1/2$ (1 when $\sigma$ constant) %\medskip \begin{block}{$\dt$-expansion of the evolution operator} \centerequation{P_\dt \varphi = \varphi + \dt \, \mathcal A_1 \varphi + \dt^2 \mathcal A_2 \varphi + \dots + \dt^{p+1} \mathcal A_{p+1} \varphi + \dt^{p+2} r_{\varphi,\dt}} \end{block} \smallskip \bu {\red Weak order}~$p$ when $\mathcal A_k = \mathcal L^k/k!$ for $1 \leq k \leq p$ \end{frame} \begin{frame}\frametitle{Example: Euler-Maruyama, weak order~1} \medskip \bu Scheme $x^{n+1} = \Phi_\dt(x^n,G^n) = x^n + \dt\,b(x^n)+\sqrt{\dt}\,\sigma(x^n)\,G^n$ \bigskip \bu Note that $P_\dt \varphi(x) = \expect_G\left[ \varphi\big(\Phi_\dt(x,G)\big) \right]$ \bigskip \bu Technical tool: {\blue Taylor expansion} \vspace{-0.2cm} \[ \varphi(x + \delta) = \varphi(x) + \delta^T \nabla \varphi(x) + \frac12 \delta^T \nabla^2\varphi(x) \delta + \frac16 D^3\varphi(x):\delta^{\otimes 3} + \dots \] \medskip \bu Replace $\delta$ with $\sqrt{\dt}\, \sigma(x)\,G + \dt\,b(x)$ and {\blue gather in powers of $\dt$} \[ \begin{aligned} \varphi\big(\Phi_\dt(x,G)\big) & = \varphi(x) + \sqrt{\dt}\, \sigma(x)\,G \cdot \nabla \varphi(x) \\ & \ \ \ + \dt \left(\frac{\sigma(x)^2}{2} G^T \left[\nabla^2\varphi(x)\right]G + b(x)\cdot\nabla \varphi(x) \right) + \dots \end{aligned} \] \medskip \bu Taking {\blue expectations w.r.t. $G$} leads to \[ P_\dt\varphi(x) = \varphi(x) + \dt \underbrace{\left(\frac{\sigma(x)^2}{2} \Delta \varphi(x) + b(x)\cdot\nabla \varphi(x) \right)}_{= \mathcal{L}\varphi(x)} + \mathrm{O}(\dt^2) \] \end{frame} \begin{frame}\frametitle{Error estimates on the invariant measure (equilibrium)} \bu {\red Assumptions} on the operators in the weak-type expansion \begin{itemize} \item invariance of $\pi$ by $\mathcal A_k$ for $1 \leq k \leq p$, namely $\dps \int_\cX \mathcal A_k \varphi \, d\pi = 0$ \item $\dps \int_\cX \mathcal A_{p+1}\varphi \, d\pi = \int_\cX g_{p+1} \varphi \, d\pi$ (\textit{i.e.} $g_{p+1} = \mathcal A_{p+1}^* \mathbf{1}$) \end{itemize} \begin{block}{Error estimates on $\pi_\dt$} \centerequation{ \int_\cX \varphi \, d\pi_\dt = \int_\cX \varphi \Big(1+\dt^{p}f_{p+1}\Big) d\pi + \dt^{p+1} R_{\varphi,\dt} } \end{block} \medskip \bu In fact, $f_{p+1} = -\left( \mathcal A_1^* \right)^{-1} g_{p+1}$ \begin{itemize} \item when $\mathcal A_1 = \mathcal L$, the first order correction can be {\red estimated} by some integrated correlation function as $\dps \int_0^{+\infty} \expect\Big(\varphi(x_t)g_{p+1}(x_0)\Big) \, dt$ \item in general, first order term can be removed by Romberg extrapolation \end{itemize} \medskip \bu Error on invariant measure can be {\blue (much) smaller} than the weak error \end{frame} %----------------------------------------------------------- \begin{frame}\frametitle{Sketch of proof (1)} {\bf Step~1: Establish the error estimate for $\varphi \in \mathrm{Ran}(P_\dt-\I)$} \medskip \bu Idea: $\pi_\dt = \pi (1 + \dt^p f_{p+1} + \dots)$ \medskip \bu by definition of $\pi_\dt$ \[ \int_\cX \left[ \left(\frac{P_\dt-\I}{\dt}\right) \psi \right] d\pi_\dt = 0 \] \bu compare to first order correction to the invariant measure \[ \begin{aligned} & \int_\cX \left[ \left(\frac{P_\dt-\I}{\dt}\right)\psi\right] (1+\dt^{p}f_{p+1})\, d\pi \\ & \qquad = \dt^{p} \int_\cX \Big( \mathcal A_{p+1}\psi + (\mathcal A_1 \psi) f_{p+1} \Big) d\pi + \mathrm{O}\left(\dt^{p+1}\right) \\ & \qquad = \dt^p \int_\cX \Big( g_{p+1} + \mathcal A_1^* f_{p+1} \Big) \psi \, d\pi + \mathrm{O}\left(\dt^{p+1}\right) \end{aligned} \] \begin{block}{} Suggests $f_{p+1} = -\left( \mathcal A_1^* \right)^{-1} g_{p+1}$ \end{block} \end{frame} %----------------------------------------------------------- \begin{frame}\frametitle{Sketch of proof (2)} {\bf Step~2: Define an approximate inverse} \medskip \bu Issue: derivatives of $(\I-P_\dt)^{-1}\varphi$ are not controlled \bigskip \bu Consider $\dps \left(\Pi \frac{P_\dt-\I}{\dt} \Pi\right) Q_\dt\psi = \psi + \dt^{p+1} \widetilde{r}_{\psi,\dt}$ where \vspace{-0.2cm} \[ \Pi \varphi = \varphi - \int_\cX \varphi \, d\pi \] \bu Idea of the construction: truncate the formal series expression \[ (A + \dt \, B)^{-1} = A^{-1} - \dt \, A^{-1}B A^{-1} + \dt^{2} \, A^{-1}B A^{-1}B A^{-1} + \dots \] \bigskip {\bf Step~3: Conclusion} \medskip \bu Write the invariances with $\dps \Pi \left(\frac{P_\dt-\I}{\dt}\right) \Pi \psi$ instead of $\dps \left(\frac{P_\dt-\I}{\dt}\right) \psi$ \medskip \bu Replace $\psi$ by $Q_\dt \varphi$, and gather in~$R_{\varphi,\dt}$ all the higher order terms \end{frame} \begin{frame}\frametitle{Examples of splitting schemes for Langevin dynamics (1)} \bu Example: Langevin dynamics, discretized using a {\blue splitting} strategy \[ A = M^{-1} p \cdot \nabla_q, \quad B_\eta = \Big(-\nabla V(q) + \eta\,F\Big)\cdot \nabla_p, \quad C = -M^{-1} p \cdot \nabla_p + \frac1\beta \Delta_p \] \bu Note that $\mathcal L_\eta = A + B_\eta + \gamma C$ \medskip \bu Trotter splitting $\to$ weak order 1 \[ P^{ZYX}_\dt = \e^{\dt Z} \e^{\dt Y} \e^{\dt X} = \e^{\dt \mathcal L} + \mathrm{O}(\dt^2) \] \bu Strang splitting $\to$ {\blue weak order 2} \[ P^{ZYXYZ}_\dt = \e^{\dt Z/2} \e^{\dt Y/2} \e^{\dt X} \e^{\dt Y/2} \e^{\dt Z/2} = \e^{\dt \mathcal L} + \mathrm{O}(\dt^3) \] \bu Other category: {\red Geometric Langevin}\footnote{N.~Bou-Rabee and H.~Owhadi, {\em SIAM J. Numer. Anal.} (2010)} algorithms, \textit{e.g.} $P_\dt^{\gamma C,A,B_\eta,A}$ \\ $\to$ weak order 1 but measure preserved at order 2 in $\dt$ \end{frame} \begin{frame}\frametitle{Examples of splitting schemes for Langevin dynamics (2)} \small \bu $P_\dt^{B_\eta,A,\gamma C}$ corresponds to %\begin{equation} %\label{eq:Langevin_splitting} $\dps \left\{ \begin{aligned} \widetilde{p}^{n+1} & = p^n + \Big(-\nabla V(q^{n}) + \eta F\Big)\dt, \\ q^{n+1} & = q^n + \dt \, M^{-1} \widetilde{p}^{n+1}, \\ p^{n+1} & = \alpha_\dt \widetilde{p}^{n+1} + \sqrt{\frac{1-\alpha^2_\dt}{\beta}M} \, G^n \end{aligned} \right.$ \\[5pt] %\end{equation} where $G^n$ are i.i.d. Gaussian and $\alpha_\dt = \exp(-\gamma M^{-1} \dt)$ \bigskip \bu $P^{\gamma C,B_\eta,A,B_\eta,\gamma C}_\dt$ for %\[ $\dps \left\{ \begin{aligned} \widetilde{p}^{n+1/2} & = \alpha_{\dt/2} p^{n} + \sqrt{\frac{1-\alpha_{\dt}}{\beta}M} \, G^{n}, \\ p^{n+1/2} & = \widetilde{p}^{n+1/2} + \frac{\dt}{2} \Big( -\nabla V(q^{n})+\eta F\Big), \\ q^{n+1} & = q^n + \dt \, M^{-1} p^{n+1/2}, \\ \widetilde{p}^{n+1} & = p^{n+1/2} + \frac{\dt}{2} \Big(- \nabla V(q^{n+1}) +\eta F\Big), \\ p^{n+1} & = \alpha_{\dt/2} \widetilde{p}^{n+1} + \sqrt{\frac{1-\alpha_{\dt}}{\beta}M} \, G^{n+1/2} \end{aligned} \right.$ %\] \end{frame} \begin{frame}\frametitle{Error estimates on linear response} \begin{block}{Error estimates for nonequilibrium dynamics} There exists a function $f_{\alpha,1,\gamma} \in H^1(\mu)$ such that \vspace{-0.3cm} \[ \int_{\mathcal E} \psi \, d{\mu}_{\gamma,\eta,\dt} = \int_{\mathcal E} \psi \Big(1+ \eta f_{0,1,\gamma} + \dt^\alpha f_{\alpha,0,\gamma} + \eta \dt^\alpha f_{\alpha,1,\gamma} \Big) d{\mu} + r_{\psi,\gamma,\eta,\dt}, \] where the remainder is compatible with linear response \vspace{-0.1cm} \[ \left|r_{\psi,\gamma,\eta,\dt}\right| \leq K(\eta^2 + \dt^{\alpha+1}), \qquad \left|r_{\psi,\gamma,\eta,\dt} - r_{\psi,\gamma,0,\dt}\right| \leq K \eta (\eta + \dt^{\alpha+1}) \] \end{block} \medskip \bu Corollary: error estimates on the {\blue numerically computed mobility} \[ \begin{aligned} \rho_{F,\dt} & = \lim_{\eta \to 0} \frac{1}{\eta} \left(\int_{\mathcal E} F^T M^{-1} p \, \mu_{\gamma,\eta,\dt}(d{q}\,d{p}) - \int_{\mathcal E} F^T M^{-1} p \, \mu_{\gamma,0,\dt}(d{q}\,d{p}) \right) \\ & = \rho_{F} + \dt^\alpha \int_{\mathcal E} F^T M^{-1} p \, f_{\alpha,1,\gamma} \, d{\mu} + \dt^{\alpha+1} r_{\gamma,\dt} \end{aligned} \] \bu Results in the {\red overdamped} limit\footnote{B.~Leimkuhler, C.~Matthews and G.~Stoltz, {\em IMA J. Numer. Anal.} (2015)} \bigskip \end{frame} \begin{frame}\frametitle{Numerical results} \begin{figure} \begin{center} \includegraphics[width=6.2cm]{figures/LR.eps} \includegraphics[width=6.2cm]{figures/mobility_Langevin.eps} \end{center} \end{figure} \small {\bf Left:} Linear response of the average velocity as a function of $\eta$ for the scheme associated with $P_\dt^{\gamma C, B_\eta,A,B_\eta, \gamma C}$ and $\dt = 0.01, \gamma = 1$. \\ \smallskip {\bf Right:} Scaling of the mobility $\nu_{F,\gamma,\dt}$ for the first order scheme $P_\dt^{A,B_\eta,\gamma C}$ and the second order scheme $P_\dt^{\gamma C, B_\eta,A,B_\eta, \gamma C}$. \end{frame} %----------------------------------------------------------- \begin{frame}\frametitle{Error estimates on Green-Kubo formulas (1)} \bu For methods of {\bf weak order}~1, {\red Riemman sum} ($\phi,\varphi$ average 0 w.r.t. $\pi$) \vspace{-0.2cm} \[ \begin{aligned} & \int_0^{+\infty} \expect \Big( \phi(x_t) \varphi(x_0) \Big) d{t} = \dt \sum_{n=0}^{+\infty} \expect_\dt \left(\Pi_\dt \phi\left(x^{n}\right)\varphi\left(x^0\right)\right) + \mathrm{O}(\dt) \\[-7pt] & \mathrm{where} \ \Pi_\dt \phi = \phi - \int_\cX \phi \, d\pi_\dt \end{aligned} \] \bu Correlation approximated in practice using $K$ independent realizations %\bi %\item truncating the integration (decay estimates) %\item using empirical averages ($K$ independent realizations) \vspace{-0.2cm} \[ \begin{aligned} & \expect_\dt \left(\Pi_\dt \phi\left(x^{n}\right)\varphi\left(x^0\right)\right) \simeq \frac{1}{K}\sum_{m=1}^{K} \left( \phi(x^{n,k}) - \overline{\phi}^{n,K} \right)\left( \varphi(x^{n,k}) - \overline{\varphi}^{n,K} \right) \\[-10pt] & \mathrm{where} \ \overline{\phi}^{n,K} = \frac1K \sum_{m=1}^{K} \phi(x^{n,k}) \end{aligned} \] \bu For methods of {\bf weak order} 2, {\blue trapezoidal rule} \vspace{-0.1cm} \[ \begin{aligned} \int_0^{+\infty} \expect \Big( \phi(x_t) \varphi(x_0) \Big) d{t} & = \frac{\dt}{2} \expect_\dt \left(\Pi_\dt \phi\left(x^{0}\right)\varphi\left(x^0\right)\right) \\ & \ \ + \dt \sum_{n=1}^{+\infty} \expect_\dt \left(\Pi_\dt \phi\left(x^{n}\right)\varphi\left(x^0\right)\right) + \mathrm{O}(\dt^2) \end{aligned} \] %\bu Allows to quantify the variance $\dps \frac{\sigma^2_{A,\dt}}{N_{\rm iter}\dt} \simeq \frac{\dps 2 \int_0^{+\infty} \expect\left[\delta A(x_t)\delta A(x_0)\right] \, dt}{T}$ where $T = N_{\rm iter}\dt$ \end{frame} %----------------------------------------------------------- \begin{frame}\frametitle{Error estimates on Green-Kubo formulas (2)} \bu Error of {\red order~$\alpha$ on invariant measure}: $\dps \int_\cX \psi \, d{\pi}_\dt = \int_\cX \psi \, d{\pi} + \mathrm{O}(\dt^\alpha)$ \medskip \bu Expansion of the evolution operator ($p+1 \geq \alpha$ and $\mathcal A_1 = \mathcal L$) \[ P_\dt \varphi = \varphi + \dt \, \mathcal L \varphi + \dt^2 \mathcal A_2 \varphi + \dots + \dt^{p+1} \mathcal A_{p+1} \varphi + \dt^{p+2} r_{\varphi,\dt} \] \begin{block}{Ergodicity of the numerical scheme} \centerequation{ \forall n \in \mathbb{N}, \qquad \left\| P_\dt^n \right \|_{\mathcal{B}(L^\infty_{\Li_s,\dt})} \leq C_s \e^{-\lambda_s n\dt} } where $\mathcal{K}_s$ is a Lyapunov function ($1+|p|^{2s}$ for Langevin) and \[ L^\infty_{\Li_s,\dt} = \left\{ \frac{\varphi}{\mathcal{K}_s} \in L^\infty(\cX), \ \int_\cX \varphi \, d\pi_\dt = 0\right\} \] \end{block} \bu Proof: Lyapunov condition + uniform-in-$\dt$ minorization condition\footnote{M. Hairer and J. Mattingly, \emph{Progr. Probab.} (2011)} \end{frame} %----------------------------------------------------------- \begin{frame}\frametitle{Error estimates on Green-Kubo formulas (3)} \begin{block}{Error estimates on integrated correlation functions} Observables $\varphi,\psi$ with average~0 w.r.t. invariant measure~$\pi$ \[ \int_0^{+\infty} \expect \Big( \psi(x_t) \varphi(x_0) \Big) d{t} = \dt \sum_{n=0}^{+\infty} \expect_\dt \left(\widetilde{\psi}_{\dt,\alpha}\left(x^{n}\right)\varphi\left(x^0\right)\right) + \dt^\alpha r^{\psi,\varphi}_\dt, \] where $\expect_\dt$ denotes expectations w.r.t. initial conditions $x_0 \sim \pi_\dt$ and over all realizations of the Markov chain $(x^n)$, and \[ \widetilde{\psi}_{\dt,\alpha} = \psi_{\dt,\alpha} - \int_\cX \psi_{\dt,\alpha} \, d\pi_\dt\] with $\dps \psi_{\dt,\alpha} = \Big(\I + \dt \,\mathcal A_2 \mathcal L^{-1} + \dots + \dt^{\alpha-1} \mathcal A_{\alpha}\mathcal L^{-1} \Big)\psi$ \end{block} \bu Useful when $\mathcal A_k \mathcal L^{-1}$ can be computed, \emph{e.g.} $\mathcal A_k = a_k \mathcal L^{k}$ \medskip \bu Reduces to trapezoidal rule for second order schemes \end{frame} %----------------------------------------------------------- \begin{frame}\frametitle{Sketch of proof (1)} \bu Define $\dps \Pi_\dt \varphi = \varphi - \int_\cX \varphi \, d\pi_\dt$ \smallskip \bu Since $\mathcal L^{-1}\psi$ has average~0 w.r.t.~$\pi$, introduce $\pi_\dt$ as \begin{align*} \int_\cX \left(-\mathcal L^{-1} \psi\right) \varphi \, d{\pi} & = \int_\cX \left(-\mathcal L^{-1} \psi\right) \Pi_\dt \varphi \, d{\pi} \nonumber \\ %& = \int_\cX \left(-\mathcal L^{-1} \psi\right) \Pi_\dt \varphi \, d{\pi}_\dt + \dt^\alpha r^{\psi,\varphi}_\dt, \nonumber \\ & = \int_\cX \Pi_\dt \left(-\mathcal L^{-1} \psi\right) \Pi_\dt \varphi \, d{\pi}_\dt + \dt^\alpha r^{\psi,\varphi}_\dt, \end{align*} \bu Rewrite $-\Pi_\dt \mathcal L^{-1}$ in terms of $P_\dt$ as \[ \begin{aligned} & -\Pi_\dt \mathcal L^{-1} \psi = -\Pi_\dt \left(\dt\sum_{n=0}^{+\infty} P_\dt^n \right) \Pi_\dt \left(\frac{\I - P_\dt}{\dt}\right) \mathcal L^{-1} \psi \\ & \ \ = \dt \left(\sum_{n=0}^{+\infty} \left[ \Pi_\dt P_\dt \Pi_\dt\right]^n \right) \left(\mathcal L + \dots + \dt^{\alpha-1} S_{\alpha-1} + \dt^\alpha \widetilde{R}_{\alpha,\dt}\right) \mathcal L^{-1} \psi, \\ & \ \ = \dt \sum_{n=0}^{+\infty} \left[ \Pi_\dt P_\dt \Pi_\dt\right]^n \widetilde{\psi}_{\dt,\alpha} + \dt^\alpha \left(\frac{\I - P_\dt}{\dt}\right)^{-1} \Pi_\dt \widetilde{R}_{\alpha,\dt} \mathcal L^{-1} \psi. \end{aligned} \] \end{frame} %----------------------------------------------------------- \begin{frame}\frametitle{Sketch of proof (2)} \bu Uniform resolvent bounds $\dps \left\| \left(\frac{\I - P_\dt}{\dt}\right)^{-1} \right \|_{\mathcal{B}(L^\infty_{\Li_s,\dt})} \leq \frac{C_s}{\lambda_s}$ \medskip \bu Coming back to the initial equality, \[ \int_\cX \left(-\mathcal L^{-1} \psi\right) \varphi \, d{\pi} = \dt \int_\cX \sum_{n=0}^{+\infty} \left( \Pi_\dt P_\dt^n \widetilde{\psi}_{\dt,\alpha} \right) \left( \Pi_\dt \varphi \right) d{\pi}_\dt + \mathrm{O}\left(\dt^\alpha\right) \] \bu Rewrite finally \[ \begin{aligned} \int_\cX \sum_{n=0}^{+\infty} \left( \Pi_\dt P_\dt^n \widetilde{\psi}_{\dt,\alpha} \right)\left( \Pi_\dt \varphi \right) d{\pi}_\dt & = \int_\cX \sum_{n=0}^{+\infty} \left(P_\dt^n \widetilde{\psi}_{\dt,\alpha} \right) \varphi \, d{\pi}_\dt \\ & = \sum_{n=0}^{+\infty} \expect_\dt \left(\widetilde{\psi}_{\dt,\alpha}\left(q^{n},p^{n}\right)\varphi\left(q^0,p^0\right)\right) \end{aligned} \] \end{frame} \begin{frame}\frametitle{Numerical results} \vspace{-0.5cm} \begin{figure} \begin{center} \includegraphics[width=11.8cm]{figures/error_diffusion.eps} %\includegraphics[width=8.2cm]{figures/error_diffusion_zoom.eps} \end{center} \end{figure} \end{frame} \begin{frame} \begin{center} \Huge{Conclusion and perspectives} \end{center} \end{frame} \begin{frame} {Main points: recall the outline!} \bu {\bf Definition and examples of nonequilibrium systems} \bigskip \bu {\bf Computation of transport coefficients} \begin{itemize} \item a survey of computational techniques \item linear response theory \item relationship with Green-Kubo formulas \end{itemize} \bigskip \bu {\bf Elements of numerical analysis} \begin{itemize} \item estimation of biases due to timestep discretization \item {\blue (largely) open issue: variance reduction} \item {\red (not discussed) use of non-reversible dynamics to enhance sampling} \end{itemize} \end{frame} \begin{frame} {Variance reduction techniques?} \bu {\blue Importance sampling?} Invariant probability measures $\psi_\infty$, $\psi_\infty^A$ for \[ dq_t = b(q_t) \, dt + \sigma dW_t, \qquad dq_t = \Big( b(q_t) + \nabla A(q_t) \Big) dt + \sigma dW_t \] In general {\red $\psi_\infty^A \neq Z^{-1} \psi_\infty \mathrm{e}^{A}$} (consider $b(q) = F$ and $A = \widetilde{V}$) \bigskip \bu {\blue Stratification?} (as in TI...) Consider $q \in \mathbb{T}^2$, $\psi_\infty = \mathbf{1}_{\mathbb{T}^2}$ \[ \left \{ \begin{aligned} dq^1_t & = \partial_{q_2}U(q^1_t,q^2_t) + \sqrt{2} \, dW_t^1 \\ dq^2_t & = - \partial_{q_1}U(q^1_t,q^2_t) + \sqrt{2} \, dW_t^2 \end{aligned} \right. \] Constraint $\xi(q) = q_2$, {\red constrained dynamics} \[ dq^1_t = f(q^1_t) \, dt + \sqrt{2} \, dW_t^1, \qquad f(q^1) = \partial_{q_2}U(q^1,0). \] Then $\dps \psi_\infty(q^1) = Z^{-1} \int_0^{1} \e^{V(q^1+y)-V(q^1)-Fy} \, dy \neq \mathbf{1}_{\mathbb{T}}(q^1)$ \\ where $\dps F = \int_0^1 f$ and $\dps V(q^1) = -\int_0^{q^1} (f(s) - F) \, ds$ \end{frame} \fi \end{document} % vim: ts=4 sw=4