Brian Farrell

Hi Trevor, You can certainly assume equilibrium with an Al mineral to solve for the concentration of dissolved Al species. Just swap in your mineral (microcline) for Al+++ from the Basis pane and specify its mass (though the actual amount will not affect your results). Generally for minerals you specify a free mass. Hope this helps, Brian Farrell Aqueous Solutions LLC

Hi DL, I didn't get your message. As for your question, you should decouple redox pairs when you believe your system is not in redox equilibrium (low temperature natural waters, perhaps), or when you have analytical data for the same element in several redox states. If you set Eh (or pe, f O2(g), O2(aq)) and leave Fe2+/Fe3+ coupled then the equilibrium between Fe2+/Fe3+ is determined by that Eh value. However, if you decouple Fe2+/Fe3+, then the redox equilibrium between them is set by constraining each component separately (Fe2+ and Fe3+). In this case, the system Eh has no effect on Fe speciation. If you decouple the pair and only specify Fe2+, then no Fe3+ species will form (similarly decoupling the pair and only constraining Fe3+ will prevent Fe2+ from forming), so be careful with this. You should check out Section 7.3 of the GWB Essentials Guide (Release 8) or Chapter 7 of Craig Bethke's Geochemical and Biogeochemical Reaction Modeling text for more on redox disequilibrium. Hope this helps, Brian

From Craig Bethke's Geochemical and Biogeochemical Reaction Modeling text, Section 8.1.2: "The B-dot equation...is considered reasonably accurate in predicting the activities of Na+ and Cl- ions to concentrations as large as several molal, and of other species to ionic strenghts up to about 0.3 to 1 molal." You should take a look at Chapter 8 of that text, or Section 7.4 of the GWB Essentials Guide for more on the implementation of activity coefficients and activity models in GWB. As for your specific problem, I would use thermo.dat. Hope this helps, Brian

Hi, For this example, I would stick with the thermo.dat dataset. Your initial solution does not have a very high ionic strength and you are only evaporating about a third of the water. Running this example, the highest ionic strength reached is about 0.05 molal, so the B-dot equation (extended form of Debye-Huckel) works just fine. Thermo.dat also has the advantage of being parameterized from 0-300 C, unlike the phrqpitz dataset, which is defined at 25 C, though it has a limited (and largely untested) provision for calculations at other temperatures. Unfortunately you cannot just add thermodynamic data for Al and Si to the phrqpitz dataset (which uses the virial activity model), without adding appropriate virial coefficients to the bottom section of the dataset. And if you are editing a database, be sure to follow the conventions outlined in the Appendix to the GWB Reference Manual. (Elements added in alphabetical order by name (Oxygen is always first, however), the list of Basis species must begin with H2O, etc.) Hope this helps, Brian Farrell Aqueous Solutions

