https://mooseframework.inl.gov
MinimalGeochemicalSystemTest.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #include "gtest/gtest.h"
11 
12 #include "PertinentGeochemicalSystem.h"
13 
14 const Real eps = 1E-12; // accounts for precision loss when substituting reactions
15 
23 TEST(PertinentGeochemicalSystemTest, basisExceptions)
24 {
25  GeochemicalDatabaseReader database("database/moose_testdb.json");
26 
27  try
28  {
29  PertinentGeochemicalSystem model(
30  database, {"H2O", "Ca++", "(O-phth)--", "H2O", "Na+"}, {}, {}, {}, {}, {}, "O2(aq)", "e-");
31  FAIL() << "Missing expected exception.";
32  }
33  catch (const std::exception & e)
34  {
35  std::string msg(e.what());
36  ASSERT_TRUE(msg.find("H2O exists more than once in the basis species list") !=
37  std::string::npos)
38  << "Failed with unexpected error message: " << msg;
39  }
40 
41  try
42  {
43  PertinentGeochemicalSystem model(database,
44  {"Ca++", "H2O", "H+", "(O-phth)--", "Am++++"},
45  {},
46  {},
47  {},
48  {},
49  {},
50  "O2(aq)",
51  "e-");
52  FAIL() << "Missing expected exception.";
53  }
54  catch (const std::exception & e)
55  {
56  std::string msg(e.what());
57  ASSERT_TRUE(msg.find("First member of basis species list must be H2O") != std::string::npos)
58  << "Failed with unexpected error message: " << msg;
59  }
60 
61  try
62  {
63  PertinentGeochemicalSystem model(
64  database,
65  {"H2O", "Ca++", "H+", "(O-phth)--", "Am++++", "does_not_exist"},
66  {},
67  {},
68  {},
69  {},
70  {},
71  "O2(aq)",
72  "e-");
73  FAIL() << "Missing expected exception.";
74  }
75  catch (const std::exception & e)
76  {
77  std::string msg(e.what());
78  ASSERT_TRUE(msg.find("does_not_exist does not exist in the basis species or redox species in "
79  "database/moose_testdb.json") != std::string::npos)
80  << "Failed with unexpected error message: " << msg;
81  }
82 }
83 
93 TEST(PertinentGeochemicalSystemTest, mineralExceptions)
94 {
95  GeochemicalDatabaseReader database("database/moose_testdb.json");
96 
97  try
98  {
99  PertinentGeochemicalSystem model(database,
100  {"H2O", "Ca++", "H+", "(O-phth)--", "HCO3-"},
101  {"Calcite", "Calcite_asdf", "Calcite"},
102  {},
103  {},
104  {},
105  {},
106  "O2(aq)",
107  "e-");
108  FAIL() << "Missing expected exception.";
109  }
110  catch (const std::exception & e)
111  {
112  std::string msg(e.what());
113  ASSERT_TRUE(msg.find("Calcite exists more than once in the minerals list") != std::string::npos)
114  << "Failed with unexpected error message: " << msg;
115  }
116 
117  try
118  {
119  PertinentGeochemicalSystem model(database,
120  {"H2O", "Ca++", "H+", "Am++++", "HCO3-"},
121  {"Calcite", "Calcite_asdf", "does_not_exist"},
122  {},
123  {},
124  {},
125  {},
126  "O2(aq)",
127  "e-");
128  FAIL() << "Missing expected exception.";
129  }
130  catch (const std::exception & e)
131  {
132  std::string msg(e.what());
133  ASSERT_TRUE(msg.find("does_not_exist does not exist in database database/moose_testdb.json") !=
134  std::string::npos)
135  << "Failed with unexpected error message: " << msg;
136  }
137 
138  try
139  {
140  PertinentGeochemicalSystem model(database,
141  {"H2O", "H+", "(O-phth)--", "HCO3-"},
142  {"Calcite"},
143  {},
144  {},
145  {},
146  {},
147  "O2(aq)",
148  "e-");
149  FAIL() << "Missing expected exception.";
150  }
151  catch (const std::exception & e)
152  {
153  std::string msg(e.what());
154  ASSERT_TRUE(msg.find("The reaction for Calcite depends on Ca++ which is not reducable to a set "
155  "of basis species") != std::string::npos)
156  << "Failed with unexpected error message: " << msg;
157  }
158 
159  try
160  {
161  PertinentGeochemicalSystem model(
162  database, {"H2O", "H+", "O2(aq)"}, {"Fe(OH)3(ppd)"}, {}, {}, {}, {}, "O2(aq)", "e-");
163  FAIL() << "Missing expected exception.";
164  }
165  catch (const std::exception & e)
166  {
167  std::string msg(e.what());
168  ASSERT_TRUE(msg.find("The reaction for Fe(OH)3(ppd) depends on Fe+++ which is not reducable to "
169  "a set of basis species") != std::string::npos)
170  << "Failed with unexpected error message: " << msg;
171  }
172 
173  try
174  {
175  PertinentGeochemicalSystem model(database,
176  {"H2O", "H+", "O2(aq)", "Fe++", ">(s)FeOH"},
177  {"Fe(OH)3(ppd)"},
178  {},
179  {},
180  {},
181  {},
182  "O2(aq)",
183  "e-");
184  FAIL() << "Missing expected exception.";
185  }
186  catch (const std::exception & e)
187  {
188  std::string msg(e.what());
189  ASSERT_TRUE(msg.find("The sorbing sites for Fe(OH)3(ppd) include >(w)FeOH which is not in the "
190  "basis_species list") != std::string::npos)
191  << "Failed with unexpected error message: " << msg;
192  }
193 
194  // This should pass since the basis species can make Fe+++, which then makes Fe(OH)3(ppd), and the
195  // surface sites for Fe(OH)3(ppd) are included
196  PertinentGeochemicalSystem model(database,
197  {"H2O", "H+", "O2(aq)", "Fe++", ">(s)FeOH", ">(w)FeOH"},
198  {"Fe(OH)3(ppd)"},
199  {},
200  {},
201  {},
202  {},
203  "O2(aq)",
204  "e-");
205 }
206 
214 TEST(PertinentGeochemicalSystemTest, gasExceptions)
215 {
216  GeochemicalDatabaseReader database("database/moose_testdb.json");
217 
218  try
219  {
220  PertinentGeochemicalSystem model(database,
221  {"H2O", "Ca++", "H+", "HCO3-", "CH4(aq)"},
222  {"Calcite"},
223  {"CH4(g)", "CH4(g)"},
224  {},
225  {},
226  {},
227  "O2(aq)",
228  "e-");
229  FAIL() << "Missing expected exception.";
230  }
231  catch (const std::exception & e)
232  {
233  std::string msg(e.what());
234  ASSERT_TRUE(msg.find("CH4(g) exists more than once in the gases list") != std::string::npos)
235  << "Failed with unexpected error message: " << msg;
236  }
237 
238  try
239  {
240  PertinentGeochemicalSystem model(database,
241  {"H2O", "Ca++", "H+", "HCO3-", "CH4(aq)"},
242  {"Calcite"},
243  {"CH4(g)", "does_not_exist"},
244  {},
245  {},
246  {},
247  "O2(aq)",
248  "e-");
249  FAIL() << "Missing expected exception.";
250  }
251  catch (const std::exception & e)
252  {
253  std::string msg(e.what());
254  ASSERT_TRUE(msg.find("does_not_exist does not exist in database database/moose_testdb.json") !=
255  std::string::npos)
256  << "Failed with unexpected error message: " << msg;
257  }
258 
259  try
260  {
261  PertinentGeochemicalSystem model(database,
262  {"H2O", "Ca++", "H+", "HCO3-"},
263  {"Calcite"},
264  {"CH4(g)"},
265  {},
266  {},
267  {},
268  "O2(aq)",
269  "e-");
270  FAIL() << "Missing expected exception.";
271  }
272  catch (const std::exception & e)
273  {
274  std::string msg(e.what());
275  ASSERT_TRUE(msg.find("The reaction for CH4(g) depends on CH4(aq) which is not reducable to a "
276  "set of basis species") != std::string::npos)
277  << "Failed with unexpected error message: " << msg;
278  }
279 
280  // this should pass: CH4(g) depends on the redox couple CH4(aq) which in turn depends on the basis
281  // species provided
282  PertinentGeochemicalSystem model(database,
283  {"H2O", "Ca++", "H+", "HCO3-", "O2(aq)"},
284  {"Calcite"},
285  {"CH4(g)"},
286  {},
287  {},
288  {},
289  "O2(aq)",
290  "e-");
291 }
292 
303 TEST(PertinentGeochemicalSystemTest, kineticMineralExceptions)
304 {
305  GeochemicalDatabaseReader database("database/moose_testdb.json");
306 
307  try
308  {
309  PertinentGeochemicalSystem model(
310  database, {"H2O"}, {}, {}, {"Calcite", "Calcite_asdf", "Calcite"}, {}, {}, "O2(aq)", "e-");
311  FAIL() << "Missing expected exception.";
312  }
313  catch (const std::exception & e)
314  {
315  std::string msg(e.what());
316  ASSERT_TRUE(msg.find("Calcite exists more than once in the kinetic_minerals list") !=
317  std::string::npos)
318  << "Failed with unexpected error message: " << msg;
319  }
320 
321  try
322  {
323  PertinentGeochemicalSystem model(database,
324  {"H2O", "Ca++", "HCO3-", "H+"},
325  {"Calcite_asdf"},
326  {},
327  {"Calcite", "Calcite_asdf"},
328  {},
329  {},
330  "O2(aq)",
331  "e-");
332  FAIL() << "Missing expected exception.";
333  }
334  catch (const std::exception & e)
335  {
336  std::string msg(e.what());
337  ASSERT_TRUE(msg.find("Calcite_asdf exists in both the minerals and kinetic_minerals lists") !=
338  std::string::npos)
339  << "Failed with unexpected error message: " << msg;
340  }
341 
342  try
343  {
344  PertinentGeochemicalSystem model(database,
345  {"H2O", "Ca++", "HCO3-", "H+"},
346  {"Calcite_asdf"},
347  {},
348  {"Calcite", "does_not_exist"},
349  {},
350  {},
351  "O2(aq)",
352  "e-");
353  FAIL() << "Missing expected exception.";
354  }
355  catch (const std::exception & e)
356  {
357  std::string msg(e.what());
358  ASSERT_TRUE(msg.find("does_not_exist does not exist in database database/moose_testdb.json") !=
359  std::string::npos)
360  << "Failed with unexpected error message: " << msg;
361  }
362 
363  try
364  {
365  PertinentGeochemicalSystem model(
366  database, {"H2O", "O2(aq)", "H+"}, {}, {}, {"Fe(OH)3(ppd)"}, {}, {}, "O2(aq)", "e-");
367  FAIL() << "Missing expected exception.";
368  }
369  catch (const std::exception & e)
370  {
371  std::string msg(e.what());
372  ASSERT_TRUE(msg.find("The reaction for Fe(OH)3(ppd) depends on Fe+++ which is not reducable to "
373  "a set of basis species") != std::string::npos)
374  << "Failed with unexpected error message: " << msg;
375  }
376 
377  try
378  {
379  PertinentGeochemicalSystem model(database,
380  {"H2O", "H+", "O2(aq)", "Fe++", ">(s)FeOH"},
381  {},
382  {},
383  {"Fe(OH)3(ppd)"},
384  {},
385  {},
386  "O2(aq)",
387  "e-");
388  FAIL() << "Missing expected exception.";
389  }
390  catch (const std::exception & e)
391  {
392  std::string msg(e.what());
393  ASSERT_TRUE(msg.find("The sorbing sites for Fe(OH)3(ppd) include >(w)FeOH which is not in the "
394  "basis_species list") != std::string::npos)
395  << "Failed with unexpected error message: " << msg;
396  }
397 
398  // This should pass since the basis species can make Fe+++, which then makes Fe(OH)3(ppd), and the
399  // surface sites for Fe(OH)3(ppd) are included
400  PertinentGeochemicalSystem model(database,
401  {"H2O", "H+", "O2(aq)", "Fe++", ">(s)FeOH", ">(w)FeOH"},
402  {},
403  {},
404  {"Fe(OH)3(ppd)"},
405  {},
406  {},
407  "O2(aq)",
408  "e-");
409 }
410 
419 TEST(PertinentGeochemicalSystemTest, kineticRedoxExceptions)
420 {
421  GeochemicalDatabaseReader database("database/moose_testdb.json");
422 
423  try
424  {
425  PertinentGeochemicalSystem model(
426  database, {"H2O"}, {}, {}, {}, {"(O-phth)--", "Am++++", "(O-phth)--"}, {}, "O2(aq)", "e-");
427  FAIL() << "Missing expected exception.";
428  }
429  catch (const std::exception & e)
430  {
431  std::string msg(e.what());
432  ASSERT_TRUE(msg.find("(O-phth)-- exists more than once in the kinetic_redox list") !=
433  std::string::npos)
434  << "Failed with unexpected error message: " << msg;
435  }
436 
437  try
438  {
439  PertinentGeochemicalSystem model(
440  database, {"H2O", "(O-phth)--"}, {}, {}, {}, {"Am++++", "(O-phth)--"}, {}, "O2(aq)", "e-");
441  FAIL() << "Missing expected exception.";
442  }
443  catch (const std::exception & e)
444  {
445  std::string msg(e.what());
446  ASSERT_TRUE(msg.find("(O-phth)-- exists in both the basis_species and kinetic_redox lists") !=
447  std::string::npos)
448  << "Failed with unexpected error message: " << msg;
449  }
450 
451  try
452  {
453  PertinentGeochemicalSystem model(
454  database, {"H2O"}, {}, {}, {}, {"(O-phth)--", "does_not_exist"}, {}, "O2(aq)", "e-");
455  FAIL() << "Missing expected exception.";
456  }
457  catch (const std::exception & e)
458  {
459  std::string msg(e.what());
460  ASSERT_TRUE(msg.find("does_not_exist does not exist in database database/moose_testdb.json") !=
461  std::string::npos)
462  << "Failed with unexpected error message: " << msg;
463  }
464 
465  try
466  {
467  PertinentGeochemicalSystem model(
468  database, {"H2O", "HCO3-", "H+"}, {}, {}, {}, {"(O-phth)--"}, {}, "O2(aq)", "e-");
469  FAIL() << "Missing expected exception.";
470  }
471  catch (const std::exception & e)
472  {
473  std::string msg(e.what());
474  ASSERT_TRUE(msg.find("The reaction for (O-phth)-- depends on O2(aq) which is not reducable to "
475  "a set of basis species") != std::string::npos)
476  << "Failed with unexpected error message: " << msg;
477  }
478 
479  // this should pass
480  PertinentGeochemicalSystem model(
481  database, {"H2O", "HCO3-", "H+", "O2(aq)"}, {}, {}, {}, {"(O-phth)--"}, {}, "O2(aq)", "e-");
482 }
483 
491 TEST(PertinentGeochemicalSystemTest, kineticSurfaceExceptions)
492 {
493  GeochemicalDatabaseReader database("database/moose_testdb.json");
494 
495  try
496  {
497  PertinentGeochemicalSystem model(
498  database, {"H2O"}, {}, {}, {}, {}, {">(s)FeO-", ">(s)FeOCa+", ">(s)FeO-"}, "O2(aq)", "e-");
499  FAIL() << "Missing expected exception.";
500  }
501  catch (const std::exception & e)
502  {
503  std::string msg(e.what());
504  ASSERT_TRUE(msg.find(">(s)FeO- exists more than once in the kinetic_surface_species list") !=
505  std::string::npos)
506  << "Failed with unexpected error message: " << msg;
507  }
508 
509  try
510  {
511  PertinentGeochemicalSystem model(database,
512  {"H2O"},
513  {},
514  {},
515  {},
516  {},
517  {">(s)FeO-", ">(s)FeOCa+", "does_not_exist"},
518  "O2(aq)",
519  "e-");
520  FAIL() << "Missing expected exception.";
521  }
522  catch (const std::exception & e)
523  {
524  std::string msg(e.what());
525  ASSERT_TRUE(msg.find("does_not_exist does not exist in database database/moose_testdb.json") !=
526  std::string::npos)
527  << "Failed with unexpected error message: " << msg;
528  }
529 
530  try
531  {
532  PertinentGeochemicalSystem model(
533  database, {"H2O", "H+"}, {}, {}, {}, {}, {">(s)FeO-"}, "O2(aq)", "e-");
534  FAIL() << "Missing expected exception.";
535  }
536  catch (const std::exception & e)
537  {
538  std::string msg(e.what());
539  ASSERT_TRUE(msg.find("The reaction for >(s)FeO- depends on >(s)FeOH which is not reducable to "
540  "a set of basis species") != std::string::npos)
541  << "Failed with unexpected error message: " << msg;
542  }
543 
544  // this should pass
545  PertinentGeochemicalSystem model(
546  database, {"H2O", "H+", ">(s)FeOH"}, {}, {}, {}, {}, {">(s)FeO-"}, "O2(aq)", "e-");
547 }
548 
550 TEST(PertinentGeochemicalSystemTest, temperatures)
551 {
552  GeochemicalDatabaseReader database("database/moose_testdb.json");
553 
554  PertinentGeochemicalSystem model(database, {"H2O"}, {}, {}, {}, {}, {}, "O2(aq)", "e-");
555  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
556 
557  ASSERT_EQ(mgd.original_database->getTemperatures(), database.getTemperatures());
558 }
559 
565 TEST(PertinentGeochemicalSystemTest, names1)
566 {
567  GeochemicalDatabaseReader database("database/moose_testdb.json");
568 
569  // The following system has secondary species: CO2(aq), CO3--, CaCO3, CaOH+, OH-, (O-phth)--,
570  // >(s)FeO-
571  PertinentGeochemicalSystem model(database,
572  {"H2O", "H+", "HCO3-", "O2(aq)", "Ca++", ">(s)FeOH"},
573  {"Calcite"},
574  {"O2(g)"},
575  {"Calcite_asdf"},
576  {"CH4(aq)"},
577  {">(s)FeOCa+"},
578  "O2(aq)",
579  "e-");
580  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
581 
582  ASSERT_EQ(mgd.basis_species_index.size(), (std::size_t)6);
583  for (const auto & sp : mgd.basis_species_index)
584  ASSERT_EQ(mgd.basis_species_name[sp.second], sp.first);
585 
586  ASSERT_EQ(mgd.eqm_species_index.size(), (std::size_t)9);
587  for (const auto & sp : mgd.eqm_species_index)
588  ASSERT_EQ(mgd.eqm_species_name[sp.second], sp.first);
589 
590  ASSERT_EQ(mgd.kin_species_index.size(), (std::size_t)3);
591  for (const auto & sp : mgd.kin_species_index)
592  ASSERT_EQ(mgd.kin_species_name[sp.second], sp.first);
593 }
594 
600 TEST(PertinentGeochemicalSystemTest, names2)
601 {
602  GeochemicalDatabaseReader database("database/moose_testdb.json");
603 
604  // The following system has secondary species: CO2(aq), CO3--, CaCO3, CaOH+, OH-, (O-phth)--
605  PertinentGeochemicalSystem model(database,
606  {"H2O", "H+", "HCO3-", "O2(aq)", "Ca++", ">(s)FeOH"},
607  {},
608  {"O2(g)"},
609  {"Calcite_asdf", "Calcite"},
610  {"CH4(aq)"},
611  {">(s)FeOCa+", ">(s)FeO-"},
612  "O2(aq)",
613  "e-");
614  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
615 
616  ASSERT_EQ(mgd.basis_species_index.size(), (std::size_t)6);
617  for (const auto & sp : mgd.basis_species_index)
618  ASSERT_EQ(mgd.basis_species_name[sp.second], sp.first);
619 
620  ASSERT_EQ(mgd.eqm_species_index.size(), (std::size_t)7);
621  for (const auto & sp : mgd.eqm_species_index)
622  ASSERT_EQ(mgd.eqm_species_name[sp.second], sp.first);
623 
624  ASSERT_EQ(mgd.kin_species_index.size(), (std::size_t)5);
625  for (const auto & sp : mgd.kin_species_index)
626  ASSERT_EQ(mgd.kin_species_name[sp.second], sp.first);
627 }
628 
630 TEST(PertinentGeochemicalSystemTest, charge)
631 {
632  GeochemicalDatabaseReader database("database/moose_testdb.json");
633 
634  // The following system has secondary species: CO2(aq), CO3--, CaCO3, CaOH+, OH-, (O-phth)--,
635  // >(s)FeO-, e-
636  PertinentGeochemicalSystem model(database,
637  {"H2O", "H+", "HCO3-", "O2(aq)", "Ca++", ">(s)FeOH"},
638  {"Calcite"},
639  {},
640  {"Calcite_asdf"},
641  {"CH4(aq)"},
642  {">(s)FeOCa+"},
643  "O2(aq)",
644  "e-");
645  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
646 
647  std::map<std::string, Real> charge_gold;
648  charge_gold["H2O"] = 0.0;
649  charge_gold["H+"] = 1.0;
650  charge_gold["HCO3-"] = -1.0;
651  charge_gold["O2(aq)"] = 0.0;
652  charge_gold["Ca++"] = 2.0;
653  charge_gold[">(s)FeOH"] = 0.0;
654  charge_gold["CO2(aq)"] = 0.0;
655  charge_gold["CO3--"] = -2.0;
656  charge_gold["CaCO3"] = 0.0;
657  charge_gold["CaOH+"] = 1.0;
658  charge_gold["OH-"] = -1.0;
659  charge_gold["(O-phth)--"] = -2.0;
660  charge_gold[">(s)FeO-"] = -1.0;
661  charge_gold[">(s)FeOCa+"] = -1.0;
662  charge_gold["Calcite"] = 0.0;
663  charge_gold["e-"] = -1.0;
664  charge_gold["Calcite_asdf"] = 0.0;
665  charge_gold["CH4(aq)"] = 0.0;
666  for (const auto & sp : mgd.basis_species_index)
667  ASSERT_EQ(mgd.basis_species_charge[sp.second], charge_gold[sp.first]);
668  for (const auto & sp : mgd.eqm_species_index)
669  ASSERT_EQ(mgd.eqm_species_charge[sp.second], charge_gold[sp.first]);
670  for (const auto & sp : mgd.kin_species_index)
671  ASSERT_EQ(mgd.kin_species_charge[sp.second], charge_gold[sp.first]);
672 }
673 
675 TEST(PertinentGeochemicalSystemTest, radius)
676 {
677  GeochemicalDatabaseReader database("database/moose_testdb.json");
678 
679  // The following system has secondary species: CO2(aq), CO3--, CaCO3, CaOH+, OH-, (O-phth)--,
680  // >(s)FeO-, >(s)FeOCa+, e-
681  PertinentGeochemicalSystem model(database,
682  {"H2O", "H+", "HCO3-", "O2(aq)", "Ca++", ">(s)FeOH"},
683  {"Calcite"},
684  {},
685  {},
686  {"CH4(aq)"},
687  {},
688  "O2(aq)",
689  "e-");
690  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
691 
692  std::map<std::string, Real> radius_gold;
693  radius_gold["H2O"] = 0.0;
694  radius_gold["H+"] = 9.0;
695  radius_gold["HCO3-"] = 4.5;
696  radius_gold["O2(aq)"] = -0.5;
697  radius_gold["Ca++"] = 6.0;
698  radius_gold[">(s)FeOH"] = 0.0;
699  radius_gold["CO2(aq)"] = 4.0;
700  radius_gold["CO3--"] = 4.5;
701  radius_gold["CaCO3"] = 4.0;
702  radius_gold["CaOH+"] = 4.0;
703  radius_gold["OH-"] = 3.5;
704  radius_gold["(O-phth)--"] = 4.0;
705  radius_gold[">(s)FeO-"] = -1.5;
706  radius_gold[">(s)FeOCa+"] = -1.5;
707  radius_gold["Calcite"] = 0.0;
708  radius_gold["e-"] = 0.0;
709  for (const auto & sp : mgd.basis_species_index)
710  ASSERT_EQ(mgd.basis_species_radius[sp.second], radius_gold[sp.first]);
711  for (const auto & sp : mgd.eqm_species_index)
712  ASSERT_EQ(mgd.eqm_species_radius[sp.second], radius_gold[sp.first]);
713 }
714 
716 TEST(PertinentGeochemicalSystemTest, molecular_weight)
717 {
718  GeochemicalDatabaseReader database("database/moose_testdb.json");
719 
720  // The following system has secondary species: CO2(aq), CO3--, CaCO3, CaOH+, OH-, (O-phth)--,
721  // >(s)FeO-, >(s)FeOCa+, e-
722  PertinentGeochemicalSystem model(database,
723  {"H2O", "H+", "HCO3-", "O2(aq)", "Ca++", ">(s)FeOH"},
724  {"Calcite"},
725  {},
726  {"Calcite_asdf"},
727  {"CH4(aq)"},
728  {">(s)FeOCa+"},
729  "O2(aq)",
730  "e-");
731  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
732 
733  std::map<std::string, Real> molecular_weight_gold;
734  molecular_weight_gold["H2O"] = 18.0152;
735  molecular_weight_gold["H+"] = 1.0079;
736  molecular_weight_gold["HCO3-"] = 61.0171;
737  molecular_weight_gold["O2(aq)"] = 31.9988;
738  molecular_weight_gold["Ca++"] = 40.0800;
739  molecular_weight_gold[">(s)FeOH"] = 72.8543;
740  molecular_weight_gold["CO2(aq)"] = 44.0098;
741  molecular_weight_gold["CO3--"] = 60.0092;
742  molecular_weight_gold["CaCO3"] = 100.0892;
743  molecular_weight_gold["CaOH+"] = 57.0873;
744  molecular_weight_gold["OH-"] = 17.0073;
745  molecular_weight_gold["(O-phth)--"] = 164.1172;
746  molecular_weight_gold[">(s)FeO-"] = 71.8464;
747  molecular_weight_gold[">(s)FeOCa+"] = 279.7144;
748  molecular_weight_gold["Calcite"] = 100.0892;
749  molecular_weight_gold["e-"] = 0.0;
750  molecular_weight_gold["Calcite_asdf"] = 111.0892;
751  molecular_weight_gold["CH4(aq)"] = 16.0426;
752  molecular_weight_gold[">(s)FeOCa+"] = 279.7144;
753  for (const auto & sp : mgd.basis_species_index)
754  ASSERT_EQ(mgd.basis_species_molecular_weight[sp.second], molecular_weight_gold[sp.first]);
755  for (const auto & sp : mgd.eqm_species_index)
756  ASSERT_EQ(mgd.eqm_species_molecular_weight[sp.second], molecular_weight_gold[sp.first]);
757  for (const auto & sp : mgd.kin_species_index)
758  ASSERT_EQ(mgd.kin_species_molecular_weight[sp.second], molecular_weight_gold[sp.first]);
759 }
760 
762 TEST(PertinentGeochemicalSystemTest, molecular_volume)
763 {
764  GeochemicalDatabaseReader database("database/moose_testdb.json");
765 
766  // The following system has secondary species: CO2(aq), CO3--, CaCO3, CaOH+, OH-, (O-phth)--,
767  // >(s)FeO-, >(s)FeOCa+, e-
768  PertinentGeochemicalSystem model(database,
769  {"H2O", "H+", "HCO3-", "O2(aq)", "Ca++", ">(s)FeOH"},
770  {"Calcite"},
771  {},
772  {"Calcite_asdf"},
773  {"CH4(aq)"},
774  {">(s)FeOCa+"},
775  "O2(aq)",
776  "e-");
777  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
778 
779  std::map<std::string, Real> molecular_volume_gold;
780  molecular_volume_gold["H2O"] = 0.0;
781  molecular_volume_gold["H+"] = 0.0;
782  molecular_volume_gold["HCO3-"] = 0.0;
783  molecular_volume_gold["O2(aq)"] = 0.0;
784  molecular_volume_gold["Ca++"] = 0.0;
785  molecular_volume_gold[">(s)FeOH"] = 0.0;
786  molecular_volume_gold["CO2(aq)"] = 0.0;
787  molecular_volume_gold["CO3--"] = 0.0;
788  molecular_volume_gold["CaCO3"] = 0.0;
789  molecular_volume_gold["CaOH+"] = 0.0;
790  molecular_volume_gold["OH-"] = 0.0;
791  molecular_volume_gold["(O-phth)--"] = 0.0;
792  molecular_volume_gold[">(s)FeO-"] = 0.0;
793  molecular_volume_gold[">(s)FeOCa+"] = 0.0;
794  molecular_volume_gold["Calcite"] = 36.9340;
795  molecular_volume_gold["e-"] = 0.0;
796  molecular_volume_gold["Calcite_asdf"] = 136.9340;
797  molecular_volume_gold["CH4(aq)"] = 0.0;
798  molecular_volume_gold[">(s)FeOCa+"] = 0.0;
799  for (const auto & sp : mgd.basis_species_index)
800  ASSERT_EQ(mgd.basis_species_molecular_volume[sp.second], molecular_volume_gold[sp.first]);
801  for (const auto & sp : mgd.eqm_species_index)
802  ASSERT_EQ(mgd.eqm_species_molecular_volume[sp.second], molecular_volume_gold[sp.first]);
803  for (const auto & sp : mgd.kin_species_index)
804  ASSERT_EQ(mgd.kin_species_molecular_volume[sp.second], molecular_volume_gold[sp.first]);
805 }
806 
808 TEST(PertinentGeochemicalSystemTest, surfaceComplexationInfo)
809 {
810  GeochemicalDatabaseReader database("database/moose_testdb.json");
811 
812  // The following system has secondary species: CO2(aq), CO3--, CaCO3, CaOH+, OH-, (O-phth)--,
813  // >(s)FeO-, >(s)FeOCa+, e-
814  PertinentGeochemicalSystem model1(database,
815  {"H2O", "H+", "HCO3-", "O2(aq)", "Ca++", ">(s)FeOH"},
816  {"Calcite"},
817  {},
818  {},
819  {"CH4(aq)"},
820  {},
821  "O2(aq)",
822  "e-");
823  ModelGeochemicalDatabase mgd1 = model1.modelGeochemicalDatabase();
824 
825  ASSERT_EQ(mgd1.surface_complexation_info.size(), (std::size_t)0);
826 
827  // The following system has secondary species: CO2(aq), CO3--, CaCO3, CaOH+, OH-, (O-phth)--,
828  // >(s)FeO-, >(s)FeOCa+, e-
829  PertinentGeochemicalSystem model2(database,
830  {"H2O", "H+", "O2(aq)", "Fe++", ">(s)FeOH", ">(w)FeOH"},
831  {"Fe(OH)3(ppd)"},
832  {},
833  {},
834  {},
835  {},
836  "O2(aq)",
837  "e-");
838  ModelGeochemicalDatabase mgd2 = model2.modelGeochemicalDatabase();
839 
840  ASSERT_EQ(mgd2.surface_complexation_info.count("Fe(OH)3(ppd)"), (std::size_t)1);
841  ASSERT_EQ(mgd2.surface_complexation_info["Fe(OH)3(ppd)"].surface_area, 600.0);
842  ASSERT_EQ(mgd2.surface_complexation_info["Fe(OH)3(ppd)"].sorption_sites[">(s)FeOH"], 0.005);
843  ASSERT_EQ(mgd2.surface_complexation_info["Fe(OH)3(ppd)"].sorption_sites[">(w)FeOH"], 0.2);
844 
845  // The following system has secondary species: CO2(aq), CO3--, CaCO3, CaOH+, OH-, (O-phth)--,
846  // >(s)FeO-, >(s)FeOCa+, e-
847  PertinentGeochemicalSystem model3(database,
848  {"H2O", "H+", "O2(aq)", "Fe++", ">(s)FeOH", ">(w)FeOH"},
849  {},
850  {},
851  {"Goethite"},
852  {},
853  {},
854  "O2(aq)",
855  "e-");
856  ModelGeochemicalDatabase mgd3 = model3.modelGeochemicalDatabase();
857 
858  ASSERT_EQ(mgd3.surface_complexation_info.count("Goethite"), (std::size_t)1);
859  ASSERT_EQ(mgd3.surface_complexation_info["Goethite"].surface_area, 60.0);
860  ASSERT_EQ(mgd3.surface_complexation_info["Goethite"].sorption_sites[">(s)FeOH"], 0.05);
861  ASSERT_EQ(mgd3.surface_complexation_info["Goethite"].sorption_sites[">(w)FeOH"], 0.222);
862 }
863 
865 TEST(PertinentGeochemicalSystemTest, surfaceComplexationRepeatedException)
866 {
867  GeochemicalDatabaseReader database("database/moose_testdb.json");
868 
869  try
870  {
871  PertinentGeochemicalSystem model3(database,
872  {"H2O", "H+", "O2(aq)", "Fe++", ">(s)FeOH", ">(w)FeOH"},
873  {"Goethite", "Fe(OH)3(ppd)"},
874  {},
875  {},
876  {},
877  {},
878  "O2(aq)",
879  "e-");
880  FAIL() << "Missing expected exception.";
881  }
882  catch (const std::exception & e)
883  {
884  std::string msg(e.what());
885  ASSERT_TRUE(msg.find("The sorbing site >(s)FeOH appears in more than one sorbing mineral") !=
886  std::string::npos)
887  << "Failed with unexpected error message: " << msg;
888  }
889 }
890 
892 TEST(PertinentGeochemicalSystemTest, excessSorbingSitesException)
893 {
894  GeochemicalDatabaseReader database("database/moose_testdb.json");
895 
896  try
897  {
898  PertinentGeochemicalSystem model2(database,
899  {"H2O", "H+", "Fe+++", "sorbsite1", "sorbsite2"},
900  {"problematic_sorber"},
901  {},
902  {},
903  {},
904  {},
905  "O2(aq)",
906  "e-");
907  FAIL() << "Missing expected exception.";
908  }
909  catch (const std::exception & e)
910  {
911  std::string msg(e.what());
912  ASSERT_TRUE(msg.find("It is an error for any equilibrium species (such as problem_eqm) to have "
913  "a reaction involving more than one sorbing site") != std::string::npos)
914  << "Failed with unexpected error message: " << msg;
915  }
916 }
917 
919 TEST(PertinentGeochemicalSystemTest, surfaceSorptionBuilding)
920 {
921  GeochemicalDatabaseReader database("database/moose_testdb.json");
922 
923  PertinentGeochemicalSystem modelC(
924  database, {"H2O", "H+", "Ca++", "HCO3-"}, {"Calcite"}, {}, {}, {}, {}, "O2(aq)", "e-");
925  ModelGeochemicalDatabase mgdC = modelC.modelGeochemicalDatabase();
926 
927  EXPECT_EQ(mgdC.surface_sorption_name.size(), (std::size_t)0);
928  EXPECT_EQ(mgdC.surface_sorption_area.size(), (std::size_t)0);
929  for (unsigned j = 0; j < mgdC.eqm_species_name.size(); ++j)
930  EXPECT_EQ(mgdC.surface_sorption_related[j], false);
931 
932  PertinentGeochemicalSystem model2(database,
933  {"H2O", "H+", "Fe+++", ">(s)FeOH", ">(w)FeOH"},
934  {"Fe(OH)3(ppd)fake", "Goethite"},
935  {},
936  {},
937  {},
938  {},
939  "O2(aq)",
940  "e-");
941  ModelGeochemicalDatabase mgd2 = model2.modelGeochemicalDatabase();
942  EXPECT_EQ(mgd2.surface_complexation_info.count("Goethite"), (std::size_t)1);
943  EXPECT_EQ(mgd2.surface_sorption_name.size(), (std::size_t)1);
944  EXPECT_EQ(mgd2.surface_sorption_name[0], "Goethite");
945  EXPECT_EQ(mgd2.surface_sorption_area.size(), (std::size_t)1);
946  EXPECT_EQ(mgd2.surface_sorption_area[0], 60.0);
947  EXPECT_EQ(mgd2.surface_sorption_related.size(), mgd2.eqm_species_name.size());
948  EXPECT_EQ(mgd2.surface_sorption_number.size(), mgd2.eqm_species_name.size());
949  // Eqm species are: OH- >(s)FeO- Fe(OH)3(ppd)fake Goethite
950  EXPECT_EQ(mgd2.eqm_species_index.count(">(s)FeO-"), (std::size_t)1);
951  const unsigned posn = mgd2.eqm_species_index.at(">(s)FeO-");
952  for (unsigned j = 0; j < mgd2.eqm_species_name.size(); ++j)
953  if (j == posn)
954  EXPECT_TRUE(mgd2.surface_sorption_related[j]);
955  else
956  EXPECT_FALSE(mgd2.surface_sorption_related[j]);
957  EXPECT_EQ(mgd2.surface_sorption_number[posn], (unsigned int)0);
958 }
959 
961 TEST(PertinentGeochemicalSystemTest, GasChi)
962 {
963  GeochemicalDatabaseReader database("database/moose_testdb.json");
964 
965  // The following system has secondary species: CO2(aq), CO3--, CaCO3, CaOH+, OH-, (O-phth)--,
966  // >(s)FeO-, >(s)FeOCa+, e-
967  PertinentGeochemicalSystem model1(database,
968  {"H2O", "H+", "HCO3-", "O2(aq)", "Ca++", ">(s)FeOH"},
969  {"Calcite"},
970  {},
971  {},
972  {"CH4(aq)"},
973  {},
974  "O2(aq)",
975  "e-");
976  ModelGeochemicalDatabase mgd1 = model1.modelGeochemicalDatabase();
977 
978  ASSERT_EQ(mgd1.gas_chi.size(), (std::size_t)0);
979 
980  PertinentGeochemicalSystem model2(database,
981  {"H2O", "Ca++", "H+", "HCO3-", "O2(aq)"},
982  {"Calcite"},
983  {"CH4(g)"},
984  {},
985  {},
986  {},
987  "O2(aq)",
988  "e-");
989  ModelGeochemicalDatabase mgd2 = model2.modelGeochemicalDatabase();
990 
991  ASSERT_EQ(mgd2.gas_chi.count("CH4(g)"), (std::size_t)1);
992  ASSERT_EQ(mgd2.gas_chi["CH4(g)"][0], -537.779);
993  ASSERT_EQ(mgd2.gas_chi["CH4(g)"][1], 1.54946);
994  ASSERT_EQ(mgd2.gas_chi["CH4(g)"][2], -.000927827);
995  ASSERT_EQ(mgd2.gas_chi["CH4(g)"][3], 1.20861);
996  ASSERT_EQ(mgd2.gas_chi["CH4(g)"][4], -.00370814);
997  ASSERT_EQ(mgd2.gas_chi["CH4(g)"][5], 3.33804e-6);
998 }
999 
1003 TEST(PertinentGeochemicalSystemTest, stoichiometry1)
1004 {
1005  GeochemicalDatabaseReader database("database/moose_testdb.json");
1006 
1007  // The following system has secondary species: CO2(aq), CO3--, CaCO3, CaOH+, OH-, (O-phth)--
1008  PertinentGeochemicalSystem model(database,
1009  {"H2O", "H+", "HCO3-", "O2(aq)", "Ca++", ">(s)FeOH"},
1010  {"Calcite"},
1011  {},
1012  {"Calcite_asdf"},
1013  {"CH4(aq)"},
1014  {">(s)FeO-", ">(s)FeOCa+"},
1015  "O2(aq)",
1016  "e-");
1017  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
1018 
1019  ASSERT_EQ(mgd.eqm_species_index.size(), (std::size_t)7);
1020  for (const auto & sp : {"CO2(aq)", "CO3--", "CaCO3", "CaOH+", "OH-", "(O-phth)--", "Calcite"})
1021  ASSERT_EQ(mgd.eqm_species_index.count(sp), (std::size_t)1);
1022 
1023  std::map<std::string, DenseMatrix<Real>> stoi_gold;
1024  for (const auto & sp : {"CO2(aq)",
1025  "CO3--",
1026  "CaCO3",
1027  "CaOH+",
1028  "OH-",
1029  "(O-phth)--",
1030  ">(s)FeO-",
1031  ">(s)FeOCa+",
1032  "Calcite",
1033  "Calcite_asdf",
1034  "CH4(aq)"})
1035  stoi_gold[sp] = DenseMatrix<Real>(1, 6);
1036  // remember the order of primaries: {"H2O", "H+", "HCO3-", "O2(aq)", "Ca++", ">(s)FeOH"},
1037  stoi_gold["CO2(aq)"](0, 0) = -1;
1038  stoi_gold["CO2(aq)"](0, 1) = 1;
1039  stoi_gold["CO2(aq)"](0, 2) = 1;
1040  stoi_gold["CO3--"](0, 1) = -1;
1041  stoi_gold["CO3--"](0, 2) = 1;
1042  stoi_gold["CaCO3"](0, 4) = 1;
1043  stoi_gold["CaCO3"](0, 2) = 1;
1044  stoi_gold["CaCO3"](0, 1) = -1;
1045  stoi_gold["CaOH+"](0, 4) = 1;
1046  stoi_gold["CaOH+"](0, 0) = 1;
1047  stoi_gold["CaOH+"](0, 1) = -1;
1048  stoi_gold["OH-"](0, 0) = 1;
1049  stoi_gold["OH-"](0, 1) = -1;
1050  stoi_gold["(O-phth)--"](0, 0) = -5;
1051  stoi_gold["(O-phth)--"](0, 2) = 8;
1052  stoi_gold["(O-phth)--"](0, 1) = 6;
1053  stoi_gold["(O-phth)--"](0, 3) = -7.5;
1054  stoi_gold[">(s)FeO-"](0, 5) = 1;
1055  stoi_gold[">(s)FeO-"](0, 1) = -1;
1056  stoi_gold[">(s)FeOCa+"](0, 5) = 1;
1057  stoi_gold[">(s)FeOCa+"](0, 4) = 1;
1058  stoi_gold[">(s)FeOCa+"](0, 1) = -1;
1059  stoi_gold["Calcite"](0, 4) = 1;
1060  stoi_gold["Calcite"](0, 2) = 1;
1061  stoi_gold["Calcite"](0, 1) = -1;
1062  stoi_gold["Calcite_asdf"](0, 4) = 2;
1063  stoi_gold["Calcite_asdf"](0, 2) = 1;
1064  stoi_gold["Calcite_asdf"](0, 1) = -1;
1065  // remember the order of primaries: {"H2O", "H+", "HCO3-", "O2(aq)", "Ca++", ">(s)FeOH"},
1066  stoi_gold["CH4(aq)"](0, 0) = 1;
1067  stoi_gold["CH4(aq)"](0, 1) = 1;
1068  stoi_gold["CH4(aq)"](0, 2) = 1;
1069  stoi_gold["CH4(aq)"](0, 3) = -2;
1070 
1071  for (const auto & sp : {"CO2(aq)", "CO3--", "CaCO3", "CaOH+", "OH-", "(O-phth)--", "Calcite"})
1072  {
1073  const unsigned row = mgd.eqm_species_index[sp];
1074  ASSERT_EQ(mgd.eqm_stoichiometry.sub_matrix(row, 1, 0, 6), stoi_gold[sp]);
1075  }
1076  for (const auto & sp : {"Calcite_asdf", "CH4(aq)", ">(s)FeO-", ">(s)FeOCa+"})
1077  {
1078  const unsigned row = mgd.kin_species_index[sp];
1079  ASSERT_EQ(mgd.kin_stoichiometry.sub_matrix(row, 1, 0, 6), stoi_gold[sp]);
1080  }
1081 }
1082 
1084 TEST(PertinentGeochemicalSystemTest, log10K1)
1085 {
1086  GeochemicalDatabaseReader database("database/moose_testdb.json");
1087 
1088  // The following system has secondary species: CO2(aq), CO3--, CaCO3, CaOH+, OH-, (O-phth)--,
1089  // >(s)FeO-, >(s)FeOCa+
1090  PertinentGeochemicalSystem model(database,
1091  {"H2O", "H+", "HCO3-", "O2(aq)", "Ca++", ">(s)FeOH"},
1092  {"Calcite"},
1093  {},
1094  {},
1095  {"CH4(aq)"},
1096  {},
1097  "O2(aq)",
1098  "e-");
1099  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
1100 
1101  ASSERT_EQ(mgd.eqm_log10K(mgd.eqm_species_index["CO2(aq)"], 0), -6.5570);
1102  ASSERT_EQ(mgd.eqm_log10K(mgd.eqm_species_index["CO2(aq)"], 1), -6.3660);
1103  ASSERT_EQ(mgd.eqm_log10K(mgd.eqm_species_index["CO2(aq)"], 2), -6.3325);
1104  ASSERT_EQ(mgd.eqm_log10K(mgd.eqm_species_index["CO2(aq)"], 3), -6.4330);
1105  ASSERT_EQ(mgd.eqm_log10K(mgd.eqm_species_index["CO2(aq)"], 4), -6.7420);
1106  ASSERT_EQ(mgd.eqm_log10K(mgd.eqm_species_index["CO2(aq)"], 5), -7.1880);
1107  ASSERT_EQ(mgd.eqm_log10K(mgd.eqm_species_index["CO2(aq)"], 6), -7.7630);
1108  ASSERT_EQ(mgd.eqm_log10K(mgd.eqm_species_index["CO2(aq)"], 7), -8.4650);
1109  ASSERT_EQ(mgd.eqm_log10K(mgd.eqm_species_index["CO3--"], 0), 10.6169);
1110  ASSERT_EQ(mgd.eqm_log10K(mgd.eqm_species_index["CaCO3"], 0), 7.5520);
1111  ASSERT_EQ(mgd.eqm_log10K(mgd.eqm_species_index["CaOH+"], 0), 13.7095);
1112  ASSERT_EQ(mgd.eqm_log10K(mgd.eqm_species_index["OH-"], 0), 14.9325);
1113  ASSERT_EQ(mgd.eqm_log10K(mgd.eqm_species_index["(O-phth)--"], 0), 594.3211);
1114  ASSERT_EQ(mgd.eqm_log10K(mgd.eqm_species_index[">(s)FeO-"], 0), 8.9300);
1115  ASSERT_EQ(mgd.eqm_log10K(mgd.eqm_species_index[">(s)FeOCa+"], 0), 1.7200);
1116  ASSERT_EQ(mgd.eqm_log10K(mgd.eqm_species_index["Calcite"], 0), 2.0683);
1117 }
1118 
1124 TEST(PertinentGeochemicalSystemTest, stoichiometry2)
1125 {
1126  GeochemicalDatabaseReader database("database/moose_testdb.json");
1127 
1128  // The following system has secondary species: CO2(aq), CO3--, OH-, (O-phth)--, CH4(aq), Fe+++,
1129  // >(s)FeO-
1130  PertinentGeochemicalSystem model(database,
1131  {"H2O", "H+", ">(s)FeOH", ">(w)FeOH", "Fe++", "HCO3-", "O2(aq)"},
1132  {"Fe(OH)3(ppd)fake"},
1133  {"CH4(g)fake"},
1134  {"Fe(OH)3(ppd)"},
1135  {},
1136  {},
1137  "O2(aq)",
1138  "e-");
1139  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
1140 
1141  ASSERT_EQ(mgd.eqm_species_index.size(), (std::size_t)9);
1142  for (const auto & sp : {"CO2(aq)",
1143  "CO3--",
1144  "OH-",
1145  "(O-phth)--",
1146  "CH4(aq)",
1147  "Fe+++",
1148  ">(s)FeO-",
1149  "Fe(OH)3(ppd)fake",
1150  "CH4(g)fake"})
1151  ASSERT_EQ(mgd.eqm_species_index.count(sp), (std::size_t)1);
1152 
1153  std::map<std::string, DenseMatrix<Real>> stoi_gold;
1154  for (const auto & sp : {"CO2(aq)",
1155  "CO3--",
1156  "OH-",
1157  "(O-phth)--",
1158  "CH4(aq)",
1159  "Fe+++",
1160  ">(s)FeO-",
1161  "Fe(OH)3(ppd)fake",
1162  "Fe(OH)3(ppd)",
1163  "CH4(g)fake"})
1164  stoi_gold[sp] = DenseMatrix<Real>(1, 7);
1165  // remember the order of primaries:
1166  // {"H2O", "H+", ">(s)FeOH", ">(w)FeOH", "Fe++", "HCO3-", "O2(aq)"}
1167  stoi_gold["CO2(aq)"](0, 0) = -1;
1168  stoi_gold["CO2(aq)"](0, 1) = 1;
1169  stoi_gold["CO2(aq)"](0, 5) = 1;
1170  stoi_gold["CO3--"](0, 5) = 1;
1171  stoi_gold["CO3--"](0, 1) = -1;
1172  stoi_gold["OH-"](0, 0) = 1;
1173  stoi_gold["OH-"](0, 1) = -1;
1174  stoi_gold["(O-phth)--"](0, 0) = -5;
1175  stoi_gold["(O-phth)--"](0, 5) = 8;
1176  stoi_gold["(O-phth)--"](0, 1) = 6;
1177  stoi_gold["(O-phth)--"](0, 6) = -7.5;
1178  stoi_gold["CH4(aq)"](0, 0) = 1;
1179  stoi_gold["CH4(aq)"](0, 1) = 1;
1180  stoi_gold["CH4(aq)"](0, 5) = 1;
1181  stoi_gold["CH4(aq)"](0, 6) = -2;
1182  stoi_gold["Fe+++"](0, 0) = -0.5;
1183  stoi_gold["Fe+++"](0, 4) = 1;
1184  stoi_gold["Fe+++"](0, 1) = 1;
1185  stoi_gold["Fe+++"](0, 6) = 0.25;
1186  stoi_gold[">(s)FeO-"](0, 2) = 1;
1187  stoi_gold[">(s)FeO-"](0, 1) = -1;
1188  stoi_gold["Fe(OH)3(ppd)fake"](0, 1) = -1;
1189  stoi_gold["Fe(OH)3(ppd)fake"](0, 0) = 2;
1190  stoi_gold["Fe(OH)3(ppd)fake"](0, 4) = 2;
1191  stoi_gold["Fe(OH)3(ppd)fake"](0, 6) = 0.5;
1192  stoi_gold["CH4(g)fake"](0, 0) = 3;
1193  stoi_gold["CH4(g)fake"](0, 4) = -2;
1194  stoi_gold["CH4(g)fake"](0, 5) = 3.5;
1195  stoi_gold["CH4(g)fake"](0, 6) = -4.5;
1196  stoi_gold["Fe(OH)3(ppd)"](0, 1) = -2;
1197  stoi_gold["Fe(OH)3(ppd)"](0, 4) = 1;
1198  stoi_gold["Fe(OH)3(ppd)"](0, 0) = 2.5;
1199  stoi_gold["Fe(OH)3(ppd)"](0, 6) = 0.25;
1200 
1201  for (const auto & sp : {"CO2(aq)",
1202  "CO3--",
1203  "OH-",
1204  "(O-phth)--",
1205  "CH4(aq)",
1206  "Fe+++",
1207  ">(s)FeO-",
1208  "Fe(OH)3(ppd)fake",
1209  "CH4(g)fake"})
1210  {
1211  const unsigned row = mgd.eqm_species_index[sp];
1212  ASSERT_EQ(mgd.eqm_stoichiometry.sub_matrix(row, 1, 0, 7), stoi_gold[sp]);
1213  }
1214  ASSERT_EQ(mgd.kin_stoichiometry.sub_matrix(0, 1, 0, 7), stoi_gold["Fe(OH)3(ppd)"]);
1215 }
1216 
1222 TEST(PertinentGeochemicalSystemTest, stoichiometry3)
1223 {
1224  GeochemicalDatabaseReader database("database/moose_testdb.json");
1225 
1226  // The following system has secondary species: CO2(aq), CO3--, OH-, CH4(aq), Fe+++
1227  PertinentGeochemicalSystem model(database,
1228  {"H2O", "H+", ">(s)FeOH", ">(w)FeOH", "Fe++", "HCO3-", "O2(aq)"},
1229  {},
1230  {"CH4(g)fake"},
1231  {"Fe(OH)3(ppd)", "Fe(OH)3(ppd)fake"},
1232  {"(O-phth)--"},
1233  {">(s)FeO-"},
1234  "O2(aq)",
1235  "e-");
1236  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
1237 
1238  ASSERT_EQ(mgd.eqm_species_index.size(), (std::size_t)6);
1239  for (const auto & sp : {"CO2(aq)", "CO3--", "OH-", "CH4(aq)", "Fe+++", "CH4(g)fake"})
1240  ASSERT_EQ(mgd.eqm_species_index.count(sp), (std::size_t)1);
1241 
1242  std::map<std::string, DenseMatrix<Real>> stoi_gold;
1243  for (const auto & sp : {"CO2(aq)",
1244  "CO3--",
1245  "OH-",
1246  "(O-phth)--",
1247  "CH4(aq)",
1248  "Fe+++",
1249  ">(s)FeO-",
1250  "Fe(OH)3(ppd)fake",
1251  "Fe(OH)3(ppd)",
1252  "CH4(g)fake"})
1253  stoi_gold[sp] = DenseMatrix<Real>(1, 7);
1254  // remember the order of primaries:
1255  // {"H2O", "H+", ">(s)FeOH", ">(w)FeOH", "Fe++", "HCO3-", "O2(aq)"}
1256  stoi_gold["CO2(aq)"](0, 0) = -1;
1257  stoi_gold["CO2(aq)"](0, 1) = 1;
1258  stoi_gold["CO2(aq)"](0, 5) = 1;
1259  stoi_gold["CO3--"](0, 5) = 1;
1260  stoi_gold["CO3--"](0, 1) = -1;
1261  stoi_gold["OH-"](0, 0) = 1;
1262  stoi_gold["OH-"](0, 1) = -1;
1263  stoi_gold["(O-phth)--"](0, 0) = -5;
1264  stoi_gold["(O-phth)--"](0, 5) = 8;
1265  stoi_gold["(O-phth)--"](0, 1) = 6;
1266  stoi_gold["(O-phth)--"](0, 6) = -7.5;
1267  stoi_gold["CH4(aq)"](0, 0) = 1;
1268  stoi_gold["CH4(aq)"](0, 1) = 1;
1269  stoi_gold["CH4(aq)"](0, 5) = 1;
1270  stoi_gold["CH4(aq)"](0, 6) = -2;
1271  stoi_gold["Fe+++"](0, 0) = -0.5;
1272  stoi_gold["Fe+++"](0, 4) = 1;
1273  stoi_gold["Fe+++"](0, 1) = 1;
1274  stoi_gold["Fe+++"](0, 6) = 0.25;
1275  stoi_gold[">(s)FeO-"](0, 2) = 1;
1276  stoi_gold[">(s)FeO-"](0, 1) = -1;
1277  stoi_gold["Fe(OH)3(ppd)fake"](0, 1) = -1;
1278  stoi_gold["Fe(OH)3(ppd)fake"](0, 0) = 2;
1279  stoi_gold["Fe(OH)3(ppd)fake"](0, 4) = 2;
1280  stoi_gold["Fe(OH)3(ppd)fake"](0, 6) = 0.5;
1281  stoi_gold["CH4(g)fake"](0, 0) = 3;
1282  stoi_gold["CH4(g)fake"](0, 4) = -2;
1283  stoi_gold["CH4(g)fake"](0, 5) = 3.5;
1284  stoi_gold["CH4(g)fake"](0, 6) = -4.5;
1285  stoi_gold["Fe(OH)3(ppd)"](0, 1) = -2;
1286  stoi_gold["Fe(OH)3(ppd)"](0, 4) = 1;
1287  stoi_gold["Fe(OH)3(ppd)"](0, 0) = 2.5;
1288  stoi_gold["Fe(OH)3(ppd)"](0, 6) = 0.25;
1289 
1290  for (const auto & sp : {"CO2(aq)", "CO3--", "OH-", "CH4(aq)", "Fe+++", "CH4(g)fake"})
1291  {
1292  const unsigned row = mgd.eqm_species_index[sp];
1293  ASSERT_EQ(mgd.eqm_stoichiometry.sub_matrix(row, 1, 0, 7), stoi_gold[sp]);
1294  }
1295 
1296  for (const auto & sp : {"Fe(OH)3(ppd)", "Fe(OH)3(ppd)fake", "(O-phth)--", ">(s)FeO-"})
1297  {
1298  const unsigned row = mgd.kin_species_index[sp];
1299  ASSERT_EQ(mgd.kin_stoichiometry.sub_matrix(row, 1, 0, 7), stoi_gold[sp]);
1300  }
1301 }
1302 
1308 TEST(PertinentGeochemicalSystemTest, log10K2)
1309 {
1310  GeochemicalDatabaseReader database("database/moose_testdb.json");
1311 
1312  // The following system has secondary species: CO2(aq), CO3--, OH-, (O-phth)--, CH4(aq), Fe+++,
1313  // >(s)FeO-
1314  PertinentGeochemicalSystem model(database,
1315  {"H2O", "H+", ">(s)FeOH", ">(w)FeOH", "Fe++", "HCO3-", "O2(aq)"},
1316  {"Fe(OH)3(ppd)fake"},
1317  {"CH4(g)fake"},
1318  {},
1319  {},
1320  {},
1321  "O2(aq)",
1322  "e-");
1323  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
1324 
1325  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index["CO2(aq)"], 0), -6.5570, eps);
1326  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index["CO3--"], 0), 10.6169, eps);
1327  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index["OH-"], 0), 14.9325, eps);
1328  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index["(O-phth)--"], 0), 594.3211, eps);
1329  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index["CH4(aq)"], 0), 157.8920, eps);
1330  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index["Fe+++"], 0), -10.0553, eps);
1331  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index[">(s)FeO-"], 0), 8.93, eps);
1332  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index[">(s)FeO-"], 1), 8.93 - 0.3 * (25 - 0), eps);
1333  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index[">(s)FeO-"], 2), 8.93 - 0.3 * (60 - 0), eps);
1334  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index[">(s)FeO-"], 3), 8.93 - 0.3 * (100 - 0), eps);
1335  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index[">(s)FeO-"], 4), 8.93 - 0.3 * (150 - 0), eps);
1336  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index[">(s)FeO-"], 5), 8.93 - 0.3 * (200 - 0), eps);
1337  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index[">(s)FeO-"], 6), 8.93 - 0.3 * (250 - 0), eps);
1338  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index[">(s)FeO-"], 7), 8.93 - 0.3 * (300 - 0), eps);
1339  ASSERT_NEAR(
1340  mgd.eqm_log10K(mgd.eqm_species_index["Fe(OH)3(ppd)fake"], 0), 6.1946 + 2 * (-10.0553), eps);
1341  ASSERT_NEAR(
1342  mgd.eqm_log10K(mgd.eqm_species_index["Fe(OH)3(ppd)fake"], 1), 4.8890 + 2 * (-8.4878), eps);
1343  ASSERT_NEAR(
1344  mgd.eqm_log10K(mgd.eqm_species_index["Fe(OH)3(ppd)fake"], 2), 3.4608 + 2 * (-6.6954), eps);
1345  ASSERT_NEAR(
1346  mgd.eqm_log10K(mgd.eqm_species_index["Fe(OH)3(ppd)fake"], 3), 2.2392 + 2 * (-5.0568), eps);
1347  ASSERT_NEAR(
1348  mgd.eqm_log10K(mgd.eqm_species_index["Fe(OH)3(ppd)fake"], 4), 1.1150 + 2 * (-3.4154), eps);
1349  ASSERT_NEAR(
1350  mgd.eqm_log10K(mgd.eqm_species_index["Fe(OH)3(ppd)fake"], 5), 0.2446 + 2 * (-2.0747), eps);
1351  ASSERT_NEAR(
1352  mgd.eqm_log10K(mgd.eqm_species_index["Fe(OH)3(ppd)fake"], 6), -0.5504 + 2 * (-0.8908), eps);
1353  ASSERT_NEAR(
1354  mgd.eqm_log10K(mgd.eqm_species_index["Fe(OH)3(ppd)fake"], 7), -1.5398 + 2 * (0.2679), eps);
1355  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index["CH4(g)fake"], 0),
1356  -2.6487 + 2 * 157.8920 - 2 * (-10.0553),
1357  eps);
1358  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index["CH4(g)fake"], 1),
1359  -2.8202 + 2 * 144.1080 - 2 * (-8.4878),
1360  eps);
1361  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index["CH4(g)fake"], 2),
1362  -2.9329 + 2 * 127.9360 - 2 * (-6.6954),
1363  eps);
1364  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index["CH4(g)fake"], 3),
1365  -2.9446 + 2 * 112.8800 - 2 * (-5.0568),
1366  eps);
1367  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index["CH4(g)fake"], 4),
1368  -2.9163 + 2 * 97.7060 - 2 * (-3.4154),
1369  eps);
1370  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index["CH4(g)fake"], 5),
1371  -2.7253 + 2 * 85.2880 - 2 * (-2.0747),
1372  eps);
1373  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index["CH4(g)fake"], 6),
1374  -2.4643 + 2 * 74.7500 - 2 * (-0.8908),
1375  eps);
1376  ASSERT_NEAR(mgd.eqm_log10K(mgd.eqm_species_index["CH4(g)fake"], 7),
1377  -2.1569 + 2 * 65.6500 - 2 * (0.2679),
1378  eps);
1379 }
1380 
1385 TEST(PertinentGeochemicalSystemTest, secondarySpecies2)
1386 {
1387  GeochemicalDatabaseReader database("database/moose_testdb.json");
1388 
1389  // The following system has secondary species: CO2(aq), CO3--, CaCO3, CaOH+, OH-, >(s)FeO-,
1390  // >(s)FeOCa+
1391  PertinentGeochemicalSystem model(
1392  database,
1393  {"H2O", "H+", "HCO3-", "O2(aq)", "Ca++", ">(s)FeOH", "(O-phth)--"},
1394  {"Calcite"},
1395  {},
1396  {},
1397  {"CH4(aq)"},
1398  {},
1399  "O2(aq)",
1400  "e-");
1401  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
1402 
1403  ASSERT_EQ(mgd.eqm_species_index.size(), (std::size_t)8);
1404  for (const auto & sp :
1405  {"CO2(aq)", "CO3--", "CaCO3", "CaOH+", "OH-", ">(s)FeO-", ">(s)FeOCa+", "Calcite"})
1406  ASSERT_EQ(mgd.eqm_species_index.count(sp), (std::size_t)1);
1407 }
1408 
1410 TEST(PertinentGeochemicalSystemTest, isMineral)
1411 {
1412  GeochemicalDatabaseReader database("database/moose_testdb.json");
1413 
1414  // The following system has secondary species: CO2(aq), CO3--, OH-, (O-phth)--, CH4(aq), Fe+++,
1415  // >(s)FeO-, e-
1416  PertinentGeochemicalSystem model(database,
1417  {"H2O", "H+", ">(s)FeOH", ">(w)FeOH", "Fe++", "HCO3-", "O2(aq)"},
1418  {"Fe(OH)3(ppd)fake"},
1419  {"CH4(g)fake"},
1420  {},
1421  {},
1422  {},
1423  "O2(aq)",
1424  "e-");
1425  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
1426 
1427  for (const auto & species : mgd.basis_species_index)
1428  ASSERT_EQ(mgd.basis_species_mineral[species.second], false);
1429  for (const auto & species : mgd.eqm_species_index)
1430  if (species.first == "Fe(OH)3(ppd)fake")
1431  ASSERT_EQ(mgd.eqm_species_mineral[species.second], true);
1432  else
1433  ASSERT_EQ(mgd.eqm_species_mineral[species.second], false);
1434 
1435  PertinentGeochemicalSystem model2(
1436  database,
1437  {"H2O", "H+", ">(s)FeOH", ">(w)FeOH", "Fe++", "HCO3-", "O2(aq)"},
1438  {},
1439  {"CH4(g)fake"},
1440  {"Fe(OH)3(ppd)fake"},
1441  {"(O-phth)--"},
1442  {">(s)FeO-"},
1443  "O2(aq)",
1444  "e-");
1445  ModelGeochemicalDatabase mgd2 = model2.modelGeochemicalDatabase();
1446 
1447  for (const auto & species : mgd2.basis_species_index)
1448  ASSERT_EQ(mgd2.basis_species_mineral[species.second], false);
1449  for (const auto & species : mgd2.eqm_species_index)
1450  ASSERT_EQ(mgd2.eqm_species_mineral[species.second], false);
1451  for (const auto & species : mgd2.kin_species_index)
1452  if (species.first == "Fe(OH)3(ppd)fake")
1453  ASSERT_EQ(mgd2.kin_species_mineral[species.second], true);
1454  else
1455  ASSERT_EQ(mgd2.kin_species_mineral[species.second], false);
1456 }
1457 
1459 TEST(PertinentGeochemicalSystemTest, isGas)
1460 {
1461  GeochemicalDatabaseReader database("database/moose_testdb.json");
1462 
1463  // The following system has secondary species: CO2(aq), CO3--, OH-, (O-phth)--, CH4(aq), Fe+++,
1464  // >(s)FeO-, e-
1465  PertinentGeochemicalSystem model(database,
1466  {"H2O", "H+", ">(s)FeOH", ">(w)FeOH", "Fe++", "HCO3-", "O2(aq)"},
1467  {"Fe(OH)3(ppd)fake"},
1468  {"CH4(g)fake", "O2(g)"},
1469  {},
1470  {},
1471  {},
1472  "O2(aq)",
1473  "e-");
1474  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
1475 
1476  for (const auto & species : mgd.basis_species_index)
1477  ASSERT_EQ(mgd.basis_species_gas[species.second], false);
1478  for (const auto & species : mgd.eqm_species_index)
1479  if (species.first == "CH4(g)fake" || species.first == "O2(g)")
1480  ASSERT_EQ(mgd.eqm_species_gas[species.second], true);
1481  else
1482  ASSERT_EQ(mgd.eqm_species_gas[species.second], false);
1483 }
1484 
1486 TEST(PertinentGeochemicalSystemTest, isTransported)
1487 {
1488  GeochemicalDatabaseReader database("database/moose_testdb.json");
1489 
1490  // The following system has secondary species: CO2(aq), CO3--, OH-, (O-phth)--, CH4(aq), Fe+++,
1491  // >(s)FeO-, e-
1492  PertinentGeochemicalSystem model(database,
1493  {"H2O", "H+", ">(s)FeOH", ">(w)FeOH", "Fe++", "HCO3-", "O2(aq)"},
1494  {"Fe(OH)3(ppd)"},
1495  {"CH4(g)fake"},
1496  {},
1497  {},
1498  {},
1499  "O2(aq)",
1500  "e-");
1501  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
1502 
1503  for (const auto & species : mgd.basis_species_index)
1504  if (species.first == ">(s)FeOH" || species.first == ">(w)FeOH")
1505  EXPECT_EQ(mgd.basis_species_transported[species.second], false);
1506  else
1507  EXPECT_EQ(mgd.basis_species_transported[species.second], true);
1508  for (const auto & species : mgd.eqm_species_index)
1509  if (species.first == "Fe(OH)3(ppd)" || species.first == ">(s)FeO-")
1510  EXPECT_EQ(mgd.eqm_species_transported[species.second], false);
1511  else
1512  EXPECT_EQ(mgd.eqm_species_transported[species.second], true);
1513 
1514  PertinentGeochemicalSystem model2(
1515  database,
1516  {"H2O", "H+", ">(s)FeOH", ">(w)FeOH", "Fe++", "HCO3-", "O2(aq)"},
1517  {},
1518  {"CH4(g)fake"},
1519  {"Fe(OH)3(ppd)"},
1520  {"(O-phth)--"},
1521  {">(s)FeO-"},
1522  "O2(aq)",
1523  "e-");
1524  ModelGeochemicalDatabase mgd2 = model2.modelGeochemicalDatabase();
1525 
1526  for (const auto & species : mgd2.basis_species_index)
1527  if (species.first == ">(s)FeOH" || species.first == ">(w)FeOH")
1528  ASSERT_EQ(mgd2.basis_species_transported[species.second], false);
1529  else
1530  ASSERT_EQ(mgd2.basis_species_transported[species.second], true);
1531  for (const auto & species : mgd2.eqm_species_index)
1532  ASSERT_EQ(mgd2.eqm_species_transported[species.second], true);
1533  for (const auto & species : mgd2.kin_species_index)
1534  if (species.first == "Fe(OH)3(ppd)" || species.first == ">(s)FeO-")
1535  ASSERT_EQ(mgd2.kin_species_transported[species.second], false);
1536  else
1537  ASSERT_EQ(mgd2.kin_species_transported[species.second], true);
1538 }
1539 
1541 TEST(PertinentGeochemicalSystemTest, redoxCapture)
1542 {
1543  GeochemicalDatabaseReader database("database/moose_testdb.json");
1544 
1545  PertinentGeochemicalSystem model_no_redox(
1546  database,
1547  {"H2O", "H+", ">(s)FeOH", ">(w)FeOH", "Fe++", "HCO3-", "O2(aq)"},
1548  {"Fe(OH)3(ppd)fake"},
1549  {"CH4(g)fake", "O2(g)"},
1550  {},
1551  {},
1552  {},
1553  "O2(aq)",
1554  "e-");
1555  ModelGeochemicalDatabase mgd_no_redox = model_no_redox.modelGeochemicalDatabase();
1556 
1557  EXPECT_EQ(mgd_no_redox.redox_stoichiometry.m(), (unsigned)1);
1558  EXPECT_EQ(mgd_no_redox.redox_log10K.m(), (unsigned)1);
1559 
1560  PertinentGeochemicalSystem model_redox(
1561  database,
1562  {"H2O", "H+", "HCO3-", "O2(aq)", "(O-phth)--", "CH4(aq)", "StoiCheckRedox", "Fe+++"},
1563  {"Fe(OH)3(ppd)fake"},
1564  {"CH4(g)fake", "O2(g)"},
1565  {},
1566  {},
1567  {},
1568  "O2(aq)",
1569  "e-");
1570  ModelGeochemicalDatabase mgd_redox = model_redox.modelGeochemicalDatabase();
1571 
1572  EXPECT_EQ(mgd_redox.redox_lhs, "e-");
1573 
1574  // StoiCheckRedox is not expressed in terms of O2(aq), and there is no Fe++ so Fe+++ does not
1575  // have a pair
1576  EXPECT_EQ(mgd_redox.redox_stoichiometry.m(), (unsigned)3);
1577  EXPECT_EQ(mgd_redox.redox_log10K.m(), (unsigned)3);
1578 
1579  // not sure which order the redox has been ordered in. The reactions are:
1580  // e- = (1/4/7.5)(O-phth)-- + (1/2 + 5/4/7.5)H2O + (-1 - 6/4/7.5)H+ - 8/4/7.5HCO3-
1581  // e- = (1/8)CH4(aq) + (1/2 - 1/8)H2O - (1+1/8)H+ - (1/8)HCO3-
1582  const bool ophth_is_slot_one = (mgd_redox.redox_stoichiometry(1, 4) > 1.0E-6);
1583  const unsigned ophth_slot = (ophth_is_slot_one ? 1 : 2);
1584  const unsigned ch4_slot = (ophth_is_slot_one ? 2 : 1);
1585 
1586  // e- = (1/4/7.5)(O-phth)-- + (1/2 + 5/4/7.5)H2O + (-1 - 6/4/7.5)H+ - 8/4/7.5HCO3-
1587  Real boa = 1.0 / 4.0 / 7.5;
1588  EXPECT_EQ(mgd_redox.redox_stoichiometry(ophth_slot, 0), 0.5 + 5.0 * boa);
1589  EXPECT_EQ(mgd_redox.redox_stoichiometry(ophth_slot, 1), -1 - 6.0 * boa);
1590  EXPECT_EQ(mgd_redox.redox_stoichiometry(ophth_slot, 2), -8.0 * boa);
1591  EXPECT_EQ(mgd_redox.redox_stoichiometry(ophth_slot, 3), 0.0);
1592  EXPECT_EQ(mgd_redox.redox_stoichiometry(ophth_slot, 4), boa);
1593  EXPECT_EQ(mgd_redox.redox_stoichiometry(ophth_slot, 5), 0.0);
1594  EXPECT_EQ(mgd_redox.redox_stoichiometry(ophth_slot, 6), 0.0);
1595  EXPECT_EQ(mgd_redox.redox_stoichiometry(ophth_slot, 7), 0.0);
1596  EXPECT_NEAR(
1597  mgd_redox.redox_log10K(ophth_slot, 0), -boa * 594.3211 + 22.76135 - 0.25 * (-2.6610), 1E-8);
1598  EXPECT_NEAR(
1599  mgd_redox.redox_log10K(ophth_slot, 1), -boa * 542.8292 + 20.7757 - 0.25 * (-2.8990), 1E-8);
1600  EXPECT_NEAR(
1601  mgd_redox.redox_log10K(ophth_slot, 2), -boa * 482.3612 + 18.513025 - 0.25 * (-3.0580), 1E-8);
1602  EXPECT_NEAR(
1603  mgd_redox.redox_log10K(ophth_slot, 3), -boa * 425.9738 + 16.4658 - 0.25 * (-3.1250), 1E-8);
1604  EXPECT_NEAR(
1605  mgd_redox.redox_log10K(ophth_slot, 4), -boa * 368.7004 + 14.473225 - 0.25 * (-3.0630), 1E-8);
1606  EXPECT_NEAR(
1607  mgd_redox.redox_log10K(ophth_slot, 5), -boa * 321.8658 + 12.92125 - 0.25 * (-2.9140), 1E-8);
1608  EXPECT_NEAR(
1609  mgd_redox.redox_log10K(ophth_slot, 6), -boa * 281.8216 + 11.68165 - 0.25 * (-2.6600), 1E-8);
1610  EXPECT_NEAR(
1611  mgd_redox.redox_log10K(ophth_slot, 7), -boa * 246.4849 + 10.67105 - 0.25 * (-2.4100), 1E-8);
1612 
1613  // e- = (1/8)CH4(aq) + (1/2 - 1/8)H2O - (1+1/8)H+ - (1/8)HCO3-
1614  boa = 1.0 / 8.0;
1615  EXPECT_EQ(mgd_redox.redox_stoichiometry(ch4_slot, 0), 0.5 - boa);
1616  EXPECT_EQ(mgd_redox.redox_stoichiometry(ch4_slot, 1), -1 - boa);
1617  EXPECT_EQ(mgd_redox.redox_stoichiometry(ch4_slot, 2), -boa);
1618  EXPECT_EQ(mgd_redox.redox_stoichiometry(ch4_slot, 3), 0.0);
1619  EXPECT_EQ(mgd_redox.redox_stoichiometry(ch4_slot, 4), 0.0);
1620  EXPECT_EQ(mgd_redox.redox_stoichiometry(ch4_slot, 5), boa);
1621  EXPECT_EQ(mgd_redox.redox_stoichiometry(ch4_slot, 6), 0.0);
1622  EXPECT_EQ(mgd_redox.redox_stoichiometry(ch4_slot, 7), 0.0);
1623  EXPECT_NEAR(
1624  mgd_redox.redox_log10K(ch4_slot, 0), -boa * 157.8920 + 22.76135 - 0.25 * (-2.6610), 1E-8);
1625  EXPECT_NEAR(
1626  mgd_redox.redox_log10K(ch4_slot, 1), -boa * 144.1080 + 20.7757 - 0.25 * (-2.8990), 1E-8);
1627  EXPECT_NEAR(
1628  mgd_redox.redox_log10K(ch4_slot, 2), -boa * 127.9360 + 18.513025 - 0.25 * (-3.0580), 1E-8);
1629  EXPECT_NEAR(
1630  mgd_redox.redox_log10K(ch4_slot, 3), -boa * 112.8800 + 16.4658 - 0.25 * (-3.1250), 1E-8);
1631  EXPECT_NEAR(
1632  mgd_redox.redox_log10K(ch4_slot, 4), -boa * 97.7060 + 14.473225 - 0.25 * (-3.0630), 1E-8);
1633  EXPECT_NEAR(
1634  mgd_redox.redox_log10K(ch4_slot, 5), -boa * 85.2880 + 12.92125 - 0.25 * (-2.9140), 1E-8);
1635  EXPECT_NEAR(
1636  mgd_redox.redox_log10K(ch4_slot, 6), -boa * 74.7500 + 11.68165 - 0.25 * (-2.6600), 1E-8);
1637  EXPECT_NEAR(
1638  mgd_redox.redox_log10K(ch4_slot, 7), -boa * 65.6500 + 10.67105 - 0.25 * (-2.4100), 1E-8);
1639 }
1640 
1642 TEST(PertinentGeochemicalSystemTest, redoxCapture_db_strange)
1643 {
1644  GeochemicalDatabaseReader database("database/moose_testdb_e_strange.json");
1645 
1646  PertinentGeochemicalSystem model_redox(
1647  database, {"H2O", "H+", "HCO3-", "O2(aq)", "(O-phth)--"}, {}, {}, {}, {}, {}, "O2(aq)", "e-");
1648  ModelGeochemicalDatabase mgd_redox = model_redox.modelGeochemicalDatabase();
1649 
1650  EXPECT_EQ(mgd_redox.redox_lhs, "e-");
1651 
1652  EXPECT_EQ(mgd_redox.redox_stoichiometry.m(), (unsigned)1);
1653  EXPECT_EQ(mgd_redox.redox_log10K.m(), (unsigned)1);
1654 
1655  // Check the reaction:
1656  // e- = H2O + 3O2(aq) + OH- = 2H2O + 3O2(aq) - H+
1657  EXPECT_EQ(mgd_redox.redox_stoichiometry(0, 0), 2.0);
1658  EXPECT_EQ(mgd_redox.redox_stoichiometry(0, 1), -1.0);
1659  EXPECT_EQ(mgd_redox.redox_stoichiometry(0, 2), 0.0);
1660  EXPECT_EQ(mgd_redox.redox_stoichiometry(0, 3), 3.0);
1661  EXPECT_EQ(mgd_redox.redox_stoichiometry(0, 4), 0.0);
1662  EXPECT_NEAR(mgd_redox.redox_log10K(0, 0), 22.76135 + 14.9325, 1E-8);
1663  EXPECT_NEAR(mgd_redox.redox_log10K(0, 1), 20.7757 + 13.9868, 1E-8);
1664  EXPECT_NEAR(mgd_redox.redox_log10K(0, 2), 18.513025 + 13.0199, 1E-8);
1665  EXPECT_NEAR(mgd_redox.redox_log10K(0, 3), 16.4658 + 12.2403, 1E-8);
1666  EXPECT_NEAR(mgd_redox.redox_log10K(0, 4), 14.473225 + 11.5940, 1E-8);
1667  EXPECT_NEAR(mgd_redox.redox_log10K(0, 5), 12.92125 + 11.2191, 1E-8);
1668  EXPECT_NEAR(mgd_redox.redox_log10K(0, 6), 11.68165 + 11.0880, 1E-8);
1669  EXPECT_NEAR(mgd_redox.redox_log10K(0, 7), 10.67105 + 1001.2844, 1E-8);
1670 }
1671 
1673 TEST(PertinentGeochemicalSystemTest, redoxCapture_redox_o_strange)
1674 {
1675  GeochemicalDatabaseReader database("database/moose_testdb_e_strange.json");
1676 
1677  PertinentGeochemicalSystem model_redox(
1678  database, {"H2O", "H+", "HCO3-", "O2(aq)", "(O-phth)_0"}, {}, {}, {}, {}, {}, "O2(aq)", "e-");
1679  ModelGeochemicalDatabase mgd_redox = model_redox.modelGeochemicalDatabase();
1680 
1681  EXPECT_EQ(mgd_redox.redox_lhs, "e-");
1682 
1683  EXPECT_EQ(mgd_redox.redox_stoichiometry.m(), (unsigned)1);
1684  EXPECT_EQ(mgd_redox.redox_log10K.m(), (unsigned)1);
1685 
1686  // Check the reaction:
1687  // e- = H2O + 3O2(aq) + OH- = 2H2O + 3O2(aq) - H+
1688  EXPECT_EQ(mgd_redox.redox_stoichiometry(0, 0), 2.0);
1689  EXPECT_EQ(mgd_redox.redox_stoichiometry(0, 1), -1.0);
1690  EXPECT_EQ(mgd_redox.redox_stoichiometry(0, 2), 0.0);
1691  EXPECT_EQ(mgd_redox.redox_stoichiometry(0, 3), 3.0);
1692  EXPECT_EQ(mgd_redox.redox_stoichiometry(0, 4), 0.0);
1693  EXPECT_NEAR(mgd_redox.redox_log10K(0, 0), 22.76135 + 14.9325, 1E-8);
1694  EXPECT_NEAR(mgd_redox.redox_log10K(0, 1), 20.7757 + 13.9868, 1E-8);
1695  EXPECT_NEAR(mgd_redox.redox_log10K(0, 2), 18.513025 + 13.0199, 1E-8);
1696  EXPECT_NEAR(mgd_redox.redox_log10K(0, 3), 16.4658 + 12.2403, 1E-8);
1697  EXPECT_NEAR(mgd_redox.redox_log10K(0, 4), 14.473225 + 11.5940, 1E-8);
1698  EXPECT_NEAR(mgd_redox.redox_log10K(0, 5), 12.92125 + 11.2191, 1E-8);
1699  EXPECT_NEAR(mgd_redox.redox_log10K(0, 6), 11.68165 + 11.0880, 1E-8);
1700  EXPECT_NEAR(mgd_redox.redox_log10K(0, 7), 10.67105 + 1001.2844, 1E-8);
1701 }
1702 
1704 TEST(PertinentGeochemicalSystemTest, addKineticRateExceptions)
1705 {
1706  GeochemicalDatabaseReader database("database/moose_testdb.json");
1707  PertinentGeochemicalSystem model(database,
1708  {"H2O", "H+", "HCO3-", "O2(aq)", "Ca++"},
1709  {},
1710  {},
1711  {"Calcite"},
1712  {"CH4(aq)"},
1713  {},
1714  "O2(aq)",
1715  "e-");
1716 
1717  try
1718  {
1719  KineticRateUserDescription rate("Ca++",
1720  1.0,
1721  2.0,
1722  true,
1723  0.0,
1724  0.0,
1725  0.0,
1726  {"H2O"},
1727  {3.0},
1728  {0.0},
1729  {0.0},
1730  4.0,
1731  5.0,
1732  6.0,
1733  7.0,
1734  DirectionChoiceEnum::BOTH,
1735  "H2O",
1736  0.0,
1737  -1.0,
1738  0.0);
1739  model.addKineticRate(rate);
1740  FAIL() << "Missing expected exception.";
1741  }
1742  catch (const std::exception & e)
1743  {
1744  std::string msg(e.what());
1745  ASSERT_TRUE(
1746  msg.find(
1747  "Cannot prescribe a kinetic rate to species Ca++ since it is not a kinetic species") !=
1748  std::string::npos)
1749  << "Failed with unexpected error message: " << msg;
1750  }
1751 
1752  try
1753  {
1754  KineticRateUserDescription rate("CH4(aq)",
1755  1.0,
1756  2.0,
1757  true,
1758  0.0,
1759  0.0,
1760  0.0,
1761  {"H2O", "H++"},
1762  {3.0, 1.0},
1763  {0.0, 0.0},
1764  {0.0, 0.0},
1765  4.0,
1766  5.0,
1767  6.0,
1768  7.0,
1769  DirectionChoiceEnum::BOTH,
1770  "H2O",
1771  0.0,
1772  -1.0,
1773  0.0);
1774  model.addKineticRate(rate);
1775  FAIL() << "Missing expected exception.";
1776  }
1777  catch (const std::exception & e)
1778  {
1779  std::string msg(e.what());
1780  ASSERT_TRUE(msg.find("Promoting species H++ must be a basis or a secondary species") !=
1781  std::string::npos)
1782  << "Failed with unexpected error message: " << msg;
1783  }
1784 
1785  try
1786  {
1787  KineticRateUserDescription rate("CH4(aq)",
1788  1.0,
1789  2.0,
1790  true,
1791  0.0,
1792  0.0,
1793  0.0,
1794  {},
1795  {},
1796  {},
1797  {},
1798  4.0,
1799  5.0,
1800  6.0,
1801  7.0,
1802  DirectionChoiceEnum::BOTH,
1803  "CH4(aq)",
1804  0.0,
1805  -1.0,
1806  0.0);
1807  model.addKineticRate(rate);
1808  FAIL() << "Missing expected exception.";
1809  }
1810  catch (const std::exception & e)
1811  {
1812  std::string msg(e.what());
1813  ASSERT_TRUE(msg.find("Progeny CH4(aq) must be a basis or a secondary species") !=
1814  std::string::npos)
1815  << "Failed with unexpected error message: " << msg;
1816  }
1817 }
1818 
1820 TEST(PertinentGeochemicalSystemTest, addKineticRate)
1821 {
1822  GeochemicalDatabaseReader database("database/moose_testdb.json");
1823  PertinentGeochemicalSystem model(database,
1824  {"H2O", "H+", "HCO3-", "O2(aq)", "Ca++"},
1825  {},
1826  {},
1827  {"Calcite"},
1828  {"CH4(aq)"},
1829  {},
1830  "O2(aq)",
1831  "e-");
1832  KineticRateUserDescription rate("CH4(aq)",
1833  1.0,
1834  2.0,
1835  true,
1836  0.0,
1837  0.0,
1838  0.0,
1839  {"H2O", "OH-", "O2(aq)", "CO2(aq)", "CaCO3"},
1840  {3.0, 3.1, 3.2, 3.3, 3.4},
1841  {1.0, 2.0, 3.0, 4.0, 5.0},
1842  {1.25, 1.5, 1.75, 2.25, 3.25},
1843  4.0,
1844  5.0,
1845  6.0,
1846  7.0,
1847  DirectionChoiceEnum::BOTH,
1848  "H+",
1849  2.0,
1850  -1.0,
1851  0.0);
1852  model.addKineticRate(rate);
1853 
1854  const ModelGeochemicalDatabase & mgd = model.modelGeochemicalDatabase();
1855 
1856  EXPECT_EQ(mgd.kin_rate.size(), (std::size_t)1);
1857  EXPECT_EQ(mgd.kin_rate[0].kinetic_species_index, mgd.kin_species_index.at("CH4(aq)"));
1858  EXPECT_EQ(mgd.kin_rate[0].description.intrinsic_rate_constant, 1.0);
1859  EXPECT_EQ(mgd.kin_rate[0].description.area_quantity, 2.0);
1860  EXPECT_EQ(mgd.kin_rate[0].description.multiply_by_mass, true);
1861  std::vector<Real> pi_gold(mgd.basis_species_index.size() + mgd.eqm_species_index.size(), 0.0);
1862  std::vector<Real> pmi_gold(mgd.basis_species_index.size() + mgd.eqm_species_index.size(), 0.0);
1863  std::vector<Real> pmk_gold(mgd.basis_species_index.size() + mgd.eqm_species_index.size(), 0.0);
1864  EXPECT_EQ(mgd.kin_rate[0].promoting_indices.size(), pi_gold.size());
1865  EXPECT_EQ(mgd.kin_rate[0].promoting_monod_indices.size(), pmi_gold.size());
1866  EXPECT_EQ(mgd.kin_rate[0].promoting_half_saturation.size(), pmk_gold.size());
1867  pi_gold[0] = 3.0; // H2O
1868  pmi_gold[0] = 1.0;
1869  pmk_gold[0] = 1.25;
1870  pi_gold[5 + mgd.eqm_species_index.at("OH-")] = 3.1;
1871  pmi_gold[5 + mgd.eqm_species_index.at("OH-")] = 2.0;
1872  pmk_gold[5 + mgd.eqm_species_index.at("OH-")] = 1.5;
1873  pi_gold[3] = 3.2; // O2(aq)
1874  pmi_gold[3] = 3.0;
1875  pmk_gold[3] = 1.75;
1876  pi_gold[5 + mgd.eqm_species_index.at("CO2(aq)")] = 3.3;
1877  pmi_gold[5 + mgd.eqm_species_index.at("CO2(aq)")] = 4.0;
1878  pmk_gold[5 + mgd.eqm_species_index.at("CO2(aq)")] = 2.25;
1879  pi_gold[5 + mgd.eqm_species_index.at("CaCO3")] = 3.4;
1880  pmi_gold[5 + mgd.eqm_species_index.at("CaCO3")] = 5.0;
1881  pmk_gold[5 + mgd.eqm_species_index.at("CaCO3")] = 3.25;
1882  for (unsigned i = 0; i < pi_gold.size(); ++i)
1883  {
1884  EXPECT_EQ(mgd.kin_rate[0].promoting_indices[i], pi_gold[i]);
1885  EXPECT_EQ(mgd.kin_rate[0].promoting_monod_indices[i], pmi_gold[i]);
1886  EXPECT_EQ(mgd.kin_rate[0].promoting_half_saturation[i], pmk_gold[i]);
1887  }
1888  EXPECT_EQ(mgd.kin_rate[0].description.theta, 4.0);
1889  EXPECT_EQ(mgd.kin_rate[0].description.eta, 5.0);
1890  EXPECT_EQ(mgd.kin_rate[0].description.activation_energy, 6.0);
1891  EXPECT_EQ(mgd.kin_rate[0].description.one_over_T0, 7.0);
1892  EXPECT_EQ(mgd.kin_rate[0].progeny_index, mgd.basis_species_index.at("H+"));
1893  EXPECT_EQ(mgd.kin_rate[0].description.progeny, "H+");
1894  EXPECT_EQ(mgd.kin_rate[0].description.progeny_efficiency, 2.0);
1895 
1896  KineticRateUserDescription ratec("Calcite",
1897  7.0,
1898  6.0,
1899  false,
1900  0.0,
1901  0.0,
1902  0.0,
1903  {"H+"},
1904  {-3.0},
1905  {-1.0},
1906  {-2.0},
1907  5.0,
1908  4.0,
1909  3.0,
1910  2.0,
1911  DirectionChoiceEnum::BOTH,
1912  "OH-",
1913  1.125,
1914  -1.0,
1915  0.0);
1916  model.addKineticRate(ratec);
1917 
1918  EXPECT_EQ(mgd.kin_rate.size(), (std::size_t)2);
1919  EXPECT_EQ(mgd.kin_rate[1].kinetic_species_index, mgd.kin_species_index.at("Calcite"));
1920  EXPECT_EQ(mgd.kin_rate[1].description.intrinsic_rate_constant, 7.0);
1921  EXPECT_EQ(mgd.kin_rate[1].description.area_quantity, 6.0);
1922  EXPECT_EQ(mgd.kin_rate[1].description.multiply_by_mass, false);
1923  EXPECT_EQ(mgd.kin_rate[1].promoting_indices.size(), pi_gold.size());
1924  std::fill(pi_gold.begin(), pi_gold.end(), 0.0);
1925  std::fill(pmi_gold.begin(), pmi_gold.end(), 0.0);
1926  std::fill(pmk_gold.begin(), pmk_gold.end(), 0.0);
1927  pi_gold[1] = -3.0; // H++
1928  pmi_gold[1] = -1.0;
1929  pmk_gold[1] = -2.0;
1930  for (unsigned i = 0; i < pi_gold.size(); ++i)
1931  {
1932  EXPECT_EQ(mgd.kin_rate[1].promoting_indices[i], pi_gold[i]);
1933  EXPECT_EQ(mgd.kin_rate[1].promoting_monod_indices[i], pmi_gold[i]);
1934  EXPECT_EQ(mgd.kin_rate[1].promoting_half_saturation[i], pmk_gold[i]);
1935  }
1936  EXPECT_EQ(mgd.kin_rate[1].description.theta, 5.0);
1937  EXPECT_EQ(mgd.kin_rate[1].description.eta, 4.0);
1938  EXPECT_EQ(mgd.kin_rate[1].description.activation_energy, 3.0);
1939  EXPECT_EQ(mgd.kin_rate[1].description.one_over_T0, 2.0);
1940  EXPECT_EQ(mgd.kin_rate[1].progeny_index, 5 + mgd.eqm_species_index.at("OH-"));
1941  EXPECT_EQ(mgd.kin_rate[1].description.progeny, "OH-");
1942  EXPECT_EQ(mgd.kin_rate[1].description.progeny_efficiency, 1.125);
1943 }
1944 
1950 TEST(PertinentGeochemicalSystemTest, kin_log10K1)
1951 {
1952  GeochemicalDatabaseReader database("database/moose_testdb.json");
1953 
1954  PertinentGeochemicalSystem model(database,
1955  {"H2O", "H+", ">(s)FeOH", ">(w)FeOH", "Fe++", "HCO3-", "O2(aq)"},
1956  {},
1957  {"CH4(g)fake"},
1958  {"Fe(OH)3(ppd)fake"},
1959  {"(O-phth)--"},
1960  {">(s)FeO-"},
1961  "O2(aq)",
1962  "e-");
1963  const ModelGeochemicalDatabase & mgd = model.modelGeochemicalDatabase();
1964 
1965  ASSERT_NEAR(mgd.kin_log10K(mgd.kin_species_index.at("Fe(OH)3(ppd)fake"), 0),
1966  6.1946 + 2 * (-10.0553),
1967  eps);
1968  ASSERT_NEAR(
1969  mgd.kin_log10K(mgd.kin_species_index.at("Fe(OH)3(ppd)fake"), 1), 4.8890 + 2 * (-8.4878), eps);
1970  ASSERT_NEAR(
1971  mgd.kin_log10K(mgd.kin_species_index.at("Fe(OH)3(ppd)fake"), 2), 3.4608 + 2 * (-6.6954), eps);
1972  ASSERT_NEAR(
1973  mgd.kin_log10K(mgd.kin_species_index.at("Fe(OH)3(ppd)fake"), 3), 2.2392 + 2 * (-5.0568), eps);
1974  ASSERT_NEAR(
1975  mgd.kin_log10K(mgd.kin_species_index.at("Fe(OH)3(ppd)fake"), 4), 1.1150 + 2 * (-3.4154), eps);
1976  ASSERT_NEAR(
1977  mgd.kin_log10K(mgd.kin_species_index.at("Fe(OH)3(ppd)fake"), 5), 0.2446 + 2 * (-2.0747), eps);
1978  ASSERT_NEAR(mgd.kin_log10K(mgd.kin_species_index.at("Fe(OH)3(ppd)fake"), 6),
1979  -0.5504 + 2 * (-0.8908),
1980  eps);
1981  ASSERT_NEAR(
1982  mgd.kin_log10K(mgd.kin_species_index.at("Fe(OH)3(ppd)fake"), 7), -1.5398 + 2 * (0.2679), eps);
1983  EXPECT_EQ(mgd.kin_log10K(mgd.kin_species_index.at("(O-phth)--"), 0), 594.3211);
1984  EXPECT_EQ(mgd.kin_log10K(mgd.kin_species_index.at("(O-phth)--"), 1), 542.8292);
1985  EXPECT_EQ(mgd.kin_log10K(mgd.kin_species_index.at("(O-phth)--"), 2), 482.3612);
1986  EXPECT_EQ(mgd.kin_log10K(mgd.kin_species_index.at("(O-phth)--"), 3), 425.9738);
1987  EXPECT_EQ(mgd.kin_log10K(mgd.kin_species_index.at("(O-phth)--"), 4), 368.7004);
1988  EXPECT_EQ(mgd.kin_log10K(mgd.kin_species_index.at("(O-phth)--"), 5), 321.8658);
1989  EXPECT_EQ(mgd.kin_log10K(mgd.kin_species_index.at("(O-phth)--"), 6), 281.8216);
1990  EXPECT_EQ(mgd.kin_log10K(mgd.kin_species_index.at("(O-phth)--"), 7), 246.4849);
1991  ASSERT_NEAR(mgd.kin_log10K(mgd.kin_species_index.at(">(s)FeO-"), 0), 8.93, eps);
1992  ASSERT_NEAR(mgd.kin_log10K(mgd.kin_species_index.at(">(s)FeO-"), 1), 8.93 - 0.3 * (25 - 0), eps);
1993  ASSERT_NEAR(mgd.kin_log10K(mgd.kin_species_index.at(">(s)FeO-"), 2), 8.93 - 0.3 * (60 - 0), eps);
1994  ASSERT_NEAR(mgd.kin_log10K(mgd.kin_species_index.at(">(s)FeO-"), 3), 8.93 - 0.3 * (100 - 0), eps);
1995  ASSERT_NEAR(mgd.kin_log10K(mgd.kin_species_index.at(">(s)FeO-"), 4), 8.93 - 0.3 * (150 - 0), eps);
1996  ASSERT_NEAR(mgd.kin_log10K(mgd.kin_species_index.at(">(s)FeO-"), 5), 8.93 - 0.3 * (200 - 0), eps);
1997  ASSERT_NEAR(mgd.kin_log10K(mgd.kin_species_index.at(">(s)FeO-"), 6), 8.93 - 0.3 * (250 - 0), eps);
1998  ASSERT_NEAR(mgd.kin_log10K(mgd.kin_species_index.at(">(s)FeO-"), 7), 8.93 - 0.3 * (300 - 0), eps);
1999 }
2000 
2002 TEST(PertinentGeochemicalSystemTest, initSwapToOrigBasis)
2003 {
2004  GeochemicalDatabaseReader database("database/moose_testdb.json");
2005  PertinentGeochemicalSystem model(database, {"H2O", "H+"}, {}, {}, {}, {}, {}, "O2(aq)", "e-");
2006  ModelGeochemicalDatabase mgd = model.modelGeochemicalDatabase();
2007  ASSERT_EQ(mgd.swap_to_original_basis.n(), 0.0);
2008 }
2009 
2011 TEST(PertinentGeochemicalSystemTest, getIndexOfOriginalBasisSpecies)
2012 {
2013  GeochemicalDatabaseReader database("database/moose_testdb.json");
2014  PertinentGeochemicalSystem model(database, {"H2O", "H+"}, {}, {}, {}, {}, {}, "O2(aq)", "e-");
2015  try
2016  {
2017  model.getIndexOfOriginalBasisSpecies("OH-");
2018  FAIL() << "Missing expected exception.";
2019  }
2020  catch (const std::exception & e)
2021  {
2022  std::string msg(e.what());
2023  ASSERT_TRUE(msg.find("species OH- is not in the original basis") != std::string::npos)
2024  << "Failed with unexpected error message: " << msg;
2025  };
2026  EXPECT_EQ(model.getIndexOfOriginalBasisSpecies("H2O"), (unsigned)0);
2027  EXPECT_EQ(model.getIndexOfOriginalBasisSpecies("H+"), (unsigned)1);
2028 }
2029 
2031 TEST(PertinentGeochemicalSystemTest, originalBasisNames)
2032 {
2033  GeochemicalDatabaseReader database("database/moose_testdb.json");
2034  PertinentGeochemicalSystem model(database, {"H2O", "H+"}, {}, {}, {}, {}, {}, "O2(aq)", "e-");
2035  EXPECT_EQ(model.originalBasisNames()[0], "H2O");
2036  EXPECT_EQ(model.originalBasisNames()[1], "H+");
2037 }
2038 
2040 TEST(PertinentGeochemicalSystemTest, allMinerals)
2041 {
2042  GeochemicalDatabaseReader database("database/moose_testdb_all_minerals.json");
2043  PertinentGeochemicalSystem model(database,
2044  {"H2O", "H+", ">(s)FeOH", ">(w)FeOH", "Fe++", "HCO3-", "O2(aq)"},
2045  {"*"},
2046  {},
2047  {"Fe(OH)3(ppd)fake"},
2048  {"(O-phth)--"},
2049  {">(s)FeO-"},
2050  "O2(aq)",
2051  "e-");
2052  ModelGeochemicalDatabase mgd2 = model.modelGeochemicalDatabase();
2053 
2054  for (const auto & species : mgd2.basis_species_index)
2055  ASSERT_EQ(mgd2.basis_species_mineral[species.second], false);
2056  EXPECT_EQ(mgd2.eqm_species_index.count("Calcite"), (std::size_t)0);
2057  EXPECT_EQ(mgd2.eqm_species_index.count("Calcite_asdf"), (std::size_t)0);
2058  EXPECT_EQ(mgd2.eqm_species_index.count("Fe(OH)3(ppd)fake"), (std::size_t)0);
2059  EXPECT_EQ(mgd2.eqm_species_index.count("Fe(OH)3(ppd)"), (std::size_t)1);
2060  EXPECT_EQ(mgd2.eqm_species_index.count("Goethite"), (std::size_t)1);
2061  EXPECT_EQ(mgd2.eqm_species_index.count("Something"), (std::size_t)1);
2062  for (const auto & species : mgd2.kin_species_index)
2063  if (species.first == "Fe(OH)3(ppd)fake")
2064  ASSERT_EQ(mgd2.kin_species_mineral[species.second], true);
2065  else
2066  ASSERT_EQ(mgd2.kin_species_mineral[species.second], false);
2067 }
ModelGeochemicalDatabase::eqm_species_molecular_weight
std::vector< Real > eqm_species_molecular_weight
all quantities have a molecular weight (g)
Definition: PertinentGeochemicalSystem.h:211
ModelGeochemicalDatabase::basis_species_transported
std::vector< bool > basis_species_transported
basis_species_transported[j] = true iff the j^th basis species is transported in reactive-transport s...
Definition: PertinentGeochemicalSystem.h:170
PertinentGeochemicalSystem::getIndexOfOriginalBasisSpecies
unsigned getIndexOfOriginalBasisSpecies(const std::string &name) const
Definition: PertinentGeochemicalSystem.C:69
ModelGeochemicalDatabase::kin_species_charge
std::vector< Real > kin_species_charge
all kinetic quantities have a charge (mineral charge = 0)
Definition: PertinentGeochemicalSystem.h:312
ModelGeochemicalDatabase::eqm_species_gas
std::vector< bool > eqm_species_gas
eqm_species_gas[i] = true iff the i^th equilibrium species is a gas
Definition: PertinentGeochemicalSystem.h:199
ModelGeochemicalDatabase::basis_species_index
std::unordered_map< std::string, unsigned > basis_species_index
basis_species_index[name] = index of the basis species, within all ModelGeochemicalDatabase internal ...
Definition: PertinentGeochemicalSystem.h:158
mgd
const ModelGeochemicalDatabase mgd
Definition: GeochemistryActivityCoefficientsDebyeHuckelTest.C:28
ModelGeochemicalDatabase::kin_stoichiometry
DenseMatrix< Real > kin_stoichiometry
kin_stoichiometry(i, j) = stoichiometric coefficient for kinetic species "i" in terms of the basis sp...
Definition: PertinentGeochemicalSystem.h:330
ModelGeochemicalDatabase::kin_species_molecular_weight
std::vector< Real > kin_species_molecular_weight
all quantities have a molecular weight (g/mol)
Definition: PertinentGeochemicalSystem.h:315
ModelGeochemicalDatabase::kin_species_molecular_volume
std::vector< Real > kin_species_molecular_volume
all quantities have a molecular volume (cm^3/mol) (only nonzero for minerals, however) ...
Definition: PertinentGeochemicalSystem.h:318
ModelGeochemicalDatabase::eqm_species_index
std::unordered_map< std::string, unsigned > eqm_species_index
eqm_species_index[name] = index of the equilibrium species (secondary aqueous species, redox couples in equilibrium with the aqueous solution, minerals in equilibrium with the aqueous solution, gases in equilibrium with the aqueous solution) within all ModelGeochemicalDatabase internal datastrcutres, with given name
Definition: PertinentGeochemicalSystem.h:190
ModelGeochemicalDatabase::swap_to_original_basis
DenseMatrix< Real > swap_to_original_basis
Swap matrix that allows expression in terms of the original basis.
Definition: PertinentGeochemicalSystem.h:367
database
const GeochemicalDatabaseReader database("database/moose_testdb.json", true, true, false)
ModelGeochemicalDatabase::kin_species_index
std::unordered_map< std::string, unsigned > kin_species_index
kin_species_index[name] = index of the kinetic species, within all ModelGeochemicalDatabase internal ...
Definition: PertinentGeochemicalSystem.h:300
ModelGeochemicalDatabase::eqm_stoichiometry
DenseMatrix< Real > eqm_stoichiometry
eqm_stoichiometry(i, j) = stoichiometric coefficient for equilibrium species "i" in terms of the basi...
Definition: PertinentGeochemicalSystem.h:220
GeochemicalDatabaseReader::getTemperatures
const std::vector< Real > & getTemperatures() const
Get the temperature points that the equilibrium constant is defined at.
Definition: GeochemicalDatabaseReader.C:139
ModelGeochemicalDatabase::basis_species_radius
std::vector< Real > basis_species_radius
all quantities have an ionic radius (Angstrom) for computing activity (mineral radius = 0...
Definition: PertinentGeochemicalSystem.h:176
PertinentGeochemicalSystem::modelGeochemicalDatabase
const ModelGeochemicalDatabase & modelGeochemicalDatabase() const
Return a reference to the ModelGeochemicalDatabase structure.
Definition: PertinentGeochemicalSystem.C:850
eps
const Real eps
Definition: MinimalGeochemicalSystemTest.C:14
PertinentGeochemicalSystem::addKineticRate
void addKineticRate(const KineticRateUserDescription &description)
Adds a rate description for kinetic_species.
Definition: PertinentGeochemicalSystem.C:856
PertinentGeochemicalSystem
Constructs and stores a minimal amount of information that is pertinent to the user-defined geochemic...
Definition: PertinentGeochemicalSystem.h:382
ModelGeochemicalDatabase::kin_log10K
DenseMatrix< Real > kin_log10K
kin_log10K(i, j) = log10(equilibrium constant for the i^th kinetic species at the j^th temperature po...
Definition: PertinentGeochemicalSystem.h:324
model
const PertinentGeochemicalSystem model(database, {"H2O", "H+", "HCO3-", "O2(aq)", "Ca++", ">(s)FeOH", "radius_neg1", "radius_neg1.5"}, {"Calcite"}, {}, {"Calcite_asdf"}, {"CH4(aq)"}, {">(s)FeOCa+"}, "O2(aq)", "e-")
ModelGeochemicalDatabase::eqm_species_radius
std::vector< Real > eqm_species_radius
all quantities have an ionic radius (Angstrom) for computing activity (mineral radius = 0...
Definition: PertinentGeochemicalSystem.h:208
ModelGeochemicalDatabase::eqm_log10K
DenseMatrix< Real > eqm_log10K
eqm_log10K(i, j) = log10(equilibrium constant) for i^th equilibrium species at the j^th temperature p...
Definition: PertinentGeochemicalSystem.h:226
TEST
TEST(PertinentGeochemicalSystemTest, basisExceptions)
Test that:
Definition: MinimalGeochemicalSystemTest.C:23
ModelGeochemicalDatabase::basis_species_molecular_weight
std::vector< Real > basis_species_molecular_weight
all quantities have a molecular weight (g)
Definition: PertinentGeochemicalSystem.h:179
ModelGeochemicalDatabase::original_database
const GeochemicalDatabaseReader * original_database
a pointer to the original database used to build this ModelGeochemicalDatabase
Definition: PertinentGeochemicalSystem.h:149
Real
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
ModelGeochemicalDatabase::eqm_species_transported
std::vector< bool > eqm_species_transported
eqm_species_transported[i] = true iff the i^th eqm species is transported in reactive-transport sims ...
Definition: PertinentGeochemicalSystem.h:202
ModelGeochemicalDatabase::eqm_species_molecular_volume
std::vector< Real > eqm_species_molecular_volume
all quantities have a molecular volume (cm^3) (only nonzero for minerals, however) ...
Definition: PertinentGeochemicalSystem.h:214
PertinentGeochemicalSystem::originalBasisNames
std::vector< std::string > originalBasisNames() const
Definition: PertinentGeochemicalSystem.C:86
ModelGeochemicalDatabase::basis_species_charge
std::vector< Real > basis_species_charge
all quantities have a charge (mineral charge = 0, gas charge = 0, oxide charge = 0) ...
Definition: PertinentGeochemicalSystem.h:173
KineticRateUserDescription
Holds a user-specified description of a kinetic rate.
Definition: GeochemistryKineticRateCalculator.h:75
ModelGeochemicalDatabase::basis_species_gas
std::vector< bool > basis_species_gas
basis_species_gas[j] = true iff the j^th basis species is a gas
Definition: PertinentGeochemicalSystem.h:167
ModelGeochemicalDatabase::basis_species_mineral
std::vector< bool > basis_species_mineral
basis_species_mineral[j] = true iff the j^th basis species is a mineral
Definition: PertinentGeochemicalSystem.h:164
ModelGeochemicalDatabase::eqm_species_charge
std::vector< Real > eqm_species_charge
all quantities have a charge (mineral charge = 0, gas charge = 0, oxide charge = 0) ...
Definition: PertinentGeochemicalSystem.h:205
ModelGeochemicalDatabase::kin_species_name
std::vector< std::string > kin_species_name
kin_species_name[j] = name of the j^th kinetic species
Definition: PertinentGeochemicalSystem.h:303
ModelGeochemicalDatabase
Data structure to hold all relevant information from the database file.
Definition: PertinentGeochemicalSystem.h:99
EM::j
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
ModelGeochemicalDatabase::eqm_species_mineral
std::vector< bool > eqm_species_mineral
eqm_species_mineral[i] = true iff the i^th equilibrium species is a mineral
Definition: PertinentGeochemicalSystem.h:196
ModelGeochemicalDatabase::basis_species_molecular_volume
std::vector< Real > basis_species_molecular_volume
all quantities have a molecular volume (cm^3) (only nonzero for minerals, however) ...
Definition: PertinentGeochemicalSystem.h:182
DenseMatrix< Real >::n
unsigned int n() const
DenseMatrix< Real >::sub_matrix
DenseMatrix sub_matrix(unsigned int row_id, unsigned int row_size, unsigned int col_id, unsigned int col_size) const
ModelGeochemicalDatabase::eqm_species_name
std::vector< std::string > eqm_species_name
eqm_species_name[i] = name of the i^th eqm species
Definition: PertinentGeochemicalSystem.h:193
ModelGeochemicalDatabase::basis_species_name
std::vector< std::string > basis_species_name
basis_species_name[j] = name of the j^th basis species
Definition: PertinentGeochemicalSystem.h:161
GeochemicalDatabaseReader
Class for reading geochemical reactions from a MOOSE geochemical database.
Definition: GeochemicalDatabaseReader.h:213
ModelGeochemicalDatabase::kin_rate
std::vector< KineticRateDefinition > kin_rate
rates given to kinetic species.
Definition: PertinentGeochemicalSystem.h:337
Generated on Wed Aug 20 2025 08:34:06 for https://mooseframework.inl.gov by   doxygen 1.8.14