The methods described above have been coded in PETSC/SLEPc environments resorting to the matrix function (MFN) object [1]}, [2]}, [3]}, [4]}, [5]}, [6]}. MFN object also provides the restarted Krylov subspace projection method (R-KSP) [7]} to compute the product of common matrix functions and a vector.
R-KSP was used to compute reference values to check the correctness of the faster methods described above. R-KSP is a robust and well-established KSP solver with built-in convergence criteria [7]}, [3]} wherein the Krylov basis is restarted until a convergence criterion is fulfilled. We refer the reader to  for details. In our case, the application of vector \(\mathbf {b}\)  on the matrix exponential is computed for a predetermined set of times \(\Big \lbrace  t_n \Big \rbrace _{0 \le n \le L}\) :
\(\mathbf {x} = \exp ({-\mathbf {A}t_n})\mathbf {b}.\) 
