phreeqc实例

长年泉 6.8 .410 .260 .071 .259 .040 .895 .025 .030

NETPATH (Plummer等，1991, 1994)和PHREEQC都具有能够进行反向模拟运算的能力。NETPATH相对PHREEQC而言有两个明显的优点：（1）NETPATH提供了完整的处理同位素的方法，包括同位素的摩尔平衡，同位素分馏，以及C-14示踪，然而PHREEQC仅具有同位素摩尔平衡的能力，和（2）NETPATH提供了数据输入和模拟进展的完全交互环境，然而PHREEQC（版本2）主要是批反应计算的程序。另一方面，在PHREEQC版本1中也是可应用完全的图形用户界面，不过它缺少同位素摩尔平衡能力（Charlton等，1997），PHREEQC版本2也可应用完全的图形用户界面（PHREEQC for Windows, V.E.A. Post, written commun., 1999, http://www.geo.vu.nl/users/posv/phreeqc.html）。PHREEＱＣ相对于NETPATH而言，进行反向模拟的主要优点是，反向模拟计算中包括分析数据不确定度的计算。这种能力使反向模拟具有更强的生命力，也就是说，在输入数据中的微小的变化不会导致模拟计算的摩尔转换产生大的变化。PHREEQC的另一种优点是所设定的任何元素都可以包含在反向模拟计算中，然而，NETPATH则局限于所选取的一组元素，这相对要全面一些。 表47. —由Garrels and Mackenzie (1967)给出的反应物组分和摩尔转换

[摩尔转换单位为毫摩尔/千克水，正数表示溶解，负数表示沉淀]

反应物

“岩盐”

“硬石膏” 高岭石 Ca-蒙脱石 CO2气体 方解石 硅 黑云母 斜长石

组分

NaCl

CaSO4?2H2O Al2Si2O5(OH)4

Ca0.17Al2.33Si3.67O10(OH)2 CO2 CaCO3 SiO2

KMg3AlSi3O10(OH)2 Na0.62Ca0.38Al1.38Si2.62O8

摩尔转移 0.016 -0.015 -0.033 -0.81 0.427 0.115 0.0 0.014 0.175

这两种泉水的分析数据见表46。由Garrels和Mackenzie (1967)给出的假定参与反应的矿物和气体的化学组分以及它们的摩尔转换见表47。反应相的特性和组分的选择是反向模拟最难的部分。一般而言，这部分是根据流动系统和沿流动路径的矿物学的知识得出的；含水层物质成分的微观和化学分析、水和矿物的同位素组分为这部分的反应物的选择提供了附加的信息。没有必要准确地知道是哪一种矿物参加了反应，但却有必要全面知道潜在的反应物。

这个例子的输入数据组见表48。SOLUTION_SPREAD数据块用来定义这两种泉水。

INVERSE_MODELING数据块用来定义所有反向模拟计算的特性，包括所用到的溶液和相，摩尔平衡方程，不确定范围，以及是否所有的或仅是“最小数量” 的模拟将会被打印，是否与不确定限一致的范围内的摩尔转换都将会被计算。一系列的标识符（在连字符之后的子关键字）指定了反向模拟的特性。

表48. —例16的数据输入

TITLE Example 16.--Inverse modeling of Sierra springs SOLUTION_SPREAD -units mmol/L Number pH Si Ca Mg Na K Alkalinity S(6) Cl

1 6.2 0.273 0.078 0.029 0.134 0.028 0.328 0.01 0.014 2 6.8 0.41 0.26 0.071 0.259 0.04 0.895 0.025 0.03

INVERSE_MODELING 1 -solutions 1 2 -uncertainty 0.025 -balances

Ca 0.05 0.025 -phases

Halite Gypsum

Kaolinite precip Ca-montmorillonite precip CO2(g) Calcite

Chalcedony precip Biotite dissolve Plagioclase dissolve -range

PHASES Biotite

KMg3AlSi3O10(OH)2 + 6H+ + 4H2O = K+ + 3Mg+2 + Al(OH)4- + 3H4SiO4 log_k 0.0 # No log_k, Inverse modeling only Plagioclase