Hi DL, SpecE8 solves for the equilibrium distribution of aqueous species in a fluid and on surfaces, as well as mineral saturation state and the fugacity of dissolved gases. If you want to know how pH affects your equilibrium system, just change the pH value a small amount, rerun your SpecE8 model, and repeat again and again. This is where React makes things easy, as it is a program for modeling reaction paths. What you need to do is start with an initial system just like in SpecE8, then go to the reactants pane. You can titrate in aqueous species like HCl or NaOH (as simple reactants) to vary the pH, or you can use a sliding pH path. You will do the same thing for Fe and P (sliding paths are calculated in terms of activity, but if you're interested in concentation just add in Fe and P as simple reactants. The GWB Reaction Modeling Guide has sections on Titration paths (section 3.1 in Release 8.0) and Sliding activity and fugacity paths (3.6) with examples that should help you figure out how to set up your calculations. You might also try diagramming Fe in Act2, choosing pH and Fe or P activity as axes. This will give you a simple predominance diagram, rather than a complete analysis that you would get from a SpecE8 or React simulation. Hope this helps, Brian Farrell

Hi Stefan, The Thermo Datasets Appendix to the GWB Reference Manual has a section on Virial coefficients. From section A.1.10 (Release 8), "Any omitted entries for the polynmial are treated as zero values." So it looks like you don't have to worry about manually entering zero values for these coefficients. Hope this helps, Brian Farrell Aqueous Solutions LLC

You might initially set porosity to 10% and try setting 20% of the aquifer volume to be some mineral, for which you would set a dissolution kinetic rate law. As the mineral dissolves, porosity would increase and so would the water content. Complete dissolution would mean porosity reaches 30%. For the mineral, you might make up some mineral and a species that it dissolves to form that would not react with any other species in your system (fake(aq) and fake(mineral)?). I'm not sure if this is what you are looking for - it depends on the purpose of your model. Let me know if that seems ok. Brian Farrell

Glad you made some progress. I suppose what determines the quality of a model is how useful it is to you. Does this satisfy most of what you were trying to accomplish? Do you have experimental data that you are trying to fit? Brian

Hi Godfrey, React (GWB Standard) can certainly be applied to problems of CO2 sequestration. React can calculate reaction paths accounting for varying CO2 fugacity, pH, temperature, ionic strength, etc. You might simulate dissolution of bedrock or precipitation of carbonate minerals, for example. You should look into this earlier post detailing one user's experience modeling CO2 sequestration with React. https://www.rockware.com/forum/index.php?showtopic=1784 Keep in mind that GWB, even the Professional version, does not currently account for multiphase flow (migration of supercritical CO2). Still, it has a lot to offer, depending on your needs. Hope to help, Brian Farrell Aqueous Solutions LLC

Hi GWB users, Thanks for the reply Mike, you brought up a number of good points about Eh-pH diagrams and geochemical modeling in general. I'll expand just a little bit. Keep in mind that activity diagrams like this assume chemical equilibrium. When a mineral is thermodynamically favored to form but kinetically limited, you can remove it from the database manually, or more easily (and less permanently) by suppressing it. Along those lines, thermodynamic databases are not perfect, and pertinent species or minerals may need to be added, equilibrium constants adjusted, etc. An Eh-pH diagram is a simplified but powerful visual representation of a system. A speciation analysis can easily be done and should be viewed in conjunction with Eh-pH diagrams to better understand your system (concentration and activity of all species, not just the most predominant, as well as mineral saturation states). Using other diagrams, like the Pb++ solubility diagram Mike created, helps explain why certain minerals may or may not form under given conditions. As for constraining the activity of your species, the best place to start is with the total analytical concentration in molal units. In dilute solutions the activity and molal concentration are similar, though this will not be true in higher ionic strength systems. After you get an idea of the important species in your system, you might choose to account for activity coefficient corrections. Best, Brian Farrell Aqueous Solutions LLC

Hi, Simulation times in X2t can vary widely, ranging form seconds or minutes to hours or days, depending on the complexity. The spatial discretization, degree of heterogeneity of the domain (physical or chemical), the amount of species loaded and reactions to consider, surface reactions with electrostatic corrections, kinetic reactions, etc. will all have a big impact on your run times. Nx and Ny are the number of nodal blocks in the x and y directions, used in the finite difference approximation to the transport equation. You'll either set Nx and overall domain length, or Nx and deltax, the size of an individual block, to define your domain. When starting a new model, begin with the simplest possible setup. Only load what you think are the most important species, use a coarse grid, etc. When a run is over, you can look at the Run pane to see how much time was spent solving the transport vs. reaction equations. Could you let me know which version of GWB you are running, and attach scripts for the two runs you were talking about? I'll try to take a look and see if I can figure anything out. Hope to help, Brian

Hi Johan, I think what you'll need to do is modify your thermodynamic dataset to include a fictitious "CuS(aq)" species. A reaction like Covellite = CuS(aq) (which replaces the current reaction for Covellite, Covellite + H+ = Cu++ + HS-) with an arbitrarily low stability (logK = -9, perhaps) will ensure that this aqueous "species" will only exist in very small quantities. You'll need to add CuS(aq) to the database under the redox species section by summing a few reactions. There are a few ways of setting this up and I believe any will work, so here's an example: Covellite + H+ = Cu++ + HS- (logK 25 = -22.8280) HS- + 2 O2(aq) = SO4-- + H+ (logK 25 = 138.3443) CuS(aq) = Covellite (logK 25 = 9.0000) -> CuS(aq) + 2 O2(aq) = Cu++ + SO4-- (logK 25 = 124.5160) Now you should have three copper redox states (Cu+, Cu++, and CuS(aq)). Decouple each of these redox couples, and specify two (balanced of course) kinetic rate laws. 1) Bornite -> Covellite + Cu+ 2) Covellite -> Cu++ + SO4-- If you want to force Covellite to form despite being undersaturated, you'll need to set order2 to 0 when you're defining variables for your kinetic rate law (order 2 is the power of the (1-Q/K term)). This should force Covellite to form until all of the Bornite is consumed. Obviously you should be careful when ignoring the thermodynamic limitations to a reaction. See the example in Geochemical and Biogeochemical Reaction Modeling Section 17.5 about creating "fake" redox species. The GWB Reaction Modeling Guide (4.2.2 and 4.4.2) discusses nonlinear rate laws and the "order" terms. Hope this helps, Brian Farrell Aqueous Solutions Makers of The Geochemist's Workbench

The "inert" command is used in SpecE8 and React to specify unreactive mineral volume. Using the standard rock density of 2.65 g/cm3, SpecE8 and React will calculate the sorbed concentration using the total mass of equilibrium, kinetic, and inert minerals. Thus it is important to specify carefully the mass of minerals in the system. In X1t and X2t, the inert command is not used. Instead, the volume of rock is calculated in one of two ways: 1) If you do not set the porosity explicitly the program figures the porosity at each node from the volumes of fluid and “rock”, the latter being the sum of the volumes of the minerals in the nodal block. The minerals include both the equilibrium minerals swapped into the basis and the reactant minerals – those reacting according to kinetic rate laws and any set as simple reactants. 2) For simulations in which you set an explicit value for porosity, the program when setting the initial conditions figures the difference at each node between the rock volume corresponding to the porosity value and the volume of minerals in the block. This volume is then carried through the simulation as “inert” or non-reactive volume. So X1t and X2t will calculate sorbed concentrations like in SpecE8 and React (sum of equilibrium, kinetic, and inert mineral volume), but the inert volume is calculated rather than being specified. If you set porosity to 30%, for example, and do not add any minerals, then 70% of the volume is rock, and using density = 2.65 g/cm3 the sorbed concentration can be determined. If porosity is set to 30% and Quartz is set to 5%, then the inert volume is 65% and Quartz volume is 5%, leading to a total volume of 70% which goes into the calculation.

