[How to add both fixed atmospheric CO2 and HCO3- in the system .]

Hello,

You're welcome. It sounds like you should set up a boundary fluid that has the CO2(g) fugacity swapped in. In this case, if you want to set the carbonate component in equilibrium with the CO2(g) in the atmosphere, then this is what you should do for the inlet fluid that you are allowing to enter into your system.

In your Initial pane, you should set the composition of the fluid in your domain. You can see an example of weathering in the Reactive Transport Modeling guide in section 3.8 of the Reactive Transport Modeling Guide. Your case isn't going to be exactly the same as this model, you can see that the CO2(g) is swapped into rainwater in the Fluids pane to set the atmospheric CO2 condition. In this example, CO2(g) fugacity is swapped into the basis and set to a higher partial pressure value to represent a fluid in equilibrium with a soil gas. I recommend playing with the example and parameters to see how the results change.

There are a plethora of other examples of reactive transport models, showing various types of reactions can be set up. In general, you set the initial fluid composition in your domain in the Initial pane. The Fluids pane is where you define the composition of your fluid entering your system and the Interval pane is where you assign them to the inlet.

If you want to troubleshoot the chemistry of the fluid composition, I would recommend doing that in React, a single node system application, before trying to put it in a full reactive transport model. Try to always start simple and build up the complexity in your model.

Hope this helps,
Jia

[Running X1t from Python]

Hello Kush,

You're welcome. You can use the suffix command to set a unique name for your output files. You can look up more information on that in the same Command Reference mentioned above. 

Best regards,
Jia

[Running X1t from Python]

Hello Kush,

You can change the output interval of the plot file by adjusting the dxplot parameter in Config -> Output dialog. For more information on dxplot, please see the GWB Command Reference.

Best,
Jia

[Tips for helping x1t to converge]

Hello Karen,

It is difficult to say what might exactly be the issue here. To troubleshoot, you can turn on "explain steps" in the Stepping dialog and check "Follow output" on the Results pane before running the simulation. This will report the rate limiting factor constraining your reaction and you can investigate further.

You can adjust some of the parameters but I would advise that you familiarize yourself with how they work before doing so. Parameters like the Courant and Xstable are stability criteria that should never be set greater than 1. The epsilon parameter sets the convergence criterion for the newton-raphson iteration, which is how closely solutions in subsequent iterations must match. Other parameters like delxi and dx_init adjusts the maximum step size and the initial step size in the simulation. I would recommend looking in the Command Reference to see how each parameter affects the simulation. I would recommend troubleshooting the input before playing with these parameters. You can see equations for mass transport in section 3.3 in the Reactive Transport Modeling Guide. If you have a copy of the Geochemical and Biogeochemical Reaction Modeling text, you can find details regarding the transport equations in chapter 23 Transport in Flowing Groundwater.

Hope this helps,
Jia Wang
Aqueous Solutions LLC

[Equation syntax for surface area.]

Hello,

You can set the specific surface area to vary with time and many other properties, like volume or mineral mass using a custom equation, script, script file, or even a compiled C++ function . The custom method and accepted syntax for each is the same as that of setting a custom kinetic rate law, which you can find detailed information and examples of in Chapter 5 Custom Rate Laws in the GWB Reaction Modeling Guide. Table 5.1 and 5.2 shows the list of internal parameters and helper functions that can be used in your customized input.

The GWB does not keep track of particle size of minerals but it does keep track of mineral volume as precipitation and dissolution occurs. You won't be able to use the particle size minerals directly in the equation. The program takes in a mineral mass or volume and a total surface area is calculated based on the specific surface area (area/g). If you have a measurement of particle size, you can plug it in the equation and set the SSA as a constant value for the simulation. If you also have a correlation between the particle size and one of the other properties (e.g. mineral volume), you can use one of the custom options mentioned above.

Hope this helps,
Jia Wang
Aqueous Solutions LLC

[Eh pH Diagrams with limits of no species formation]

You're welcome Scott. I want to follow up with some additional information about Act2 and Phase2.

In Act2, boundaries between aqueous species represent equal activity (or if stoichiometry differs, this is accounted for). In other words, they display the predominance of aqueous species. Aqueous species never “do not form at all”. Minerals fields show stability and outside their fields, they do not form. You can play around with equilibrium reaction equations in Rxn to consider species activity ratios besides the standard 50:50. Although the calculation is different, if you want more insight, you can overlay additional data on a P2plot diagram (e.g. using contours above), or make a speciation diagram in React (or cross-section in Phase2/P2plot).

