Best workflow to convert H2(g) partial pressure / mole fraction into absolute moles

Hi everyone,

I am currently post-processing results from a reactive transport simulation and need to convert 

𝑯𝟐

(

𝒈

)

 from its reported mole fraction or partial pressure into absolute moles for specific grid blocks over time.

Looking at my visualization files (gas.tec and time.out), I see the individual gas partial pressures (e.g., P_H2(g)). To get the true number of moles, my current planned approach is to calculate total gas moles (

ntotal

𝑛total) per grid block using the Ideal Gas Law:

1. Calculate total gas phase pressure (

   Ptotal

   𝑷𝐭𝐨𝐭𝐚𝐥) by summing all gas partial pressure columns or pulling total pressure from flow.out (where gas saturation 

   𝑆𝑔

   >0

   ).
2. Calculate gas volume (

   Vgas

   𝑽𝐠𝐚𝐬) using: 

   𝑉block

   ×porosity

   ×𝑆𝑔

   .
3. Solve for total moles (

   𝑛total

   =𝑃total⋅𝑉gas𝑅⋅𝑇

   ) and multiply by the 

   H2

   𝐻2 mole fraction (

   𝜒𝐻2

   =𝑃𝐻2𝑃total

   ).

Before I build a Python/Pandas parser to automate this over my output datasets, I wanted to check with the community:

•  
• Is this the correct physical workflow for extracting molar quantities from TOUGHREACT outputs, or is there an edge case regarding phase volumes or compressibility factors I am overlooking?
• Is there an alternative output flag or pre-existing post-processing tool/script (like toughio) that can output mass or moles of a gaseous species directly?
•  

Any advice, formulas, or example script snippets would be greatly appreciated!

Thanks in advance,