Best approach for adding chemistry to an existing THM CO2 injection model (THMC)?
Hi all,
I’m building a CO2 injection / geological storage model for my final year project, using OGS 6.5.7. I have a working THERMO_HYDRO_MECHANICS model (Volve field dataset, tetrahedral mesh, three permeability zones) and now want to extend it to a full THMC model with reactive transport (calcite dissolution driven by injected CO2, tracking C(4)/Ca/H via PHREEQC).
From reading the docs, changelogs, and forum, my understanding so far is:
The native <global_process_coupling>/staggered scheme is for staggering sub-equations within one process (e.g. pressure vs. concentration inside ComponentTransport, or H–M inside HM_PHASE_FIELD) — not for combining two different process types like THERMO_HYDRO_MECHANICS and ComponentTransport in one project file.
ComponentTransport has <chemically_induced_porosity_change> + a KozenyCarman permeability model, which handles the chemistry→porosity/permeability feedback natively — but this appears to be scoped to ComponentTransport only, with no equivalent hook in HYDRO_MECHANICS/THERMO_HYDRO_MECHANICS.
So my current plan is to run the THM and reactive-transport models as two separate .prj files on the same mesh, staggered externally via a Python driver, exchanging pressure/temperature (THM → C) and porosity/permeability (C → THM) between coupling intervals.
Questions:
Is this external two-file staggering actually the recommended approach for THMC with the current OGS version, or is there a more direct mechanism I’m missing?
Is there existing tooling (ogs6py, ogstools, or something in the benchmark suite) built for exactly this kind of two-process external staggering, rather than me hand-rolling the driver?
Has anyone done chemistry↔mechanics feedback in OGS before (porosity/permeability changes from calcite dissolution affecting a mechanical/stress field), and if so, what did that coupling loop look like in practice?
Happy to share my .prj files if useful. Thanks for any pointers!
Hi,
you are right, that for now we don’t have a fully coupled THMC process in OGS. However, this is on our agenda for the near future.
As you mentioned, you could try the workaround by exchanging data between two separate process runs. We did that recently with one of our models having only a one-way coupling THM->C. We had a full THM calculation and we took the timesteps by parsing the pvd file and constructing a time dependent heterogeneous parameter in the component transport project file which was used as Dirichlet condition for the pressure field. In theory, similar things could be done to achieve a two-way coupling with some kind of a Python interface (like FINAM FINAM · GitHub ).
There are some caveats, I think parts of the mass balance equation are plugged-in in the transport equation so there might be some inconsistencies, as the ComponetTransport process assumes a different mass balance than THM does. So you should check the velocity field carefully.
So in theory, your approach should work. But I can’t give you a guarantee that everything will be consistent at the end. If this matters really matters you need to wait until we have a working prototype. You are free to share your experiences and I can share the MR once it is in a working state.
One practical detail for the two-way loop: feed back porosity/permeability at the same time level as the pressure/temperature you hand to ComponentTransport, otherwise the staggering is effectively explicit and the two fields can drift apart. If the Kozeny–Carman feedback is strong, under-relax the porosity update or do a couple of inner iterations per coupling interval, since porosity→permeability→velocity→dissolution is a stiff loop. Given Jörg’s mass-balance caveat, a cheap check is to run the transport model passively (no chemistry) on the THM flux field and compare velocities before adding calcite dissolution.