The PHIKS algorithm computes actions of phi-functions for Kronecker sums of arbitrary matrices using a mu-mode approach. It is based on the approximation of the integral representation of the phi-functions by Gaussian quadrature formulas, combined with a scaling and squaring technique. For more details, see the reference manuscript
This GitHub repository contains a MATLAB implementation of PHIKS (fully compatible with GNU Octave) and all the functions and the scripts needed to reproduce the numerical examples contained in the reference manuscript.