Influence of inflow's water mass in Reactive transport

[Haohan]

H, I am coupling phreeqcrm with my transport code. To verify my code, I am comparing my result with Phreeqc's transport on Example 11 (Cation exchange).

Here is the phreeqc input with `TRANSPORT` block (using Phreeqc's transport). Its result is used to compare with my code's result.

```
DATABASE C:\Program Files (x86)\USGS\Phreeqc Interactive 3.7.3-15968\database\phreeqc.dat
TITLE Example 11.--Transport and cation exchange.
SOLUTION 0  CaCl2 inflow
        units            mmol/kgw
        temp             25.0
        pH               7.0     charge
        pe               12.5    O2(g)   -0.68
        Ca               0.6
        Cl               1.2
        -water           5
SOLUTION 1-40  Initial solution
        units            mmol/kgw
        temp             25.0
        pH               7.0     charge
        pe               12.5    O2(g)   -0.68
        Na               1.0
        K                0.2
        N(5)             1.2
END
EXCHANGE 1-40
        -equilibrate 1
        X                0.0011
COPY cell 1 101

TRANSPORT
        -cells           40
        -lengths         0.002
        -shifts          100
        -time_step       720.0
        -flow_direction  forward
        -boundary_conditions   flux  flux
        -diffusion_coefficient 0.0
        -dispersivities  0.002
        -correct_disp    true
        -punch_cells     40
        -punch_frequency 1
        -print_cells     40
        -print_frequency 20

(omit the rest ...)

```


I realized that ` -water` parameter in `SOLUTION 0` can affect the output, since different water mass gives different results. However, my code, the water mass doesn't seem to change anything. So I assume that there is something wrong in my code.

So my question is: How does phreeqc's own transport take account of the water mass in the inflow? How does it affect the whole system?


Thanks!

[dlparkhurst]

SOLUTION 0 is advected into cell 1 every 720 seconds. Logically, solution 0 should have the same volume (approximated by the mass of water) as the cells of the column, so that the volume of water in each cell remains the same throughout the transport calculation.

You have defined solution 0 with 5 kgw of water, so the volume of water in each cell will increase from 1 kgw to 5 kgw as the pulse of water moves through the column. After 40 shifts, all cells will have 5 kgw.

Under most cases, you should define SOLUTION 0 with 1 kgw, in which case, the flux of water into the column is 1 L/720 s.

[Haohan]

Thanks you for your reply. I am unsure how to take the water mass into account when coupling phreeqc with my transport code.

My code is similar to "PhreeqcRM BMI Python Example" at https://github.com/usgs-coupled/phreeqcrm/blob/master/swig/python/ex11-advect/ex11-advect.ipynb which couples simple advection transport code and phreeqcrm. I compare its results (component concentration) with the results from the following phreeqc input, to verify the transport code. Both have the same setup.

"""
    TITLE Example 11.--Transport and cation exchange.
    SOLUTION 0  CaCl2
            units            mmol/kgw
            temp             25.0
            pH               7.0     charge
            pe               12.5    O2(g)   -0.68
            Ca               0.6
            Cl               1.2
            -water         1 #kg         
    SOLUTION 1-40  Initial solution for column
            units            mmol/kgw
            temp             25.0
            pH               7.0     charge
            pe               12.5    O2(g)   -0.68
            Na               1.0
            K                0.2
            N(5)             1.2
    EXCHANGE 1-40
            -equilibrate 1
            X                0.0011
    COPY cell 1 101
    TRANSPORT
            -cells           40
            -lengths         0.002
            -shifts          100
            -time_step       720.0
            -flow_direction  forward
            -boundary_conditions   flux  flux
            -diffusion_coefficient 0
            -dispersivities  0.000
            -correct_disp    true
            -punch_cells     40
            -punch_frequency 1
            -print_cells     40
            -print_frequency 20

    SELECTED_OUTPUT
            -file            ex11trn.csv
            -reset           false
            -state           true
            -pH              true
            -pe              true
            -solution        true
            -step            true
            -high_precision  true
            -water           true
            -totals          Na Cl K Ca N(5)
    END
"""

I noticed that when I have "-water 1" in "SOLUTION 0", both gives the same results. But when water mass is not 1 (for example 5kg), phreeqc's transport gives different results, but "PhreeqcRM BMI Python Example" coupling code gives the same result as when water is 1kg. So I think there is something missing in the example coupling code.

How should I change the transport code in this example so that the water mass is taken into account?

[dlparkhurst]

The units for solution are set to mol/L for the Python example with SetUnitsSolution.

The units for exchange are set to mol/L of water. Look at the description of the methods SetUnitsExchange for how and EXCHANGE definition in the initial Phreeqc instance is transferred to the cells of the module.

For the boundary condition (SOLUTION 0), the method InitialPhreeqc2Concentrations is used, which returns mol/L because of the choice of SetUnitsSolution. The concentrations of solution 0 are independent of the mass of water you define. So, you will get the same concentrations whether you define -water 1 or -water 5.

If you want to adjust the flux of water, that is determined by the implementation of boundary conditions in your transport model. PhreeqcRM provides the concentrations, but you will have to account for the volumes of water in your transport code. The simple advection only accounts for concentration. If you want the volume of fluid to change, you probably want to change the saturation according to results from your transport code. (Note that for confined flow, it does not make sense to increase the volume of water in a cell by a factor of 5.) As a hack, I think you could sequentially increase the saturation in cells to 5 as the pulse moves through the column.

[Haohan]

Thanks. Set saturation seems to work.

I am still confused why, assuming SOLUTION 0 and SOLUTION 1-40 have the same amount of water mass, the water mass affects the molalities(concentration) of each element and pH. I tried water mass 0.5 and 1, which gives different results, despite having the same initial concentration input. Shouldn't chemical reaction and ARD transport (https://wwwbrr.cr.usgs.gov/projects/GWC_coupled/phreeqc/html/final-22.html#pgfId-252121) function depends on the concentration instead of total moles of each component?

Here is the PHREEQC input with water mass as variant.

"""
TITLE Example 11.--Transport and cation exchange.
SOLUTION 0  CaCl2
        units            mmol/kgw
        temp             25.0
        pH               7.0     charge
        pe               12.5    O2(g)   -0.68
        Ca               0.6
        Cl               1.2
        -water         ?           
SOLUTION 1-40  Initial solution for column
        units            mmol/kgw
        temp             25.0
        pH               7.0     charge
        pe               12.5    O2(g)   -0.68
        Na               1.0
        K                0.2
        N(5)             1.2
        -water           ?
EXCHANGE 1-40
        -equilibrate 1
        X                0.0011
COPY cell 1 101
TRANSPORT
        -cells           40
        -lengths         0.002
        -shifts          100
        -time_step       720.0
        -flow_direction  forward
        -boundary_conditions  flux flux
        -diffusion_coefficient 0
        -dispersivities  0.002
        -correct_disp    true
        -punch_cells     40
        -punch_frequency 1
        -print_cells     40
        -print_frequency 20

SELECTED_OUTPUT
        -file            ex11trn.csv
        -time            true
        -reset           false
        -state           true
        -pH              true
        -pe              true
        -solution        true
        -step            true
        -high_precision  true
        -water           true
        -totals          Na Cl K Ca N(5)
END
"""

[dlparkhurst]

Reactants in PHREEQC and PhreeqcRM are not concentrations, but simply moles. So, if you keep the number of exchange sites the same, and halve the amount of water, you will get different results.

In PhreeqcRM, the amount of water is determined by the representative volume, the porosity, and the saturation. The initial moles of reactants is determined by the units--SetUnitsExchange, for example--and the number of moles or reactants defined in the InitialPhreeqc instance by RunFile or RunString.