Hi, You can certainly implement adsorption in your reactive transport models using X1t/X2t. Currently available sorption models include the distribution coefficient approach (Kd), Freundlich and Langmuir isotherms, ion-exchange, and the Double Layer surface complexation model (with constant capacitance and non-electrostatic options). See chapters 9, 10, and 32 in Craig Bethke's Geochemical and Biogeochemical Reaction Modeling text for descriptions of the various sorption models and their implementation in speciation, reaction, and reactive transport models. The reference manuals which come with the GWB software will also come in handy. http://www.rockware.com/product/documentation.php?id=132 For future reference, this is the archive of old posts which have been moved from the main page. When starting a new thread, you should post on the main page to ensure a prompt response. Hope to help, Brian Farrell Aqueous Solutions LLC Makers of The Geochemist's Workbench

Tam, I'm not sure why H2(aq) and H2(g) are not in the Minteq database. Using GWB8 I can still make an Eh-pH diagram, though without H2(g) loaded the stability in the top right and bottom left corners continue to the edges of the diagram instead of stopping. Depending on your version of GWB, you may need to turn the water limit option off to even generate a plot. You can certainly add H2(aq) and H2(g) to themo_minteq.dat (I would make a copy and call it something else). Be sure to follow the instructions in the Appendix to the Reference Manual. If you're going to paste the data from thermo.dat (0 - 300 C), keep in mind that the principal temperatures are different in thermo_minteq.dat (0 - 100 C). Using thermo.dat, you can go to Rxn and write the reaction for H2(aq) + 0.5 O2(aq) = H2O . Specify the temperature as one of the new principal temperatures (ie 12.5 C) and Rxn will spit out the logK value for the reaction at that temperature. Repeat this for every temperature you need and add it to your new file. If you follow the way the reactions are written in thermo.dat, you'll need a reaction for H2(aq) added to the redox species, and for H2(g) added to the gases. Remember to add to the counters for redox species and gases, and you'll want to change O2 (g) to O2(g) so that Act2 can read the upper stability limit. Everyone is a beginner at some point. Hope this helps, Brian

