In this work, we propose the Prony fitting decomposition (PFD) as an accurate and efficient exponential series method, applicable to arbitrary interacting bath correlation functions. The resulting hierarchical equations of motion (HEOM) formalism is greatly optimized, especially in extremely low temperature regimes that would be inaccessible with other methods. For demonstration, we calibrate the present PFD against the celebrated Pad\'e spectrum decomposition method, followed by converged HEOM evaluations on the single-impurity Anderson model system.
Comment: 4 pages, 3 figures