<p>The study of Markov-modulated fluid models with upward jumps and phase transitions shows that the joint distribution is governed by a non-homogeneous linear differential system with specific boundary conditions. Spectral analysis technique was used to solve this system in which the unique solution is expressed by components that are calculated through the resolution of a linear system. The purpose of this paper is to propose a novel methodology that computes numerically the moments for these models. Our approach is recursive: the moment is obtained from the preceding moment via a linear system which in turn is expressed via the components derived for the joint distribution of the fluid level and the modulating process. Numerical illustrations are presented and analyzed.</p>