Magma migration is a key mechanism governing a wide range of geological and thermal phenomena, such as dyke propagation, volcanic activity, the genesis of ore bodies, and the development of geothermal systems. Because magma movement cannot be directly observed in most geological settings, our understanding of these processes remains constrained. As a result, constructing a robust numerical framework that can properly bridge the relevant spatial and temporal scales is essential for realistically modeling magma migration within lithospheric rocks [1]. This work is grounded in the principles of flow mechanics in deformable porous media [2] and introduces a new multiphase computational framework to simulate magma migration within lithospheric rocks. The interaction of the solid rock and fluid magma is represented by the phase-averaging technique. The mechanical behavior of the lithosphere is governed by a visco-elasto-viscoplastic constitutive law that includes both strain-softening effects and permeability enhancement driven by plastic deformation, allowing the progressive evolution of the rock to be modeled realistically. An Eulerian–Lagrangian discretization scheme [3] is used to discretize the resulting coupled and nonlinear set of equations. A series of numerical simulations is performed to demonstrate the ability of the proposed method to model magma percolation through porous lithospheric rocks. Validation of the results is achieved through comparison with well-known benchmark problems.