Hi Tam, The default database included with GWB, thermo.dat, does not include any Antimony species. If you know the logKs for the relevant species and minerals, or can find them in the literature, you can create your own database by adding these to thermo.dat (See GWB Reference Manual, Appendix: Thermo Datasets). Alternatively, you can use the thermo.com.v8.r6+.dat or thermo_minteq.dat datasets, which include Sb. Take a good look at both of these to see which, if any, contains the species you need over the temperature range of interest to you. Keep in mind that your modeling results will only be as good as the thermodynamic data upon which they are based. Hope this helps, Brian

Hi Godfrey, A few things. I deleted all the entries for which no analytical concentration was available. Your Eh value is entered as 442.2 V - this should surely be 442.2 mV. Just with these changes I was able to get a solution in SpecE8, although with a rather large charge imbalance. This is possibly due to missing analytical data or the way the data were input. I believe where you have Cr+++ listed in mg/l as Cr, and Mn++ in mg/l as Mn, that this is unnecessary. Generally you do this for the polyatomic molecules like B(OH)3 or SO4-- when you have a measurement for Boron or Sulfur, not the whole molecule. Make sure you are being careful when entering in values for redox sensitive elements like Iron or Copper. Fe++ and Cu+ are the default basis species, and you should define the Iron and Copper components in this way if you believe your system is at equilibrium. In GSS you chose Fe+++ and Cu++, hence decoupling the Fe++/Fe+++ and Cu+/Cu++ redox pairs, indicating your system is in redox disequilibrium. Generally you choose this option if you have analytical measurements for both Fe++ and Fe+++. By decouple Fe+++ and listing no input for Fe++, you prevent all ferrous iron species or minerals from forming. See the chapter on Redox disequilibrium chapter in Craig Bethke's Geochemical and Biogeochemical Reaction Modeling, or go to the SpecE8 section of the GWB Essentials Guide and look for the Redox disequilibrium section. Finally, you can calculate species concentrations and mineral saturations in GSS or SpecE8. In GSS you simply go to data - calculate - mineral saturations. I'll post both my GSS spreadsheets and SpecE8 input files. Hope to help, Brian Farrell Aqueous Solutions Makers of The Geochemist's Workbench PS I am not sure about our current plans for short courses. We will certainly keep South Africa in mind. Godfrey.sp8 Godfrey_GSS.gss

Hi, from the GWB Reaction Modeling Guide, Release 8.0: Polythermal Reaction Paths: "In a polythermal mixing model, you set temperatures for the reactants and the initial system. Over the path, the program calculates temperature from the mass of reactant that has been added to the system, using the approximate heat capacities for fluids and minerals stored in variables cpw and cpr (on the Config - Variables dialog)." So basically you set a heat capacity for the fluid (default 1 cal/g/C) and "rock" (default 0.2 cal/g/C) which is used throughout the calculation. Because the mass of the fluid (water plus solutes) is what goes into the calculation, density is not important. Hope this helps, Brian

Pressures in hydrothermal systems aren't as big a concern as they are in say metamorphic petrology, because the deviation from the steam saturation curve isn't that large, in an absolute sense. Hydrothermal chemists are commonly more worried about error in measuring species and mineral stability (the log Ks) and in estimating activity coefficients as a function of fluid composition (the gammas). If you do make an accounting of the effect of pressure on chemical stability, you generally do so by adjusting the log Ks in the thermo dataset, rather than in the activity coefficients. To my knowledge, PHREEQC, EQ3/6, etc. have no special provision for making pressure corrections to the thermo data under hydrothermal conditions. (Remember, we are talking about confining pressure here, not the partial pressures of gases, which you set in the GWB as fugacity values.) If you feel that a pressure correction is important enough to your modeling to make the effort worthwhile, you can proceed along the U2TK route, or you can post a message to the forum inquiring about what databases others may have prepared for the conditions you're studying. Additionally there are a few old topics in the archive with questions about pressure - you might look there for additional sources of information. But in any case keep in mind that the uncertainties in the log Ks and gammas are considerable and may well overshadow the pressure corrections you end up estimating.

