Model insensetive when running pflotran from 2000 yr checkpoint file

[Bisesh Joshi]

Hi all,

I'm trying to reach an equilibrium state for mineralization of SOM that produces NH4+ and DOC at first, and then run a suite of biogeochemical reactions. Model results seem to be constant at around 2000 years. However,  if I use that point as a starting point, the model is insensitive to any changes in parameters such as logK values, mineral rate constant, SOM volume fraction  and permeability of soil layers. Also, I checked SOM rates and VF and they are zero, however, SI values are non-zero. 

The problem is I want to increase NH4+ concentration in the system to calibrate it to observed values in the middle of 2-D transect. But I'm unable to get it as the model is insensitive now. Also, I'm not sure if I can use this approach for scenario analysis, as the model is insensitive. 

Please let me know if anyone has any clue how to solve this problem.

Thank you.

Best regards,

Bisesh

--

Bisesh Joshi

PhD student

Water Science and Policy

Graduate College

University of Delaware

901-297-6761

[Hammond, Glenn E]

Bisesh,

Please send a brief description of the key mass action equations (A + B -> C) and mathematical rate expressions (rate = k[A][B]) so that we can better understand the conceptual model. Or you could send a SIMPLE input file. It sounds like a mineral should be formed through transition state theory, but it is not.

Glenn

From: pflotra...@googlegroups.com <pflotra...@googlegroups.com> on behalf of Bisesh Joshi <bjo...@udel.edu>
Date: Monday, October 20, 2025 at 6:20 AM
To: pflotra...@googlegroups.com <pflotra...@googlegroups.com>
Subject: [pflotran-users: 8596] Model insensetive when running pflotran from 2000 yr checkpoint file

Check twice before you click! This email originated from outside PNNL.

--
You received this message because you are subscribed to the Google Groups "pflotran-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to pflotran-user...@googlegroups.com.
To view this discussion visit https://groups.google.com/d/msgid/pflotran-users/CAEvi1yQK0Kqu%2BvM90dXGa%3DTgCrP43GAXyn2FkkVnMEZ7gmf4zw%40mail.gmail.com.

[Bisesh Joshi]

Hi Glenn, 

Thank you for your response.

Here are my simple input file with its other input blocks.

I look forward to hear from you soon.

Best regards,

Bisesh

[Bisesh Joshi]

Hi all,

Also, a follow-up question: Model spinning up for 2000 years to reach equilibrium is consuming all the mineral mass. Is there a way to reach an equilibrium before that, or give a constant mineral source instead of a limited volume fraction (except the source-sink method)?

Thanks in advance!

Best regards,

Bisesh

To view this discussion visit https://groups.google.com/d/msgid/pflotran-users/2c285b76-79ea-4475-9abc-4ff527992f7en%40googlegroups.com.

[Hammond, Glenn E]

 Bisesh,

As you may have already assessed, the insensitivity of your spun-up runs to log K, volume fraction (VF), and other variables is primarily a result of the absence of the targeted mineral along with the aqueous solution being undersaturated concerning this mineral, preventing any precipitation from occurring.

The primary challenge lies in the dynamic nature of your conceptual model. It is uncertain whether a steady-state, spun-up solution can be established under these circumstances. Were I to approach this problem myself, I would begin by simplifying the chemical framework to focus on a single set of soil organic matter (SOM). This would allow for a clearer pathway to appropriately spinning up the model. Once this simpler model is successfully stabilized, I would then proceed to introduce greater chemical complexity incrementally.

In addition, one method to achieve swift "equilibration" between aqueous chemistry and associated minerals during the spin-up process is to use larger reaction rate constants initially. Once the simulation is restarted, these rate constants can be reduced for more realistic scenarios, provided that the names and numbers of species within the system remain consistent across simulations. Notably, you can implement this adjustment using the following undocumented keyword:

CHEMISTRY 

    … 

    MINERAL_RATE_SCALING_FACTOR <float> 

