  In the era of advanced electromagnetic and gravitational wave detectors, it
has become increasingly important to effectively combine and study the impact
of stellar evolution on binaries and dynamical systems of stars. Systematic
studies dedicated to exploring uncertain parameters in stellar evolution are
required to account for the recent observations of the stellar populations. We
present a new approach to the commonly used Single-Star Evolution (SSE) fitting
formulae, one that is more adaptable: Method of Interpolation for Single Star
Evolution (METISSE). It makes use of interpolation between sets of pre-computed
stellar tracks to approximate evolution parameters for a population of stars.
We have used METISSE with detailed stellar tracks computed by the Modules for
Experiments in Stellar Astrophysics (MESA), Bonn Evolutionary Code (BEC) and
Cambridge STARS code. METISSE better reproduces stellar tracks computed using
the STARS code compared to SSE, and is on average three times faster. Using
stellar tracks computed with MESA and BEC, we apply METISSE to explore the
differences in the remnant masses, the maximum radial expansion, and the
main-sequence lifetime of massive stars. We find that different physical
ingredients used in the evolution of stars, such as the treatment of radiation
dominated envelopes, can impact their evolutionary outcome. For stars in the
mass range 9 to 100 M$_\odot$, the predictions of remnant masses can vary by up
to 20 M$_\odot$, while the maximum radial expansion achieved by a star can
differ by an order of magnitude between different stellar models.
