A Flux-Constrained Generalized Finite Difference Method for Fractured Reservoir Simulation

Abstract Meshless methods based on generalized finite difference formulations (GFDM) offer compelling advantages for reservoir simulation by eliminating global mesh generation and enabling arbitrary point placement. However, standard GFDM implementations fail to enforce local mass conservation. We develop a flux-constrained GFDM (FC-GFDM) that achieves local mass conservation throughconstrained meshless operators without sacrificing geometric flexibility, and extend it to fractured reservoirs via the embedded discrete fracture modeling (EDFM) framework. For each node, the method defines a local support domain and performs Delaunay triangulation within this region to construct control volumes as dual cells connecting circumcenters or centroids of sharing tetrahedra.This local construction requires no global mesh connectivity, enabling arbitrary point addition or removal without re-meshing. The perpendicular bisector (PEBI) nature of circumcenter-based dual cells provides natural orthogonality for efficient two-point flux approximations. Meshless operators are formulated through constrained Moving Least Squares where stencil coefficients simultaneously satisfy Taylor series consistency and discrete flux balance over the control volume. For matrix-fracture interaction, the method employs EDFM to represent fractures as 2D surfaces embedded within the 3D matrix point distribution. The rotation-invariant nature of meshless operators eliminates directional bias in coupling calculations and significantly simplifies handling of arbitrarily oriented and intersecting fractures. Multi-point schemes addressing non-orthogonal geometries and full-tensor anisotropic permeability are presented. Numerical tests are carefully designed to showcase the conservation property and spatial discretiza- tion accuracy across irregular point distributions. The results indicate that our FC-GFDM can be practical for a diverse set of simulation problems.