END

This keyword facilitates global scaling of mineral reaction rate constants without necessitating individual modifications for each rate value.

Regards,

Glenn

To view this discussion visit https://groups.google.com/d/msgid/pflotran-users/CAEvi1yRPE6b-fKX37DrY_%2BVGw%3DJKZ1KJbcHo9NkxZjofEgvNwQ%40mail.gmail.com.

[Bisesh Joshi]

Hi Glenn,

Thanks for the suggestion. I will try that.

The mineral volume fraction is fully depleted  (i.e., it is changed from 0.075 to 0.0) with logK that I'm using (-31.01). However, lower logK values produce less NH4+ than its intended concentration.

Is there any keyword or code block for a constant mineral source in the system?

Also, can new minerals be introduced to the spun-up model? For example, I initially spin up a model for 2000 years with SOM in the system, but now I want to run different scenarios on the model using other minerals like pyrite, calcite etc?

Thanks again!

Best regards,

Bisesh

To view this discussion visit https://groups.google.com/d/msgid/pflotran-users/PH0PR09MB783683D2DF4259076E7900829AFCA%40PH0PR09MB7836.namprd09.prod.outlook.com.

[Hammond, Glenn E]

From: pflotra...@googlegroups.com <pflotra...@googlegroups.com> on behalf of Bisesh Joshi <bjo...@udel.edu>

Date: Tuesday, October 28, 2025 at 6:17 AM
To: pflotra...@googlegroups.com <pflotra...@googlegroups.com>
Subject: Re: [pflotran-users: 8602] Model insensetive when running pflotran from 2000 yr checkpoint file

Hi Glenn,

Thanks for the suggestion. I will try that.

The mineral volume fraction is fully depleted  (i.e., it is changed from 0.075 to 0.0) with logK that I'm using (-31.01). However, lower logK values produce less NH4+ than its intended concentration.

I do not recommend modifying logKs in an attempt to preserve equilibrium in a dynamic system.

Is there any keyword or code block for a constant mineral source in the system?

Try adding OVERRIDE_MINERAL_MASS_ACTION to the CHEMISTRY block and set the stoichiometry of the mineral to  zero. E.g.,

CHEMISTRY

  ...

  OVERRIDE_MINERAL_MASS_ACTION

    Calcite

      REACTION 0. Calcite + H+ <-> Ca++ + 1. HCO3-

      LOGK 2.2257 1.8487 1.3330 0.7743 0.0999 -0.5838 -1.3262 -2.2154

    /

  /

END

https://documentation.pflotran.org/user_guide/cards/subsurface/chemistry/override_mineral_mass_action_card.html#override-mineral-mass-action-card

My hope is that this trick will work. But note that you will have mass imbalance in your scenario as mineral will precipitate/dissolve with no change in mineral volume fraction. I also note that the code may crash (I have not tested zero stoichiometries).

Also, can new minerals be introduced to the spun-up model? For example, I initially spin up a model for 2000 years with SOM in the system, but now I want to run different scenarios on the model using other minerals like pyrite, calcite etc?

