Power system harmonics can increase system loses, excite resonant frequency, and cause equipment damage. It is, therefore, essential to develop accurate harmonic models for harmonic sources, such as power converters. During the last decades, renewable energy sources have become an important part of the worldwide concern with clean power generation. In response to the growing demand for medium and high power trends, multilevel converters (MCs) have been attracting growing considerations. As they become more spread out, their harmonic impact on the system could be significant. In this paper, a fast time-domain method (FTDM) suitable for modeling a MC is proposed. The FTDM essentially models the harmonic of interest as a harmonic state and solves it together with the system differential equations. Since the solution is analytical, the FTDM is very accurate. Nevertheless, to model a MC, the FTDM requires to solve a very large exponential matrix. Hence, a Krylov subspace method is proposed to solve for the exponential matrix, resulted from the FTDM. Two case studies on two different modulations (in-phase disposition and phase-shifted modulations) are performed to validate the proposed method with PSCAD/EMTDC. The computation time of the proposed method is much shorter compared to PSCAD/EMTDC and other solvers.