MarkusLenz Posted June 1, 2015 Share Posted June 1, 2015 Hello, I am working on the following problem. We have been looking at As mobilization from soil when reducing conditions are induced. We see a decrease in As concentration in liquid when reducing conditions (i.e. Mn, Fe and sulfate reducing) persist for an extended time and suspect secondary mineral formation of either As containing minerals; or precipitation of further minerals (Fe, else) that can sorb / co-precipitate As. This I would like to model now. I am intending to use react, feed it with both soil bulk element concentrations and medium composition and slide Eh from oxidizing (400 mV) to reducing (-400 mV) conditions (simulating subsequent Mn, Fe, SR conditions). Alternatively, I could use dissolved element concentrations instead of soil bulk concentrations. The soil elemental composition is [mg / kg] (only few selected elements); we have used 100 g in 1 L volume Mn 1200 Fe 53000 Zn 400 As 2010 as Scorodite Pb 140 The medium was composed of different salts that add up to (all in mol / L) Na 0.0319 K 0.0054 Cl 0.0361 Mg 0.0020 Ca 0.0010 PO4 0.0015 CO3 0.0071 SO4 0.0021 The pH was controlled by acid base dosing at 8.2. Attached is the react file. When I press run, the program gives the following error: Step 83, Xi = .7758 (15 iterations) Charge balance: Cl- molality adjusted from .03237 to .03192 Checking basis N-R didn't converge after 400 its., maximum residual = 2.84e-11, Xi = 0.7776 Cutting step size to find solution Checking basis N-R didn't converge after 999 its., maximum residual = 2.84e-11, Xi = 0.7758 -- Didn't wake up, abandoning path Seems like some of the paths are calculated, but then the software stops. I can plot some mineral phases, but only until ~-200 mV. Do you have an advice on what to change ? Thanks a lot, Best regards, Markus soil As experiment.rea Quote Link to comment Share on other sites More sharing options...
Brian Farrell Posted June 2, 2015 Share Posted June 2, 2015 Hi Markus, If I increase the convergence criterion (Config - Iteration menu, see 6.28 epsilon in the GWB Reference Manual) I can get the model to proceed further. I also noticed that all other things held constant, decreasing the step size (Config - Stepping menu, 6.19 delxi in the Reference Manual) changes the minerals that precipitate. This didn't affect convergence at all, but typically you don't want your results to change much as you take smaller steps in your simulation. Before playing around with the parameters for the numerical solution, though, I'd make sure what you've entered into the program matches your system as best as possible. I'm slightly confused by your description. For your Na, K, etc., are all these the dissolved analytical concentrations? If so, I think you've entered them properly. But for the Mn, Fe, etc., are these concentrations of elements in the solid? If you know the mineralogy, you'll want to swap minerals in equilibrium with the system into the Basis. Section 7.2, Equilibrium models, in the GWB Essentials Guide describes swapping. You might, for example, swap Scorodite in for H2AsO4-. If you're starting with oxidized conditions, you might have an iron mineral like Goethite in the system. The database you're using doesn't have Scorodite, but you can switch to one that does (thermo.tdat) or copy the mineral into the dataset you're using. I also noticed that you included Pb in your description but Pd in your model. Finally, and this may be something you're getting to eventually, but if you want to model sorption of As then you'll need to include a description of the surface chemistry. You might set up a two-layer surface complexation model for sorption to a ferric hydroxide mineral, for example, or maybe something simpler like a Langmuir isotherm. Hope this helps, Brian Farrell Aqueous Solutions Quote Link to comment Share on other sites More sharing options...
MarkusLenz Posted June 3, 2015 Author Share Posted June 3, 2015 Hi Brian, thank you very much for your quick response. To the description of the system, it is correct that the medium was fed as dissolved salts, the soil elemental composition is based on XRF bulk values. We only know minerology for the As phase (xafs), not for any other phase (Fe, Mn). The soil was aerobic, probably it would make sense to swap Fe/Mn diss. to Fe/Mn solid phases characteristic for such soils. The intial soil Na, K, Cl, Mg, Ca, PO4 etc. I realize now, I should have included as well. However, again I do not as which phase they are present, should I add the concentration to the dissolved content or make an assumption on probable solid phases? Regarding the configuration -> iteration setting, I can see how this works for the calculation. Still, I do not really understand what the convergence criterion "epsilon" stands for and what are typical values that make sense. I need to set it to 10^-8 instead of 10^-11 to have all steps calculated, is that ok? I would very much appreciate some help here. I realized somethind strange when I started playing with the plots. H2O in fluid increased considerately from ~ -100 mV onwards. Could you comment on this please? Thanks again for your help, Best regards Markus Quote Link to comment Share on other sites More sharing options...
MarkusLenz Posted June 3, 2015 Author Share Posted June 3, 2015 sorry, forgot, Pd was a type-o in the model, the soil contained 140 ppm Pb. Best Markus Quote Link to comment Share on other sites More sharing options...
Brian Farrell Posted June 3, 2015 Share Posted June 3, 2015 Hi Markus, I would use the dissolved concentration for all the components for which you have data. For those components that you don't have data, make an assumption about minerals that might be in equilibrium with the fluid. In an oxygenated system, Hematite or Goethite might be good assumptions to swap in for the Fe++ component. Since you detected Scorodite, you could swap that in for the As. You should not enter solid phase concentrations as if they were dissolved, though. In a geochemical model, the program iterates to a solution of the equations representing the distribution of chemical mass. The iterations continue until a solution is found that satisfies the equations to an acceptably small tolerance. The default value for epsilon, the convergence criterion, is 5e-11. If you were to set epsilon to a smaller value, the program would likely need to take more steps to get to an ever so slightly more accurate solution, so this is not recommended. With a larger epsilon, the tolerance is a little larger and so the program would likely converge more quickly, or converge more successfully in solving a difficult problem. If you make epsilon too large, though, you risk getting a less accurate solution. How are you actually lowering Eh in your real system? Are you adding some reduced compound, such as H2 or acetate? In that case, it might be more realistic to use a titration path (simple -> Aqueous... -> H2(aq), for example) instead of a sliding activity path (sliding Eh). The program seems to solve a slightly simplified version of your example much more easily when I use the titration path. Regards, Brian Quote Link to comment Share on other sites More sharing options...
MarkusLenz Posted June 4, 2015 Author Share Posted June 4, 2015 Hi Brian, thanks for your help. Sorry, but I am still a bit puzzled. I would swap then for most probable minerals. You write "You should not enter solid phase concentrations as if they were dissolved, though". I dont get what you mean with this. If I have 5300 ppm of Fe in solid in the reactor, I should not enter this as "goehtite 5300 free mg/kg as Fe"? How should I enter it instead? To the second point, it is now more clear to me what epsilon is representing. Is there a numerical maximum for epsilon that still represents a reasonable accurate solution (and would be acceptable for a publication)? In one of the experiments, the Eh lowered because of natural organic matter present, in the other experiment we added lactate. I will try using titration instead. Still, when using the sliding Eh, do you know why H2O in fluid increased from 55 mol to >10^6 mol at reducing Eh ? When I plot "elemental composition of system (fluid + rock)" vs. Eh, both hydrogen and oxygen are ok. Also, system parameters -> "Mass H2O" and "Mass solution" vs. Eh are ok. Maybe there is some general problem with the model/calculations causing the program to have problems in the first place. Or maybe just a plotting bug (no prob, I would anyhow not plot this) ? Thanks a lot, Best regards, Quote Link to comment Share on other sites More sharing options...
Brian Farrell Posted June 4, 2015 Share Posted June 4, 2015 Hi Markus, I think what you have in the first paragraph (goethite swapped into the Basis, goethite 5300 free mg/kg as Fe) is correct. I meant that you should not enter 5,300 mg/kg solution for Fe++. Glad that the discussion of epsilon helped. For more information on the numerical solution, try checking out Chapter 4 of the Geochemical and Biogeochemical Reaction Modeling text. Unfortunately I do not have a good maximum value for epsilon that I wouldn't recommend exceeding. That doesn't look like a bug. Rather, it illustrates the difference between a thermodynamic component (which is an abstract concept) and an actual species. You can write the reaction for the electron in terms of the GWB's default basis species (the components) as: e- + H+ + .25 O2(aq) = .5 H2O With the decreasing Eh path, you're (artificially) increasing the activity of the electron. If you plot the concentration of the H+ and O2(aq) components along with the H2O, you'll see that the H+ and O2(aq) decrease (and go to negative concentrations) while the H2O increases. Does this make sense? For more information on the difference between components and species, try checking out section 3.2.2 in the GBRM text mentioned above. Hope this helps, Brian Quote Link to comment Share on other sites More sharing options...
MarkusLenz Posted June 9, 2015 Author Share Posted June 9, 2015 Dear Brian, thank you very much for your detailed explanation, it really helped. I am still having problems swapping aqueous vs. mineral species, i am sure this is complete basic stuff and has been discussed many times. Sorry for this naive, newbie problem. I have updated the example above, one file having a sliding Eh path, one example adding H2 to decrease Eh, using the thermo.dat as database as you recommended. Since I dont know how much lactate / NOM was in fact used, i increased H2 conc. until calculation gave ~-350 mV, which is what we measured. When i now start swapping aqueous vs. mineral species, some phases work, others result in either "inital solution too supersaturated" or simply "residuals too large". For instance, if I use "H2 addition", I can swap Fe->Goethite, but not As->Scorodite; Mn-> Birnessite. In particular I do not really understand the probelm with "intial solution too supersaturated". To me it seems intuitive that the solution needs to be supersaturated (for most elements) since despite some reductive dissolution, the soil will never "dissolve" entirely. I would be glad, if you could explain where my misunderstanding is. Many thanks Best regards, Markus H2 addition.rea sliding Eh.rea Quote Link to comment Share on other sites More sharing options...
Brian Farrell Posted June 9, 2015 Share Posted June 9, 2015 Hi Markus, Before you get too involved in the reaction path model, I think it would be a good idea to verify the equilibrium state of your initial fluid so you know you're starting in the right place. For example, you know how much Ca, Cl, C, K, Mg, Na, P, and S were dissolved into the solution, and you know the pH and Eh. A quick speciation (use SpecE8, or in React disable the "precipitation" option (from Config - Iteration...) and any reaction path options) of the fluid shows that a number of minerals are supersaturated, including Hydroxyapatite, Whitlockite, Dolomite, Calcite, Magnesite, etc. Next, use React and enable precipitation. You'll probably also want to fix pH for this. The program will in this case allow supersaturated minerals to form. Take a look at the minerals that form and how precipitation affects your solution composition. As a modeler, it's up to you to allow supersaturated minerals to form or to suppress (Config - Suppress...) any that might be unlikely to form. This is commonly an iterative process. Once you're comfortable with the initial state of the simple medium, try incorporating some of the minerals likely to exist in the soil. Hope this helps, Brian Quote Link to comment Share on other sites More sharing options...
Recommended Posts
Join the conversation
You can post now and register later. If you have an account, sign in now to post with your account.