Discrete Element Method (DEM) Simulation of Coupled Thermal, Mechanical and Melt Dynamics during Formation of Plumbing Systems of Large Igneous Provinces

Linking small-scale fracture processes to lithosphere-scale magma transport remains a core challenge in understanding the development of magmatic plumbing systems in Large Igneous Provinces (LIPs). In this study, we employ a two-dimensional Discrete Element Method (DEM) to investigate the coupled thermo-hydro-mechanical evolution of plumbing systems in the continental lithosphere. Using the MatDEM framework, we simulate fracture propagation, magma migration, and heat transfer from a magma chamber located at the lithosphere-asthenosphere boundary to the upper crust. Magma transport is modeled through a pore density flow approach, allowing dynamic coupling between pore pressure, temperature, and mechanical deformation of the host rocks. Scaling principles are applied to ensure mechanical and thermal similarity between numerical models and natural systems. The initial model shows that magma overpressure and thermal expansion generate radial fractures around the magma chamber, which progressively evolve into vertically connected magma pathways (i.e., dikes). We systematically examine the influence of layering structure, pre-existing faults, lower crustal strength, crustal thickness variations, magma viscosity, and magma overpressure on plumbing system development. The existence of horizontal weak zones or mechanical boundaries, such as the Moho and intra-crustal compositional boundaries will promote sill emplacement along these horizontal boundaries prior to renewed upward magma propagation. Steeply dipping faults further localize magma ascent and control geometry and number of sub-vertical conduits. A mechanically strong lower crust acts as a barrier to vertical magma ascent, favoring magma underplating and prolonged magma storage near the Moho. Crustal thickness gradients will drive magma migration toward the thinner crust. Increasing magma viscosity reduces magma flowability and limits the extent of fracture-controlled magma networks, whereas higher magma overpressure enhances fracture opening and results in a plumbing system with wider conduit width and larger spatial distribution. Our results fit well with geological and geophysical observations of LIPs. This DEM-based approach provides a bridge between small-scale fracture processes and the large-scale magma transport and emplacement in LIPs.