\section{Computation of Fig.\ref{fig:figure_2}}\label{a:A}
To compute Fig.\ref{fig:figure_2}, we take
\begin{equation}
p(x,0)=\frac{\exp\left(-\frac{(x-1)^2}{2\sigma_i^2}\right)}{\sqrt{2\pi}\sigma_i}
\end{equation}
as initial condition for eq.\eqref{eq:FP}. In particular, $\sigma_i$ is such that $p(x,t)$ broadens inside the right well of $U(x)$ before developing a clean-cut second peak in the left well of the potential. We then solve eq.\eqref{eq:FP} with the Matlab \texttt{pdepe} function. With the solution we compute
\begin{equation}
P_0(t_R)=\int_{-\infty}^{0} p(x,t_R)
\end{equation}
for different refresh times $t_R$, and then we evaluate the failure probability
\begin{equation}
P_E=1-\left(1-P_0(t_R)\right)^\frac{t}{t_R}
\end{equation}
for different values of $t_R$ and $t\gg t_R$. As a last step, we use a spline fit of $P_E$ to sample $t_R$ for different values of $P_E$ and $t$. The results obtained in this way are plotted in Figure 2. These results are obtained with $T = 1/8$ which corresponds to $\tau_k = \num{5.3e4}\tau_w$. 
 
\section{Derivation of eq.\eqref{eq:equation3} and eq.\eqref{eq:equation2}}\label{a:B}
We derive here two important equations given in the main text, namely eq.\eqref{eq:equation3} and eq.\eqref{eq:equation2}.
We start with eq.\eqref{eq:equation2}. To derive it we assume that $T \ll 1$. This simplifies the mathematical description of the system as it implies that the intra-well relaxation mechanisms of the system are faster than the inter-well ones. If we are interested in intra-well mechanism only, then a satisfactory form for the dimensionless $p(x,t)$ is
\begin{subequations}\label{eq:twopeaks}
\begin{align}
p(x,t)=&p_0(x,t)+p_1(x,t)\\
p_0(x,t)=&P_0\frac{\exp\left(-\frac{(x+1)^2}{2\sigma(t)^2}\right)}{\sqrt{2\pi}\sigma(t)}\\
p_1(x,t)=&(1-P_0)\frac{\exp\left(-\frac{(x-1)^2}{2\sigma(t)^2}\right)}{\sqrt{2\pi}\sigma(t)}\\
\end{align}
\end{subequations}
where $P_0$ is, to all effects, constant over time. We substitute eq.\eqref{eq:twopeaks} in eq.\eqref{eq:FP} and then we approximate eq.\eqref{eq:duffing} with an harmonic potential by Taylor-expanding around $x=\pm 1$. This yields two distinct equations
\begin{subequations}\label{eq:doppia_FP}
\begin{align}
\pdv{p_0}{t}-8\left(p_0+(x+1)\pdv{p_0}{x}\right)-T\pdv{^2p_0}{x^2}&=0\\
\pdv{p_1}{t}-8\left(p_1+(x-1)\pdv{p_1}{x}\right)-T\pdv{^2p_1}{x^2}&=0\\
\end{align}
\end{subequations}
where used the fact that $p_1(x,t)$ ($p_0(x,t)$) can't affect the dynamics of the system in the left (right) well of $U(x)$ if $\Delta U \ll k_BT$. Eq.\eqref{eq:doppia_FP} can be combined into
\begin{equation}
\begin{aligned}
&\int_{-\infty}^{\infty}\left(\pdv{p_0}{t}-8\left(p_0+(x+1)\pdv{p_0}{x}\right)-T\pdv{^2p_0}{x^2} \right)(x+1)^2 \mathrm{d}x\\
&+\int_{-\infty}^{\infty} \left(\pdv{p_1}{t}-8\left(p_1+(x-1)\pdv{p_1}{x}\right)-T\pdv{^2p_1}{x^2}\right)(x-1)^2 \mathrm{d}x\\
&=0,
\end{aligned}
\end{equation}
which reduces to
\begin{equation}\label{eq:sigma_eq}
\pdv{\sigma(t)^2}{t}+16\sigma(t)^2-2T=0.
\end{equation}
Eq.\eqref{eq:sigma_eq} describes the time evolution of $\sigma(t)$ when intra-well relaxation mechanisms occurs. Its analytic solution for an initial condition $\sigma(0)=\sigma_i$ is
\begin{equation}\label{eq:sigma_t}
\sigma(t)=\sqrt{\frac{T}{8}+\exp\left(-16 t\right)\left(\sigma_i^2- \frac{T}{8}\right)}
\end{equation}
which is the dimensionless version of the eq.\eqref{eq:equation2} given the main text.