Na0.62Ca0.38Al1.38Si2.62O8 + 5.52 H+ + 2.48H2O = \\ 0.62Na+ + 0.38Ca+2 + 1.38Al+3 + 2.62H4SiO4 log_k 0.0 # No log_k, inverse modeling only END

-solutions标识符选择应用的溶液的序号，两种或是更多的溶液序号必须在其标识符之

后列出。如果仅仅是给出了两种溶液的序号，第二种溶液假定是从第一种溶液演化来的。如果给出了更多的溶液序号，最后列出的溶液假定是由前面溶液的混合物演化来的。在反向模拟中所用到溶液的定义与在PHREEQC中模拟所用到任何溶液的定义是一样的。通常分析数据是在SOLUTION或SOLUTION_SPREAD数据块中输入的，但在当前或是前面模拟中批反应计算中所定义的溶液也可使用，条件是它们须用SAVE关键字进行保存。

-uncertainty标识符设定了每一分析数据的缺省不确定限。在这个例子中，相对不确定

限为0.025(2.5%)，它适用于除pH之外的所有分析数据。在缺省状态下，pH值的不确定限为0.05个单位。pH的不确限也可以用-balances标识设置为绝对值（标准单位）。任何溶液任何数据的不确定限，能够用-balances标识符明确地设置为小数值或是绝对值（摩尔数；碱度当量）。

在缺省状态下，每个反向模拟包括包含于-phases中（除了氢和氧）任何相的每一种元

素的摩尔平衡方程。如果需要不包括在相中的元素的摩尔平衡方程，也就是说，元素没有源或汇（守恒混合），-balances标识符能够用来包含那些在反向模拟方程中的公式中的元素（见例17）。另外，-balances标识符能够用来指定每一种溶液中每一种元素的不确定限。在这个例子中为了证明这个目的，钙离子的不确定限在溶液1中设置为0.05（5%）和在溶液2中设置为0.025（2.5%）

在反向模拟计算中用到的相用-phases标识符进行定义。另外，这个标识符可以用来约

束某相仅被溶解，或是仅为沉淀。在这个例子中，高岭石，钙-蒙脱石和玉髓(SiO2)都要求仅为沉淀。这意味着在包含有高岭石相的任何模拟中高岭石都将作为沉淀生成（摩尔转移为负）；同样，也适用于钙-蒙脱石和玉髓。如果在反向模拟中出现黑云母和斜长石，则要求它们仅为溶解（摩尔转移为正）。

在反向模拟中用到所有的相必须在数据块文件中或是在输入文件中用PHASES或

EXCHANGE_SPECIES数据块来定义。这样在所有默认数据库文件中所定义的所有相，phreeqc.dat或wateq4f.dat在反向模拟中都可用。黑云母和斜长石不在缺省的数据库文件phreeqc.dat中，因此，在输入文件设置中的以PHASES数据块显性定义。为简单起见，这些相的log K’s设置为0，由于仅是应用到了矿物的化学计量，这样设置不会影响反向模拟。然而，这些相饱和指数的计算将会是不合实际的。在反向模拟中用到的相均须具有电荷平衡反应。这是因为每一种溶液都必须处于电荷平衡状态。在用最佳化方法减少目标函数的过程中（见“反向模拟方程和数值方法”），通过在不确定限范围内调整元素的浓度，保证每个模型中的每一种溶液均达到电荷平衡。如果一种溶液在给定的不确定限范围内不能调整到电荷平衡，那么在输出文件中将会标记出该溶液和找不到模型。由于在模拟过程中的所有溶液都处于电荷平衡，相也必须达到电荷平衡，或在模型中没有该相。注意，斜长石的反应（表48）有两行，由于在第一行的最后以反斜杠“\\”结束，程序将这两行解释为单逻辑行。

-range标识符表明，除了发现所有的反向模拟，在不确定限的范围之内，发现的每一种