You are correct. Their values are in MPa, so a conversion of roughly 10 is needed. By sliding the f CO2(g) to 20, 40, 60, 80, and 100 I get pH values of 4.60, 4.45, 4.36, 4.29, and 4.25, respectively, which look to be in agreement with the results of the paper. Basically your basis is the water analysis reported in Table 3, including HCO3-. On the Reactants pane, you set a sliding f CO2(g) of 20, for example, and titrate in Dolomite, Calcite, and Siderite as simple minerals. Attached is a working script, all you need is to change the CO2(g) fugacity. Hope this helps, Brian React_CO2.txt-1.rea

Hi Tam, There are many sources of uncertainty in the calculation of saturation indices, including error in the experiments from which thermodynamic data were derived and analytical error in your own water analyses. In chapter 6 of Craig Bethke's Geochemical and Biogeochemical Reaction Modeling text there is a section entitled "Interpreting saturation indices" which describes these and several additional complications. In addition, the RIMG V70 text Thermodynamics and Kinetics of Water-Rock Interaction has a chapter about thermodynamic databases which details such errors. If you could find the original source of thermodynamic data for a mineral, along with reported logK ranges, you could potentially repeat your SI calculations using the "alter log K" function with several different logK values. This would only describe one source of uncertainty, however. To the best of my knowledge, there is no rigorous way to truly define error bars in SI calculations. This is part of the uncertainty geochemists have to deal with - you have to rely on knowledge of your field area and your best judgement in interpreting saturation indices. Brian Farrell Aqueous Solutions

Hi Olivia, If the above link doesn't work it may be that something is wrong with your HASP installation, or that HASP is not starting. Try going to the search box (click on the Windows icon) and type in services, then press enter. You should see a list of various services - look for Hasp License Manager. Click on that and a small window should come up. Look for where it says Service status - ideally it should say started. You might also try running the HASP installer by itself. Paste the text below into your search box and follow the installation instructions. I think it's normal to get a message asking if the the file was installed correctly - if everything seems ok click "the file was installed correctly," and then try revisiting the link from Tom. If the link works click on the HASP keys and diagnostics buttons and send screenshots to my email. If you get another sort of error, please send that my way. C:\Program Files (x86)\Gwb\Hasp\haspdinst -i Hope this helps, Brian Farrell Aqueous Solutions (makers of the GWB)

Hi, I was able to get a pH of about 4.4 using a fucagity of 10, so it's a little closer. Scanning through the paper it looks like the equilibrium calculations were only simulated at discrete points (fugacity = 10, 8, 6, ...) and then the mineral assemblages were titrated in. So I constrained the CO2(g) fugacity directly at the beginning and then added in minerals. It seems to be closer to the paper results, except at very low fugacity. You might try contacting the authors to see if they could provide the exact script they used. Unfortunately I couldn't find any sort of data repository with such information. Best of luck. Brian React_CO2.txt.rea

Hi, could you attach a script of yours and the paper you reference so I can compare? Thanks, Brian Farrell

Hi David, To the best of my knowledge water is generally modeled in the liquid phase. In Act2 or Tact, for instance, the activity of water (which is close to 1 in dilute solutions and decreases with increasing ionic strength) can be used as a constraint. In SpecE8 or React, (which perform speciation and reaction path calculations, respectively) water activity is calculated based on the ionic strength and temperature of your system. The fugacity of steam (water vapor) is also calculated and appears in the output results. This fugacity represents the partial pressure of an external gas buffer (large reservoir) in equilibrium with your fluid, if such a gas phase were to exist. The fugacity of steam cannot be set as a constraint. Hope this helps. Brian