To compute eq.\eqref{eq:equation3} we recollect that we defined the ``refresh operation restores a non-equilibrium condition by shrinking the width of each peak of $p(x,t)$.'' without error corrections. If we assume that the refresh protocol preserves the symmetry of $U(x)$, then $p(x,t)$ can be written as eq.\eqref{eq:twopeaks} during the whole refresh  procedure. As a consequence, the sole effect of a refresh operation with duration $t_p$ is to transform
\begin{equation}\label{eq:init_cond}
p(x,t)=P_0(t)\frac{\exp\left(-\frac{(x+1)^2}{2\sigma(t)^2}\right)}{\sqrt{2\pi}\sigma(t)}+(1-P_0(t))\frac{\exp\left(-\frac{(x-1)^2}{2\sigma(t)^2}\right)}{\sqrt{2\pi}\sigma(t)}\\
\end{equation}
into
\begin{equation}\label{eq:finit_cond}
p(x,t+t_p)=P_0(t)\frac{\exp\left(-\frac{(x+1)^2}{2\sigma_i^2}\right)}{\sqrt{2\pi}\sigma_i}+(1-P_0(t))\frac{\exp\left(-\frac{(x-1)^2}{2\sigma_i^2}\right)}{\sqrt{2\pi}\sigma_i}\\
\end{equation}
where $\sigma(t)$ is given by eq.\eqref{eq:sigma_t}, $P_0(t)$ is fitted from the numerical solution of eq.\eqref{eq:FP} with eq.\eqref{eq:duffing}, and $\sigma_i=\sigma(0)$.
We now use the Gibbs entropy definition
\begin{equation}
S(t)=-k_B\int_{-\infty}^{\infty} p(x,t)\ln{p(x,t)} \mathrm{d} x
\end{equation}
to compute the entropy variation $\Delta S=S(t+t_p)-S(t)$ of the refresh protocol. Because of the $\Delta U\gg k_BT$ assumption, we have that
\begin{equation}
\begin{aligned}
&\Delta S\approx -k_B\Bigl(\int_{-\infty}^{\infty} \frac{\mathrm{e}^{-\frac{x^2}{2\sigma_i^2}}}{\sqrt{2\pi}\sigma_i}\ln\Bigl(\frac{\mathrm{e}^{-\frac{x^2}{2\sigma_i^2}}}{\sqrt{2\pi}\sigma_i} \Bigr) \mathrm{d}x \\
&-\int_{-\infty}^{\infty} \frac{\mathrm{e}^{-\frac{x^2}{2\sigma(t)^2}}}{\sqrt{2\pi}\sigma(t)}\ln\Bigl(\frac{\mathrm{e}^{-\frac{x^2}{2\sigma(t)^2}}}{\sqrt{2\pi}\sigma(t)} \Bigr) \mathrm{d}x\Bigr),
\end{aligned}
\end{equation}
which reduces to 
\begin{equation}\label{eq:DeltaS}
\Delta S\approx k_B\ln\left(\frac{\sigma_i}{\sigma(t)}\right).
\end{equation}
By using eq.\eqref{eq:DeltaS} with $t=t_R$ we obtain the eq.\eqref{eq:equation3} presented in the main text.

\section{Minimum value for $\sigma_i$}\label{a:C}
We discuss here the existence of the minimum possible value for $\sigma_i$. First of all, we observe that $\sigma_i\to0$ is a singular limit in eq.\eqref{eq:DeltaS}. This is inconsistent with the third law of thermodynamics, so there must be a minimum value for $\sigma_i$. This is given by the Heisenberg uncertainty principle. In the best case scenario this reads
\begin{equation}
\sigma_x\sigma_p=\frac{\hbar}{2}.
\end{equation}
where $\sigma_x$ ($\sigma_p$) is the uncertainty on the position $x$ (momentum $p$). According to the equipartition theorem,
\begin{equation}
\sigma_p=m{\sqrt{\langle v^2\rangle-\langle v\rangle^2}}={\sqrt{mk_B T}},
\end{equation}
so we have that
\begin{equation}\label{eq:sigma_x}
\sigma_x=\frac{\hbar}{2\sqrt{mk_BT}}
\end{equation}
Eq.\eqref{eq:sigma_x} sets the minimum possible uncertainty for $\sigma_x$. Since $\sigma_i$ describes the uncertainty of the initial $x$ value, we therefore have that $\sigma_i\geq \sigma_{iMin} =\frac{\hbar}{2\sqrt{mk_BT}}$. The existence of a $\sigma_{iMin}$ implies that, even at $t=0$, the probability of error $P_0$ is greater than zero. Clearly, this does not exclude that one can have a smaller $\sigma_i$ by accepting a larger $\sigma_p$. This would imply to operate the memory out of the thermal equilibrium, growing the dissipated energy well above the fundamental minimum.


\end{document}