模型都需要进行额外的计算以决定摩尔转移值的范围。

下列方程也包含于每一种反向模拟之中：系统中每种元素的摩尔数以及每种元素化合价

（如在-phases相中所定义的元素和在-balances中所列出的元素）的摩尔平衡方程，每一种溶液的电荷平衡，体系的碱度平衡，体系的电子平衡和体系的水平衡。在这些方程中的未知数包括相转移的摩尔数，氧化还原反应转移的摩尔数，和每一种溶液中每一种元素的不确定度（不包括氢和氧）。不确定性包括每一种溶液的碱度和pH值。优化方法解决了满足所有方程的未知数，并满足所有的不确定限，同时最小化了目标函数，这是不确定未知数的加权

和（见“反向模拟方程和数值方法”）。

该例中发现的两个反向模拟的结果见表49。模型中每一种溶液的每种成分用三列描述。

所有的列除了pH（如果有同位素，也包括在内），单位都是mol/kgw。第一栏包括溶液的初始分析数据（输入）。第二栏包括模拟中所计算的分析数据的修正(Delta)。这些修正必须在指定的不确定范围内。第三栏包含溶液修订后的分析数据，它等于原始数据加上修正值(Input+Delta)。

在溶液列表之后，将会打印出反向模拟中每一种溶液的相对比例（溶液分数）。在模型

中仅有两种溶液，通常情况下每种溶液的分数都为1.0。在反向模拟中如果超过两种溶液，除了最后一种溶液外，其它溶液所占的百分数的和通常等于1.0。这个分数实际上是来自于水的摩尔平衡，因此，如果消耗了含水矿物或是产生了一定数量的水，或是模拟蒸发（见例17），分数的和则可能不是1.0。在这个例子中，所有的分数都等于1.0；来自于石膏溶解水的量很小，不影响混合分数的四个有效数字。这个数据块的第二栏和第三栏给出的溶液分数的最小值和最大值，这是在指定的不确定限范围内取得的。只有使用了-range标识符，这两栏才为非0。

在列表中的下一个数据块包括三列，分别描述了相的摩尔转移（相摩尔转移）。第一栏

包含与打印到溶液列表中修正的浓度相一致的反向模拟。在这个例子中，修正的溶液1加上在第一栏中的摩尔转移恰好等于修正后的溶液2。摩尔转移为正表示溶解；为负表示沉淀。（注意：在多组反应计算中相集合的摩尔转移是相对于相，而不是相对于溶液的）第二和第三列的摩尔转移是在指定不确定限的范围内每种相摩尔转移的最小值和最大值。只有使用了-range标识符，这两列为非0。一般而言，这个最小值和最大值不是独立的，也就是说，取得一种相最大的摩尔转移将会影响到摩尔平衡模拟中其它相的摩尔转移。

在这个反向模拟中没有计算氧化还原摩尔转移。如果计算了某个氧化还原摩尔转移，每

种元素不同化合价之间的摩尔转移将会被打印在标题“氧化还原摩尔转移”之下。

下一个数据块打印了在某种程度上与这个模拟修正的分析数据有关的结果；如果没有进

行任何修正，打印的这三个数将都为0。首先，打印偏差的和，这是反向模拟不确定范围的总和（(??uqm?m,qm,q??q?m?q?m,qum,q)），其次，是反向模拟不确定范围加权并修正的每

种元素浓度和同位素组分的总和（

??uqm?m,qm,q，delta/uncertainty不确定范围的和）。最后，

将会打印每种溶液中任何元素或是同位素组分（元素浓度的最大部分）的最大分数修正值（元素浓度的最大修正误差）。所有的这三个值应用于打印在左栏中的模型。

对于给定的摩尔平衡模型，如果没有发现包含适当的溶液和相子集的更为简单的模型，将会为给定的模型打印“模拟包含最少的相”。

所有模拟都打印之后是简短的计算总结，列出了发现的模型数量、最小模型的数量（包

含最少相的模拟），所验证的不可实行模拟的数目，不等式计算调用的数次cl1（计算的时间