To be precise, P2plot similarly accounts for stoichiometry as necessary when displaying species predominance. In a Predominance map, minerals are treated the same as aqueous species.

Best,
Jia

[retrieval of fugacity using python]

You're welcome. I am glad that helped.

Best,
Jia

[Internal error: transform_comp]

Hello Peter,

You can get this error if different sets of species are loaded in linked ChemPlugin instances. You can check if you are working from the same or name-compatible thermo datasets in each ChemPlugin instance and use the identical “span” commands to cause each to load the same species set.

Hope this helps,
Jia Wang
Aqueous Solutions LLC

[How to add both fixed atmospheric CO2 and HCO3- in the system .]

Hello,

You can't use two parameters, in this case the CO2(g) fugacity and a bulk dissolve concentration,  to constrain one component. If you want to use the CO2(g) fugacity to set the concentration of your dissolved carbonate component, then you would want to swap it in. I would suggest that you double check whether or not that assumption reflects the equilibrium state of your situation. If you have a measured bulk concentration for the carbon component (like in a lab analysis), then it is completely valid to just use your measured value. You have to decide what describes the condition of your system accurately.

Just to clarify, when you add a concentration like mg/kg for a component, like HCO3-, this reflects the total concentration of the dissolved carbonate component. This includes the mass in HCO3-, CO2(aq), CO3--, NaCO3-. NaHCO3, etc. You can specify the free concentration of a specific species, like CO2(aq) or HCO3-, by setting the unit to 'free'. The program calculates an additional amount dissolved in equilibrium with the species concentration. For more information on free vs. bulk concentrations, please see section

Fixing gas in the Reactants pane will hold the value you supplied in the Initial pane, not in the Fluids pane. The Fluids pane is used to constrain the fluid at the boundary, to enter your system.

In terms of charge balancing errors, they are difficult to narrow down without knowing what your system looks like. It could be that you need to choose a different charge balancing ion. You typically want to choose the ion with the least certainty and the highest concentration. Other reasons can also be whether your choice of constraints accurately reflect your fluid. If you know which fluid composition is giving you errors, I would suggest that you start troubleshooting in React and see if it calculates the equilibrium state successfully. Once that is working, go back to your reactive transport modeling.

Hope this helps,
Jia Wang
Aqueous Solutions LLC

[Eh pH Diagrams with limits of no species formation]

Hello Scott,

Are you referring to the diagram generated in Phase2 or Act2?

In Act2, the boundaries represent the equilibrium reaction between species in the adjacent fields. This cannot be altered in Act2.

In a Phase2 predominance diagram, the boundary denotes nodes in the diagram where the concentration in adjacent fields are equal. This cannot be altered. However, you can add a variable to create a 2D map on top of a Predominance map. For example, you can add a contour map of mineral saturation for a specific mineral or the concentration of a species. You can see an example of contours in section 7.9 Example: Mineral solubility of the GWB Reaction Modeling Guide.

Hope this helps,
Jia Wang
Aqueous Solutions LLC

[retrieval of fugacity using python]

Hello Helge,

The results command always return an array, even if there is only one value. To retrieve the fugacity value from the array, you can use the index operator ("[0]") like you have for the pH and concentrations.

Best,
Jia

[Calculating pH of rain water from pCO2]

Hello Jaxon,

You can turn off charge balancing, either by selecting this option from the units in the GUI or in the command line using "balance off". You can perform a "go initial" run (Run menu -> Go initial) and look under Chemical Parameters in Gtplot or in the header block of the text output file for charge imbalance. The charge imbalance will inform you whether your solution is in excess of positive or negative charge. Please see chapter on Using Gtplot in the GWB Reaction Modeling User Guide if you are having issues configuring plot results.

Also note that the fugacity of SO4-- is not fixed, SO4-- is basis species and not a gas. You have swapped the HS- component with H2S(g) in your basis and then fixed the fugacity of H2S(g) in the Reactants pane.

In general, you should try to simplify your input file for troubleshooting. When I performed a "go initial" run, your input file will not converge. This typically occurs when you have specified constraints or composition that is not possible for the equilibrium calculation. I would suggest that you check the initial composition of fluid first. In my previous post, I had suggested that you input the composition of your fluid, not with CO2(g) swapped in, and then use the sliding fugacity path to arrive at an end CO2(g) fugacity desired. Do you have a pH value for your hypothetical fluid? Check that the temperature parameter and the decoupling of the HS-/SO4-- couple makes sense for your fluid.

