Isotope simulation stability

[p.m.berger]

Hello,

I am using the attached system to model the dissolution and reprecipitation of pyrite as buffered by organic material (phenol in the model). Along with this, I want to look at the change in 34S. I need to ensure I have enough phenol in the system, but when I increase the amount released into the system from 0.1 mmol/m3/yr to 1.0 the stability of the isotopic part of the model goes haywire. Can you confirm that none of the overall model stability constraints limit step size based on isotopic changes? If my only options are really lowering the Courant and step_increase then I'll probably have to drop the model as it will take too long to run.

Thanks,

Peter

iso_test.x1t

[Jia Wang]

Hello Peter,

The Xstable command controls the convergence criterion for isotopes in the Reactive Transport Modeling apps. It looks like you have already lowered the Xstable value in your model, along with the Courant constraint and the step_increase variable.

If you need to maintain a specific concentration of phenol in the solution, an alternative approach could be to create a fake phenol mineral and set your system in equilibrium with that mineral. You can set the equilibrium constant so that the dissolution reaction maintains a constant level of phenol. I decoupled HCO3-/phenol and added phenol as a basis species, swapping in the fake mineral.

This approach helped with the isotopic stability issue, and I was able to maintain it while running the simulation with a higher Courant number and step_increase value. You can experiment with these variables to reduce the run time

Hope this helps,
Jia Wang
Aqueous Solutions LLC

[Jia Wang]

Hello Peter,

Just want to follow up and make a correction to my previous post. The Xstable controls the stability of the transport equation on the isotopes and not the convergence criterion.

Best regards,
Jia

[p.m.berger]

Hello,

Thanks for looking into it. I actually was doing the solid phase organic when I first ran into problems. I switched to the simple reactant for testing and to simplify my script before posting it. I guess I'll just have to play around with it.

Peter

[Jia Wang]

Hello Peter,

 

You're welcome. Just looking at your input file again, I am thinking perhaps that transport constraint isn't the issue here but more of how the problem is set up with the chemistry. If I reset the transport constraints to their default setting and run the simulation with 1 mmol/m3/yr of phenol, I can turn on the "explain_step" command to see what is limiting the size of the step. I also shortened the simulation length to 2000 years to troubleshoot. Running X1t, it looks like the step size is limited by the Q/K of the reactant of various minerals. 

 

You can change the kinetic rates of your minerals to 0 such that you can test whether or not phenol is actually appearing in your system as expected. Just checking the concentration of your carbonate species, most of the component's concentrations are in the form of CO2(aq), HCO3-, and methane instead of phenol. When I decouple the HCO3-/phenol redox couple, phenol actually began to accumulate in your system at what seems to be a more reasonable concentration. This also helped the program run more smoothly and the issue with the isotopic composition.

 

Hope this helps,

Jia

[p.m.berger]

Hello,

As a follow up to this, I can't decouple the phenol reaction because I am trying to look at pyrite dissolution followed by reprecipitation/reduction by an organic material (phenol being the stand in for said organic material). I have gotten a 1D model to converge in a somewhat reasonable time, however, at the end it says there is no isotopic data to plot for anything, including bulk mineral and fluid, except for pyrite. I assume that just means the isotopic model broke down and xtplot is recognizing that the output is gibberish?

Thanks,

Peter

[Jia Wang]

Hello Peter,

Could you attach the new script so we can take a look?

Best regards,
Jia

[p.m.berger]

Here are the files, unfortunately, I didn't save the output.

Peter

iso2.x2t thermo.com.organ.tdat

[Jia Wang]

Hello Peter,

Thanks for attaching your file. I wasn't able to run your attached file as is in a reasonable amount of time to test. The initial calculation would begin but at the really small delxi, it would barely progress. I reset the delxi value to the default value (0.01) and shortened the length of your simulation significantly to 1 year so I can see the isotopes in the short term.

When I ran it for a short duration, I didn't encounter any output where isotopes didn't have any data. Looking at the isotope result in cross sections, I was able to see some variation in isotope that I think is attributed to the differential flow paths caused by the random porosity distribution.

Is there a case scenario that you can attach that runs in a shorter period of time that demonstrates the no data issue? It would also be really helpful if you can also provide the reduction/precipitation reactions with phenol and pyrite.

Best regards,
Jia

[p.m.berger]

Hello,

This should finish in a reasonable amount of time.

Peter

iso_bug.x1t

[Jia Wang]

Hello Peter,

Thank you for the file. Could you also provide the reaction of pyrite reprecipitation by the reduction of phenol?

Best regards,
Jia

[p.m.berger]

Well, I think the reaction would be:

Phenol(s)  + H2O  + 4 SO4--  + 2 Fe++  = 2 Pyrite  + 2 H+  + 6 HCO3-

Peter

[Jia Wang]

Hello Peter,

Thanks for attaching your file and the reaction. I think one of the issues lies in how the chemistry of your organic carbon is set up in your model. The amount of phenol in solution is so low in comparison to the other forms of carbon in this case, where a majority of the concentration is in the form of inorganic carbon species. If you decouple the Phenol and HCO3- redox pair, you can set a kinetic redox equation for the pyrite precipitation. This way, the amount of phenol can remain stable in your system. I would suggest that you start troubleshooting with the kinetic reaction in a single node system in React. When the chemistry looks correct, then you can move your reaction to the reactive transport apps.

Hope this helps,
Jia