the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
| 31 Jul 2026
𝙿𝚘𝚛𝚘𝚞𝚜𝙻𝚊𝚋 v1.1: a modular finite element framework for simulation of fractured porous media
Abstract. Numerical modeling of fractured porous media requires the combined treatment of coupled processes, nonlinear geomechanical behavior, and strong discontinuities such as fractures and faults. Although several open-source simulators are available for porous-media applications, many are oriented toward large-scale production analyses and rely on mixed discretizations or complex software architectures. This paper presents 𝙿𝚘𝚛𝚘𝚞𝚜𝙻𝚊𝚋, an open-source, object-oriented finite-element framework implemented in MATLAB for transparent formulation development, verification, and rapid prototyping of coupled porous-media models. The current release supports mechanical, single-phase flow, two-phase flow, hydro-mechanical single-phase, and hydro-mechanical two-phase formulations. Nonlinear geomechanical analyses are supported through different constitutive models and solution procedures, including quasi-static continuation schemes and nonlinear transient analyses. Pre-existing fractures and faults can be represented independently of the background mesh using recent embedded finite-element formulations, which simplify model generation by avoiding the need for mesh conformity with discontinuity geometries. The framework stands out for integrating these features within a single FEM-based discretization setting. Furthermore, the code architecture is developed around abstract OOP classes designed to support modularity and extensibility. Verification examples covering the implemented features are presented and compared with reference solutions.
- Preprint
(2009 KB) - Metadata XML
-
Supplement
(1069 KB) - BibTeX
- EndNote
Status: open (until 25 Sep 2026)
-
RC1: 'Comment on egusphere-2026-3322', Anonymous Referee #1, 11 Sep 2026
reply
Review of egusphere-2026-3322 (PorousLab 1.1)## SummaryPorousLab 1.1 is an MIT-licensed MATLAB finite element framework for two-dimensionalcoupled hydro-mechanical problems with single- and two-phase flow, nonlinear geomechanics,continuation methods and embedded strong discontinuities (E-FEM). The formulation isstandard and, as far as I checked (balance equations, weak forms, residuals) seems correct;the code is archived with example scripts and regression tests. Thecontribution is in integration and accessibility rather than in new methods: the E-FEMformulations and the two-phase treatment come from published work. Two aspects requiresubstantial revision. First, the comparison with existing frameworks is not appropriate/balanced:design choices of other codes are presented as deficiencies, in the case of PorePyI'm quite certain it is misclassified, and the dependenceon a proprietary runtime is not acknowledged.Second, the verification is qualitative throughout:no error measures or refinement studies, no fracture-network casewith a reference solution (these exist elsewhere, see below),no blocking fractures, and no example in which the embedded faultslips or opens (although fault mechanics is a main application). Two examples cannot bereproduced from the paper alone (missing material and fracture parameters).I recommend major revisions.I'm all for new frameworks, there is always something new to learn!The requests are about demonstrating the framework fairly and quantitatively, not the method itself.## General comments### Framing of the comparison with existing frameworks (Sect. 1, Table 1) is not appropriate:The introduction motivates PorousLab largely by presenting properties of other frameworks asdisadvantages: hybrid discretizations (l. 60–62, 68–69, 72–73), compiled languages andperformance orientation (l. 65–68). Most of these are deliberate design choices rather thandeficiencies so using a negative connotation is inappropriate.Also the claim are actually not substantiated, and several statements do not hold as written.(i) Example of a mistake: PorePy is classified as "Hybrid" in Table 1, and l. 72–73 state that its coupledformulation with mechanics relies on more than one discretization strategy. PorePy discretizesflow and mechanics with cell-centred finite volume schemes throughout (MPFA for flow, MPSA formechanics; see Keilegavlen et al., 2021), i.e. a single finite-volume setting. It offersRT0 mixed finite elements and a virtual element method for flow in addition, as options andthis is great.(ii) The "PDE scheme" column conflates the method family (FEM/FVM) with the number of schemesa framework offers. Example: in DuMux, geomechanics is handled by a combination ofa vertex-centered FV (for displacements) and a cell-centred TPFA/MPFA FV scheme for flow. Both are finite-volume schemes implemented in one assembly framework, and the user selects one per sub-problem;that different schemes can be combined is a feature of theseframeworks, not a limitation. The claim that hybrid schemes require "proficiency in multiplenumerical methods" (l. 69) is not substantiated in this case they just need to know finite volumesand the software abstractions hide the details.(iii) Interestingly, in this work: Appendix A introduces node-centred cells (Fig. A1),cell-wise saturation and density, and element-wise averaged mobilities in order to obtain local mass conservation. These are clearly control-volume treatments inside the finite element method(and there is a good reason for this). But this means the manuscript's own H2physics is not a pure Galerkin FEM either. The "single FEM discretization" argument(abstract l. 11, l. 75–76, 80–82) should be removed accordingly.I recommend removing the entire "PDE scheme" column. List neutrally what schemesare avaible elsewhere and rewrite l. 63–74 so that the contribution is stated positively (whatPorousLab offers: general constitutive models, continuation methods, E-FEM, a scriptedMATLAB workflow) rather than by discounting other codes.Several authors seem to be core developers ofKratos Multiphysics and of GeMA, whose physicsabstraction is adopted. A paragraph on why a separate MATLAB framework wasdeveloped, and whether PorousLab is intended as a prototyping front-end for those codes,would position the contribution more convincingly than the comparison table.### PorousLab is an open source software on a proprietary platform:PorousLab is released under the MIT licence but requires a MATLAB licence, andit is stated that it is not compatible with GNU Octave. At the same time the manuscript excludesCOMSOL and Abaqus from the comparison as commercial software and builds on"open-source" and "lower entry barrier" arguments (abstract, l. 74, 82–85, 405–407). Theexclusion is a legitimate scoping decision, but the text should say explicitly that theframework is open-source code on a proprietary runtime, and that this touches theaccessibility argument; at least in comparison with PorePy (Python, also easy to use)and MRST, which (as far as I know) actually makes some effort to run on Octave.Please also state which MATLAB featuresprevent Octave use (App Designer, specific classdef features, ...) and whether Octave supportof the non-GUI core is planned; a table note on the "Platform" column would be appropriate.### Verification against published fracture benchmarks missing:All verification cases are one-dimensional, single-fracture, or parallel-fracture problems,and the only network case (Fig. 19e–f, "random fractures") has no reference solution.The authors state that "Figure 19 shows that the pore-pressure profiles obtained withPorousLab agree closely with the reference solutions". Community benchmarks for exactly thissetting exist, with reference results from a range of codes: Flemisch et al. (2018,Adv. Water Resour. 111, 239–258) define four two-dimensional single-phase flow cases(Hydrocoin; a regular fracture network with permeable and blocking variants; a complexnetwork with mixed permeable and blocking fractures; a realistic outcrop-based network of 64fractures), with geometries, meshes and reference solutions publicly available.At least one of these cases should be added. This should (a) verify the hydraulic E-FEM on intersectingfractures rather than on isolated or parallel ones, (b) exercise blocking fractures, whichthe abstract and Sect. 1 mention but no example tests, and (c) place PorousLab on the samecomparison as some of the other codes in Table 1 already went through.The "random fractures" case should either receive areference or be labelled as a show case rather than a verification.### Quantitative verification and cost:"Good agreement" or "close agreement" is stated for every example, but no error measure,no mesh or time-step sizes for the transient cases ("adaptive time stepping" withoutparameters), and no convergence under refinement are reported. For a modeldescription paper, please give for each example the discretization parameters and DOF count,typically some error norm value with respect to the reference, and for at least one example a refinementstudy. In the same spirit, the "limited computational performance" and "small- tomedium-scale" statements should be accompanied by wall-clock times (just to get a rough ideaand be able to compare this to other codes or for someone that develops their own code lateras a reference), and DOF counts of the examples.## Specific comments- Sect. 6.1 probably MPa instead kPa for load?- also Sect. 6.1 text says "a rigid footing applies auniform vertical pressure" --- but a rigid footing imposes a uniform displacement; a uniformpressure would correspond to a flexible footing. Clarify.- Sect. 6.2: material parameters are missing.- Sect. 6.3: geometry, aperture and permeability of the parallel and random fracture sets are not given- the Jacobian in Eq. (64) omits derivatives of K and C, although Kand C depend on x in the two-phase physics (relative permeability, saturation, density);Appendix B includes these derivatives (B15–B17). Either state that Eq. (64) applies when Kand C do not depend on x, or write it consistently with Appendix B.- I think you should cite the specific version on Zenodo https://doi.org/10.5281/zenodo.20612078- As far as I know, PorePy's bulk model is linear elastic,but its fracture treatment is a frictional contact-mechanics formulation solved with asemismooth Newton method, which is more general than the penalty-based cohesive treatmentused here ("contact-penalty parameters"). "More limited in scope" seems to be a mischaracterization.- Appendix A2: Relative permeabilities are arithmetic nodal averages of saturation withoutupwinding. This is adequate for the capillary-driven McWhorter problem but is known to losemonotonicity for advection-dominated displacement; a Buckley–Leverett test would maybe show thebehaviour, otherwise the limitation should be stated.- Local mass conservation: "Guarantee local mass conservation" (l. 639, 657)"the flux terms are standard Galerkin -> please clarify what "local" refers to,and report the global mass-balance error for the transient examples.- the HM discontinuity terms: the aperture is updated from themechanical deformation and enters the cubic law (Eq. (46)) so the fracturetransmissivity depends on the jump DOFs. Appendix B lists the H2M Jacobian as I see it withoutdiscontinuities only. Is the aperture lagged? Is it really the full Jacobian.- Abstract, l. 3: "rely on mixed discretizations" suggests mixed finite elements; use "combine different discretizations" or drop the phrase.- l. 116: "include several continuation strategies" --> "includes".- l. 137: "validated on MATLAB R2022a" --> "tested on".- l. 195: "repectively" --> "respectively".- Fig. 2: the labels n, s, h, g0, t0 are not explained in the caption or text.- Eq. (43): phi_d is not defined; state its meaning (and value) for an open fracture.- l. 305: a Delaunay triangulation is only "conforming" to the interface if it isconstrained or performed per side; rephrase.- Fig. 20 caption: "The fractured is represented" --> fracture- Fig. 23 legend says "Analytical", the text says semi-analytical.- Helmig (1997) is a single-author book, remove "et al."ReplyCitation: https://doi.org/
10.5194/egusphere-2026-3322-RC1 -
CC1: 'Comment on egusphere-2026-3322', Giacomo Medici, 13 Sep 2026
reply
General comments
Very good numerical solution on fluid flow in geological media. The research is also characterized by a good degree of novelty as a consequence of the interdisciplinary angle. You can follow my specific comments to improve the manuscript.
Specific comments
Lines 16-17. “CO2 Storage, flow and transport”. You can insert recent reviews that discuss all these topics:
- Baqer, Yousef, Xiaohui Chen. 2022. A review on reactive transport model and porosity evolution in the porous media." Environmental Science and Pollution Research 29(32), 47873-47901.
- Lupi, F., Agbotui, P. Y., Medici, G. 2026. Hydraulic Conductivity in the Mesozoic units of the Umbria-Marche succession (Italy); insights towards a sustainable management of carbonate aquifers worldwide. Sustainability, 18, 9297; https://doi.org/10.3390/su18189297.
Lines 67-80. When you discuss the link between fluid flow, and geomechanics, I suggest to be more explicit and conceptual. Which are the key aspects that PorePy and PorousLab incorporate? What about the relation between fracture aperture and orientation of the stress tensor in the sub-surface?
Lines 67-80. Improve the conceptual description of the novelty of PorousLab.
Lines 67-80. Is PorousLab capable to represent exchanges matrix-fractures? Not clear in the introduction.
Line 105. Describe the specific objectives of your research by using numbers (e.g., i, ii, and iii) at the end of the introduction.
Lines 614-634. Highlight novelty of the research in the conclusion.
Figures and tables
Table 1. Are those numerical solutions for fluid flow? This is unclear from the caption. The information in the brackets is not sufficient to explain the content of the table.
Table 2. For example, the caption of Table 2 is very clear.
Figure 1. Do you need an approximate spatial scale?
Figure 2. What do you mean by “strong discontinuity”? Please, improve the caption.
Figure 6. Make the figure larger and improve the graphic resolution.
Algorithm 1. The use of bold is unclear.
Figure 15. The boundary conditions are not clear. Do you need adding labels or a legend?
Citation: https://doi.org/10.5194/egusphere-2026-3322-CC1
Model code and software
PorousLab v1.1.0 Danilo Cavalcanti https://doi.org/10.5281/zenodo.15624961
Viewed
| HTML | XML | Total | Supplement | BibTeX | EndNote | |
|---|---|---|---|---|---|---|
| 179 | 60 | 22 | 261 | 41 | 20 | 24 |
- HTML: 179
- PDF: 60
- XML: 22
- Total: 261
- Supplement: 41
- BibTeX: 20
- EndNote: 24
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1