Hope this helps,
Jia

[Calculating pH of rain water from pCO2]

Hello Jaxon,

If you swap in a fugacity value for a component, the program will use that constraint as a calculation for the initial pH concentration. From your description above, it sounds like you would like to use the sliding fugacity path to evolve the initial fluid to equilibrate at a final end point CO2 fugacity. In this case, you need to set the initial concentration of your fluid in the basis, not swap in a CO2(g) fugacity. If you have all the components needed to calculate the CO2(g) fugacity (i.e., HCO3-, H+, H2O), the program will calculate the CO2(g) fugacity as part of its equilibrium calculation. You can set up a sliding fugacity path for CO2(g) to arrive at the end result of a fluid in equilibrium with 0.003. You can pick up the end result as the starting point of a new run in React. The pickup function allows you to pick up just the fluid or the entire system, including precipitated minerals. For more information on pickup, please see the GWB Reference Manual.

Aside from that, I would consider choosing another ion for charge balancing, since your dissolved carbonate component would change according to the CO2 fugacity. Typically, you would want to choose the ion that has a high abundance and low uncertainty. For example, React is set up to use Cl- as the default ion for charge balancing, given its typically abundant presence and the ability to calculate it from charge balance based on other species concentrations measured in chemical analyses. You can also turn off the charge balance and run your calculation to see if your solution is in excess of cations or anions.

If you would like further help with troubleshooting, please provide more details regarding the assumption that you are making in the initial input. For example, you are fixing the activity of SO4-- and I am not sure why.

Hope this helps,
Jia Wang
Aqueous Solutions LLC

[Surface complexation mass balance issue.]

Hello,

By default, the program assumes the provided concentration describes the mass of dissolved solutes in fluid. At the start of the calculation, X1t calculates an additional mass adsorbed to the surface of the mineral that is in equilibrium with the solution. This is what you observe in 3metals.x1t.

Alternatively, you can set the concentration in your Initial pane as the total dissolved solute in solution and sorbed on the mineral surface. If you would like to set this constraint, go to 'Config' menu -> 'Options...' -> and check 'sorbate include'.

See Reference Manual for more details on "sorbate include". To see an example of two-layer surface complexation, see chapter 7.7.1 Surface Complexation in the GWB Essentials Guide.

Hope this helps,
Jia Wang
Aqueous Solutions LLC

[transition elements speciation in GSS file]

Hello Stanislav,

If you want the program to calculate the distribution of bulk mass across the redox states then you should pick the basis species and provide the appropriate parameters for the oxidation state.

For example, in thermo.tdat, Fe++ is the basis species and the Fe+++ is a redox species. If you only have a bulk concentration for iron, you would select Fe++. If you also provide a parameter for the oxidation state (like O2(aq) or Eh), then when you run a speciation calculation, the total concentration is going to be distributed between Fe++ and Fe+++ species.

Using the above example, if you add Fe+++ to the spreadsheet, the program will only distribute mass to just Fe+++ species during speciation. The reason for this is because the equilibrium reaction between Fe++ and Fe+++ is disabled or as we call it "decoupled" throughout the user guides. This occurs when you add a redox species to the spreadsheet.

You can open the thermodynamic dataset in the TEdit application to view all the entries available in the basis and redox species sections. If you are in a GWB application, you can always launch the thermo dataset loaded by going to "File" -> "View" and select the file that ends with the extension .tdat.

If you already have your data entered in GSS, you can launch a sample to SpecE8 to troubleshoot more closely. If you have added any redox species to your spreadsheet, the redox reaction should be automatically disabled when you go to "Config"-> "Redox Couples...".

The oxidation set by the pe, Eh, or O2(aq) concentration is going to affect all redox couples mass distribution unless they have been disabled. So if you have measurements for the separate pools of mass for Fe++ and Fe+++ species, that oxidation state parameters won't affect that pair's mass distribution.

Section 7.3 Redox disequilibrium in the Essentials Guide walks through a detailed example on decoupling redox reactions. For information on launching samples from GSS, please see section 3.5 Launching SpecE8 and React of the same user guide.

Hope this helps,
Jia

[Isotope simulation stability]

Hello Peter,

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

Best regards,
Jia

[Running X1t from Python]

Hello James,

