In this study, we propose a phase-field-based finite element model to simulate two-phase ferrofluid flows in two and three dimensions. The proposed model combines the Cahn-Hilliard equation to handle the phase field, the Poisson equation to account for magnetics, and the Navier-Stokes equation to characterize fluid flow. To efficiently handle this coupling, we present a linear, totally decoupled numerical scheme, which involves solving four separate equations independently, namely, a linear elliptic system for the phase function, a Poisson equation for the magnetic potential, a linear elliptic equation for the velocity, and a Poisson equation for the pressure. To assess the accuracy, applicability, and numerical stability of the model, we conduct simulations for several typical problems. These include investigating the deformation of a ferrofluid droplet under a two-dimensional uniform magnetic field model, the bubble coalescence in ferrofluids under a three-dimensional uniform magnetic field model, the collision of two ferrofluid droplets under two-dimensional shear flow, and the two-dimensional interfacial instability of a ferrofluid. The numerical results confirm the model's capability to robustly simulate multiphase flow problems involving high-density and high-viscosity ratios, both in two- and three-dimensional problems. Moreover, the model effectively captures fundamental phenomenological features of two-phase ferrofluid flows under large topological changes such as the Rosensweig instability.