You can overwrite any concentration or mineral volume fraction in the restart file if you chose HDF5 (.h5) as the format. I recommend using HDFView to view the values and the Python package h5py to modify the dataset within the .h5 file. There are example scripts in PFLOTRAN_DIR/src/python/* that make similar modfications to .h5 files.

Glenn

To view this discussion visit https://groups.google.com/d/msgid/pflotran-users/CAEvi1yTd54MDwJV_JKZsYLGAfU4rw70ETZNSaWeLCcfrQo3hyA%40mail.gmail.com.

[Bisesh Joshi]

Hi Glenn,

Thanks for your timely reply and suggestions.

It seems the Pflotran version installed on my local machine does not support it?

Are there any other alternatives?

Best regards,

Bisesh

You received this message because you are subscribed to a topic in the Google Groups "pflotran-users" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/pflotran-users/fa1CBDZc74A/unsubscribe.
To unsubscribe from this group and all its topics, send an email to pflotran-user...@googlegroups.com.
To view this discussion visit https://groups.google.com/d/msgid/pflotran-users/PH0PR09MB78365D9155A528F020D838F79AFDA%40PH0PR09MB7836.namprd09.prod.outlook.com.

[Hammond, Glenn E]

Update your version of PFLOTRAN.

Glenn

To view this discussion visit https://groups.google.com/d/msgid/pflotran-users/CAEvi1ySTXmLih7ptvHH4oCDN7qmpGUqTjTQqu7zqb9bWXLSZNw%40mail.gmail.com.

[Bisesh Joshi]

Hi Glenn,

Thanks for the suggestion

I have rebuilt Pflotran and tried OVERRIDE_MINERAL_MASS_ACTION. However,  it needs a unit value there:

ERROR: Non-unity reactive mineral stoichiometry in mass action override (i.e., \nu_m /= -1.): -0.000000E+00

Probably, I will have to work around it please let mw know if you have any other suggestions.

Best regards,

Bisesh 

To view this discussion visit https://groups.google.com/d/msgid/pflotran-users/PH0PR09MB7836F1BFE9C9064279019C519AFAA%40PH0PR09MB7836.namprd09.prod.outlook.com.

[Hammond, Glenn E]

Looks like I cannot outsmart myself.  This may be ill-advised, but try commenting out line 2227 at

https://bitbucket.org/pflotran/pflotran/src/723c4c72fc68c19111ac3e0f92646817039629d2/src/pflotran/reaction_database.F90#lines-2227

Glenn

To view this discussion visit https://groups.google.com/d/msgid/pflotran-users/CAEvi1ySBskmhvXArch0cN0JSbopFJqWFk38w9hSmkqOUOaDyQw%40mail.gmail.com.

[Bisesh Joshi]

Hi Glenn,

Thanks for the reply.

Unfortunately, this also did not work. 

Also, is there a Python script to extract variables from Tec block I can work on? Outputs seem to be saved in the Tecblock in the updated PFLOTRAN version, even though I specify to save output as Tec point. My script for Tec point version does not work for this.

Best regards,

Bisesh 

To view this discussion visit https://groups.google.com/d/msgid/pflotran-users/PH0PR09MB7836518B7C8D21CC250FE8299AFBA%40PH0PR09MB7836.namprd09.prod.outlook.com.

[Bisesh Joshi]

Hi Glenn, 

I hope you're doing well.

I'm able to develop the script for reading the Tecplot block outputs.

Thank you!

[Hammond, Glenn E]

The Tecplot POINT format is supported exclusively in serial (single-process) because generating such files in parallel can require extensive inter-process communication and impose a significant memory burden on process 0.

Also, I recommend extracting information from the .h5 files using the h5py Python package as the data is binary (full 64-bit precision) and the access is straightforward, though I admit that the formatting of the .h5 group time stamps is a bit tricky.

Glenn

To view this discussion visit https://groups.google.com/d/msgid/pflotran-users/CAEvi1yQ%2BNb-9ot%3Dmf7Y1oeSYdqhbBAcokSvKJzfVY8LxF21LDg%40mail.gmail.com.

[Bisesh Joshi]

Thanks for the timely reply, Glenn. 

 I'm trying to simulate the kinetic reaction of the Ferrihydrite/Goethite mineral with SOM  to produce Fe2+ + CO2 +NH4+ + H2O. But, since the log K value for iron mineral dissolution is low, it is constraining the reaction outputs.

 Is there a way to simulate this kinetically controlled  (without logK) biotic reaction? Also, can I simulate both reactions (kinetic and log K controlled ) together?

Best regards,

Bisesh

To view this discussion visit https://groups.google.com/d/msgid/pflotran-users/PH0PR09MB7836EDBF99C208CB35877E8B9AC7A%40PH0PR09MB7836.namprd09.prod.outlook.com.