I am glad to hear that you are making progress with your script and the GWB Plugin. The "result" command always returns values in an array, even when requesting a single value. If the "result" command fails for any reason, such as the requested data doesn't exist or the specified unit conversion fails, the command will return ANULL, which is -999999. In your script, the program interprets the first string "mass" as the argument for a component, element, aqueous or surface species, minerals, and end members in the fluid. So the returned -999999 is telling you that there's not a mass value for "reacted" . The right command for what you want is "mass_reacted" for your result function. You might already know this but the full list of accepted arguments can be found under the Report chapter in the GWB Reference manual.

The GWB Plugin is designed to only retrieve the result of each nodal block at the end of the simulation. If you would like to retrieve results at different time slices, you can run the simulation multiple times with the desired time level set for the end time point and retrieve the information and store it in a text file. Another option is to parse the plot output file for data at each time step. The Plugin writes out a plot dataset (with extension.xtp for an X1t run) by default in binary format. You can set the output format to write out human readable xml with the command plot = character in your X1t script.

If you are looking to create simulations that you wish to report results at each time step, you can try setting up a Chemplugin run. While GWB Plugin is designed to report results at the final time step, ChemPlugin can retrieve information at each timestep. The ChemPlugin User’s Guide describes how the software works and example scripts are provided for setting up a 1D reactive transport model.

Hope this helps,
Jia

 

[How to incorporate amount of rainwater in x1t?]

Hello,

You can specify flow within the system, but not outside. If you can assume that all the rain enters your system at 1200 mm/yr, you can set that as your discharge value. The discharge unit in the GWB is a volumetric flux, which is volume/ (cross sectional area * time). I think your rainfall measurement is already in the correct unit format but you should double check that. You can read a more detailed explanation on setting flow rates in section 3.2 in the GWB Reactive Transport Modeling guide.

You set the initial condition of your system in the Initial pane and the Fluids pane is where you set the composition of your boundary fluid chemistry. H2O is typically set as a free quantity since it is the amount of water present. In one node applications, such as SpecE8 and React, varying the water mass will change the bulk mass and volume of the system. In a reactive transport model, you set the bulk volume of your system by choosing the dimensions of your domain. 

If you choose units such as mmolal or mmolal, then the concentration set is relative to 1 kg of solvent water. Other units such as mg/kg is milligram of solute per kilogram of solution. Absolute units such as moles or grams are typically not recommended, as then you have to recalculate the total moles of a component if you change the size of your domain.

Hope this helps,
Jia Wang
Aqueous Solutions LLC

 

[transition elements speciation in GSS file]

Hello,

It's difficult to say what might be the exact problem here, but here are some information on GSS that might help.

The basis species analyte available in GSS contains all primary basis species and redox species in the thermodynamic dataset loaded.

If you select a primary basis species for that element and provide the constraints for the oxidation state, such as O2(aq), Eh, or pe, then the program will figure distribution of mass from the total concentration you specified between valence states. Since your units are in elemental concentrations, then selecting the "as element" unit is right. For example, Mn++ is a basis species in thermo.tdat, the default thermo dataset. Adding this analyte to your spreadsheet and providing the constraints for oxidation will automatically prompt the speciation calculation to account for the distribution of mass to other oxidation states, MnO4- and MnO4--.

If you select a redox species to be added to the spreadsheet, the redox coupling reaction is disabled. If you know the oxidation state of your element from an analytical measurement, then you can pick that redox species for the spreadsheet and disable the redox couple.

You can open and view the thermo dataset loaded into your GSS spreadsheet by going to the File menu -> View -> and select the file that ends with the extension .tdat. Please see section 2.4 Redox couples in the GWB Essentials Guide to see how redox couples are specified in the thermodynamic dataset. I would also suggest looking at section 7.3 Redox disequilibrium for a detailed example on decoupling redox couples. Section 7.2 gives an example for how the software calculates the equilibrium state of a fluid. This is a good starting point for you to learning the software.

If you would like further help with troubleshooting, please attach your input file and any custom thermodynamic dataset so we can take a closer look.

Hope this helps,
Jia Wang
Aqueous Solutions LLC

[something weird in Rxn logK]

I am glad that you were able to figure it out. I should've have said K-feldspar being in the wrong place in the mass action equation according to your equation and not albite.

Best regards,
Jia

[something weird in Rxn logK]

Hello,

I think you have albite in the wrong place in your mass action equation. According to your reaction, albite + K+ = k-feldspar + Na+ , the equilibrium constant is calculated as:

