The time integration of stiff systems of differential equations as \[u'(t) = F(t,u(t)), \quad u(0) = u_0,\] constitutes a heated topic in numerical analysis. In particular, exponential integrators drew a great deal of attention. In fact, similarly to implicit methods, these methods show good stability properties, allowing integration with large time steps. Each exponential integration step of length \(\tau\) (out of hundreds, thousands or even millions steps) consists of the same operation: let \(\textbf{A}\) be the linear part of \(F(t,u(t))\) (or, say, its Jacobian), one shall compute \[u^{n+1} := \varphi_0(\theta_0 \tau\textbf{A}) {v}_0 + \varphi_1(\theta_1 \tau\textbf{A}) {v}_1 + \ldots + \varphi_p(\theta_p \tau\textbf{A}) {v}_p,\] where \(\theta_0, \theta_1, \ldots, \theta_p\) are fixed scalars, \[\varphi_{\ell}( x ) := \sum_{j = 0}^\infty \frac{x^j}{(j+{\ell})!}, \quad {\ell} = 0,1,\ldots,p\] and the vectors \(v_0, v_1, \ldots, v_p \in \mathbb{C}^N\) are obtained, in a recursive fashion, as functions of linear combinations of \(\varphi\)–functions applied to vectors connected to the current state of the system \(u^{n}\).

The authors exploited this peculiarity of exponential integrators and recent advancements in numerical analysis to build a routine for computing the action of the matrix \(\varphi\)–functions arising in the exponential integration steps, called bamphi, able to recycle the information gathered through the exponential integration steps and to reach high levels of speed and accuracy. In this presentation, we outline some of bamphi’s main features and ideas.