K = [activity(K-feldspar) * activity(Na+)] / [activity(albite) * activity(K+)]

Taking log on both sides, we get:

logK = log [activity(K-feldspar) * activity(Na+)] - log[activity(albite) * activity(K+)] 

logK = log activity(K-feldspar) + log activity(Na+) - log activity(albite) - log activity(K+)

substituting the values for feldspar and albite:

logK = log (0.891) + log activity(Na+) - 0 - log activity(K+)

logK = -0.05012 + log activity(Na+) - log activity(K+)

Hope this helps,
Jia Wang
Aqueous Solutions LLC

[Running X1t from Python]

Hello James,

If you would like to run X1t input scripts directly, I think the GWB Plugin may be better suited. You can set up a X1t input file and in Python, call the GWB Plugin to run the script. You can retrieve the result from your simulation using the 'results' command. A good place to get started is the Plug-in Chapter in the GWB Reference Manual. You can find example scripts installed with your software.

ChemPlugin is a powerful software object that can be called in a client program, such as python, to configure a reactive transport model of any geometry. The user prepares a client program that spawns any number of ChemPlugin objects, called instances, that self-links into a network to represent the system being modeled. The client program at a minimum specifies the rate at which fluid flows across each link. The client can also set at each link transmissivities representing diffusion, physical mixing during flow, and heat conduction.Then, as it marches forward in time, the client prompts each instance to perform essential steps: reporting the optimum time step size, transporting mass, transferring heat, and solving the equations describing chemical reactions. While ChemPlugin is derived from the same compute engine as that of The GWB, it operates independently from the rest of the GWB modeling applications. If you are interested in learning more, visit the ChemPlugin webpage for more details regarding the software and its capabilities. You can also review the ChemPlugin Academy for example clients.  

Hope this helps,
Jia Wang
Aqueous Solutions LLC

[Optimizing X1t runs]

Hello James,

Glad to hear that helped.

Best,
Jia

[Fluid-Rock Interaction in X1t]

Hello Austin,

Thank you for attaching your scripts.

To troubleshoot, I would check if your initial system converged before any type of reactions or transport takes place. This is the most easily done by running a "go initial" command in the command pane or selecting the same option under the Run menu. This triggers the X1t to run a speciation calculation for the initial equilibrium state. When I did this, I saw that the program fails to converge. This suggests that the chemical constraints that you have set up for your initial fluid are not quite right or some critical component is missing.

Looking at the Basis and I see that most component concentrations are set at 1 mg/l except for O2(g) and Alunto glass is set in equilibrium with the initial fluid. Using your description above, it doesn’t sound like the glass is in equilibrium with the initial fluid, in that case, you should prescribe the mass in the Reactants tab, under the kinetic rate law set for the glass reaction. I unswapped Alunto glass from the Initial pane and set a small amount of SiO2(aq) in the fluid and then set 50 volume% for Alunto glass in the Reactants pane. Making these couple of changes allowed me to run the full simulation in X1t.  

Some additional information regarding setting up entries in the software. The GWB programs allow you to decouple redox pairs (e.g. Fe++ and Fe+++), which will allow you to add them as separate entries in your Basis. You can choose the redox couple to disable in the "Redox Couples…" dialog under the Config menu. I am not sure what treatment of redox species would be appropriate for your system but want to let you know that this is available. The oxidation state set by your O2(g) affects the mass distribution between all redox species unless otherwise specified. You can learn more about redox coupling in the GWB is section 2.4 Redox couples of the GWB Essentials Guide.

Hope this helps,
Jia

 

[Optimizing X1t runs]

Hello,

I couldn't run your files as they are as I don't have your modified surface datasets for Gibbsite and Goethite. Since your simulation is able to run, I would suggest that you turn on "explain step" (in the Stepping dialogs) and check the 'Follow Output' to see what is limiting the time step size. With these options enabled, you should be able to see what is the limiting constraint to your time steps.

I can't be sure how the surface complexation reactions will change your simulation here but it seems like the custom rate laws produce a much higher dissolution rate for olivine than using the Arrhenius equation and the built-in rate law. Faster kinetic reactions will typically require smaller time steps to maintain stability.

It is difficult to say how you can reduce run time since I am not really sure what the purpose behind your model. To reduce run time, it would be good to evaluate your reactions and see if they need to be constrained kinetically. Reactions that occur very fast over the course of your simulation may be prescribed as simple equilibrium reactions.

Hope this helps,
Jia Wang
Aqueous Solutions LLC