How scale/interpret results from PHREEQC to experimental data and vice-versa

[Angie]

Hi everyone!
I am trying to model metal recovery by ion exchange resins using the results from experimental data. I started with a column of 5.7cm length and an inflow velocity of 0.000095 m/s. The loading phase for the recovery of Cu-Zn runs 100-bed volumes in total and each bed volume has 10cm3.
I defined a grid of 10 cells, 1000 shifts as you can see in the attachment. I obtained a quite similar behavior for the Zn loading curve comparing with the experimental data but not for the Cu curve. I am a little bit confused about how should I scale up or down or interpret the results obtained since PHREEQC simulates 1 liter of water for each cell.
-Does my input file match with my experimental data or I am missing something.
-If it matches, how should I use these results to scale up/down in order to validate them with my experimental data?
-What parameters should I pay attention to for this purpose?

Any enlightenment would be greatly appreciated!

Angie

[MichaelZ20]

Using minteq.v4.dat or core10.dat produces graphs for both Zn and Cu.
I paid attention that the final concentrations with these databases and the default database are almost the same. But with phreeqc.dat (and also with llnl.dat) a graph for Cu is different. What is wrong with the databases?

[dlparkhurst]

I urge you to look at the output file, especially when you don't understand results. If you set the print frequencies the same as the punch frequencies (and don't -reset false), you can see the concentrations of ZnX2 and CuX2. 

As noted in the other post, the database can make a difference in your results. In phreeqc.dat, wateq4f.dat, llnl.dat, and others, an exchanger named X is defined. Your definitions replace the reactions for CuX2 and ZnX2, but other exchange reactions will be as defined in the database (CaX2, MgX2, etc), and these ions will compete for exchange sites in your simulation. If the database you are using does not contain these definitions, then only Cu and Zn will compete for exchange sites.

Note that copper concentration is about 1000 times lower concentration than Zn in SOLUTION 0. If enough pore volumes are run through the column to attain equilibrium with solution 0, then, given the similar exchange constants, Zn will ultimately be much more abundant on the exchanger than Cu. You may need to plot the Cu and Zn on separate scales to see the relative abundances on the exchanger.

As for parameters, your file matches your description of the flow, and you have defined 0.494 equivalents of exchange sites per kg (L) of water. Beyond that, log Ks and the presence or absence of competing exchange reactions will affect the fit with experimental data. EXCHANGE; -exchange_gammas could make a bit of difference.

Code: [Select]

TITLE    Column B - Teruya
SOLUTION 1-10
    temp      20
    pH        2.07 charge
    pe        6.976
    redox     pe
    units     mmol/kgw
    density   1
    Al        0.00111
    Ca        0.00773
    Cu        0.00031
    K         0.00358
    Mg        0.00123
    Na        0.04089
    S         0.005
    Si        0.00285
    -water    1 # kg
EXCHANGE_MASTER_SPECIES
X X-
EXCHANGE_SPECIES                                 #adding different exchange coefficients
X- = X-; -log_k 0
Cu+2 + 2X- = CuX2; -log_k 1.05; -gamma 6.0 0.0
Zn+2 + 2X- = ZnX2; -log_k 0.65; -gamma 5.0 0.0
EXCHANGE 1-10
 X  0.494               #[mol/kgw water], CEC=29 [g/L soil] data sheet TP220
 -equil 1               #calculate iX in equilibrium with SOLUTION 1
END
SOLUTION 0
    temp      20
    pH        3.27
    pe        6.063
    redox     pe
    units     mmol/l
    density   1
    Cu        0.01542
    S         0.651412 charge
    Zn        12.5804
    -water    1 # kg
END
PRINT
 #-reset false
 -status false
TRANSPORT
 -cells  10
 -lengths  0.0057          #0.057 [m] column
 -dispersivities  2e-3     #2mm
 -shifts  1000                 #100 BV
 -time_step 60               #v=6BV/h
 -flow_direction forward
 -boundary_conditions flux flux
 -correct_disp true
 -punch_cells  10
 -punch_frequency 10            #to select shifts for which results will be written to the selected output file
 -print_cells 10
 -print_frequency 10
USER_GRAPH
 -headings PV Cu Zn
 -axis_titles  "Pore Volumes" "Cu, mmol/kgw", "Zn, mmol/kgw"
 -axis_scale  x_axis 0 100 5
 #-axis_scale  y_axis 0 13 1
 -chart_title "Column B TP220"
 -initial_solutions  false
 -plot t
 -start
 10 graph_x  (Step_no + 0.5)/cell_no
 20 graph_y tot("Cu")*1e3
 30 graph_sy tot("Zn")*1e3
 -end
END


[Angie]

Thank you so much for all the information and hints. It was really helpful and it made me understand the simulation deeply :) I am using phreeqc.dat as a database since this one contains all the ions used in my experimental column.

 1. I would like to analyze in more detail the "Exchange Composition" OUTPUT like the "Solution Composition", so how can I get this data also as a second part of my GRID and as the second graph in my CHART. Is it possible to get both simultaneously or I should do it separately?.

2. During the loading in my experimental column I had to change the velocity of the inflow (100-bed-volumes in total): from 1 to 30-bed volume: 6 BV/h; from 31 to 50-bed volume: 1 BV/h; from 51 to 100-bed-volume: 6 BV/h. Then, to simulate the change of velocities in my column I used 3 TRANSPORT and USER_GRAPH keywords, varying only the -shifts and -time_step since the -cell_length is constant, but I cannot obtain the desired graph which should be showing the behavior of the concentrations along with the 100-bed volumes and the 3 different velocities. Instead; a)I obtain only the last third part (TRANSPORT 3-USER_GRAPH 3), Or the 3 parts but initiating from 0 the 3 of them (6 curves in total). Is it totally wrong? or there is another way to do it? 

Hope I could explain properly my inquiries.
Thanks very much again!

Angie

[dlparkhurst]

I think you need only define one USER_GRAPH, but you must separate the TRANSPORT steps with END statements that cause the TRANSPORTs to be run in sequence. (Without ENDs, only the last TRANSPORT definition will be run.)

You can plot both dissolved and exchanged on the same axis.

Code: [Select]

TITLE    Column B - Teruya
SOLUTION 1-10
    temp      20
    pH        2.07 charge
    pe        6.976
    redox     pe
    units     mmol/kgw
    density   1
    Al        0.00111
    Ca        0.00773
    Cu        0.00031
    K         0.00358
    Mg        0.00123
    Na        0.04089
    S         0.005
    Si        0.00285
    -water    1 # kg

EXCHANGE_MASTER_SPECIES
X X-
EXCHANGE_SPECIES                                 #adding different exchange coefficients
X- = X-; -log_k 0
Cu+2 + 2X- = CuX2; -log_k 1.0; -gamma 6.0 0.0
Zn+2 + 2X- = ZnX2; -log_k 0.65;-gamma 5.0 0.0

EXCHANGE 1-10
 X  0.494               #[mol/kgw], CEC=29 [g/L soil] data sheet TP220
 -equil 1               #calculate iX in equilibrium with SOLUTION 1
END

SOLUTION 0
    temp      20
    pH        3.27
    pe        6.063
    redox     pe
    units     mmol/kgw
    density   1
    Cu        0.01542
    S         0.651412 charge
    Zn        12.5804
    -water    1 # kg
END
PRINT
 #-reset false               #to see the output in more detail, exchange composition
 -status false

TRANSPORT 1
 -cells  10
 -lengths  0.0057          #0.057 [m] column
 -dispersivities  2e-3     #2mm
 -shifts  300                  #30 BV
 -time_step 60               #v=6BV/h
 -flow_direction forward
 -boundary_conditions flux flux
 -correct_disp true
 -punch_cells  10
 -punch_frequency 10            #to select shifts for which results will be written to the selected output file
 -print_cells 10
 -print_frequency 10

 USER_GRAPH 1
 -headings PV Cu CuX2 Zn ZnX2
 -axis_titles  "Pore Volumes" "Cu, mmol/kgw" "Zn, mmol/kgw"
 #-axis_scale  x_axis 0 100 5
 #-axis_scale  y_axis 0 13 1         #to have 2 "y" axis
 -chart_title "Column B TP220"
 -initial_solutions  false
 -connect_simulations true
 -plot t
 -start

 10 graph_x  total_time #(Step_no + 0.5)/cell_no
 20 graph_y tot("Cu")*1e3, MOL("CuX2")*1e3
 30 graph_sy tot("Zn")*1e3, MOL("ZnX2")*1e3
 -end
END
TRANSPORT 2
 -shifts  200                  #20 BV
 -time_step 360               #v=1BV/h
END

 TRANSPORT 3
 -shifts  500                  #50 BV
 -time_step 60               #v=6BV/h
END

[Angie]

Thanks! I could plot the graph with the 3 different velocities modifying a little bit the "10 graph_x.." as you can see in the code, but now I am struggling with 2 problems:

1. My selected data is a function of time, and not of bed-volumes anymore. I need to work with bed-volumes for comparison reasons, in case I change velocities or other parameters. Is it possible to plot the same graph as a function of bed-volumes? I couldn't modify it.
2. Since my dissolved and exchanged have a quite different axis-scale I would like to plot 2 different graphs, one for each. What am I missing here? (I attached my code)

Code: [Select]

TITLE    Column B - Teruya
SOLUTION 1-10
    temp      20
    pH        2.07 charge
    pe        6.976
    redox     pe
    units     mmol/kgw
    density   1
    Al        0.00111
    Ca        0.00773
    Cu        0.00031
    K         0.00358
    Mg        0.00123
    Na        0.04089
    S         0.005
    Si        0.00285
    -water    1 # kg

EXCHANGE_MASTER_SPECIES
X X-
EXCHANGE_SPECIES                                 #adding different exchange coefficients
X- = X-; -log_k 0
Cu+2 + 2X- = CuX2; -log_k 1.0; -gamma 6.0 0.0
Zn+2 + 2X- = ZnX2; -log_k 0.65;-gamma 5.0 0.0

EXCHANGE 1-10
 X  0.494               #[mol/kgw], CEC=29 [g/L soil] data sheet TP220
 -equil 1               #calculate iX in equilibrium with SOLUTION 1
END

SOLUTION 0
    temp      20
    pH        3.27
    pe        6.063
    redox     pe
    units     mmol/kgw
    density   1
    Cu        0.01542
    S         0.651412 charge
    Zn        12.5804
    -water    1 # kg
END

PRINT
 #-reset false               #to see the output in more detail, exchange composition
 -status false

TRANSPORT
 -cells  10
 -lengths  0.0057          #0.057 [m] column
 -dispersivities  2e-3     #2mm
 -shifts  300                  #30 BV
 -time_step 60               #v=6BV/h
 -flow_direction forward
 -boundary_conditions flux flux
 -correct_disp true
 -punch_cells  10
 -punch_frequency 10            #to select shifts for which results will be written to the selected output file
 -print_cells 10
 -print_frequency 10

 USER_GRAPH 1
 -headings t Zn Cu
 -axis_titles  "Time, in hours" "Zn, mmol/kgw" "Cu, mmol/kgw"
 -axis_scale  x_axis 0 40 5
 -axis_scale  y_axis 0 13 1
 -axis_scale secondary_y_axis 0 0.017 0.001
 -chart_title "Outflow Solution Composition - Column B TP220"
 -initial_solutions  false
 -connect_simulations false
 -plot_concentration_vs t
 -start
 10 graph_x  TOTAL_TIME/3600 + 240/3600    #(TOTAL_TIME + time step/2) in hours
 20 graph_y tot("Zn")*1e3
 30 graph_sy tot("Cu")*1e3
 -end

 USER_GRAPH 2
 -headings t ZnX2 CuX2
 -axis_scale  x_axis 0 40 5
 -axis_scale  y_axis 0 260 25
 -axis_scale secondary_y_axis 0 0.8 0.05
 -chart_title "Exchanger Composition - Column B TP220"
 -connect_simulations false
-start
 10 graph_x  TOTAL_TIME/3600 + 240/3600    #(TOTAL_TIME + time step/2) in hours
 20 graph_y  MOL("ZnX2")*1e3
 30 graph_sy MOL("CuX2")*1e3
 -end
END

TRANSPORT
 -shifts  200                  #20 BV
 -time_step 360               #v=1BV/h
END

 TRANSPORT
 -shifts  500                  #50 BV
 -time_step 60               #v=6BV/h
END

[dlparkhurst]

(1) USER_GRAPH is called every 10 shifts, so every call to USER_GRAPH 1 represents a pore volume. The PUT statement in USER_GRAPH 1 counts the number of times it has been called, and so represents the number of pore volumes.

(2) You need  -plot_concentration_vs t.

Code: [Select]

TITLE    Column B - Teruya
SOLUTION 1-10
    temp      20
    pH        2.07 charge
    pe        6.976
    redox     pe
    units     mmol/kgw
    density   1
    Al        0.00111
    Ca        0.00773
    Cu        0.00031
    K         0.00358
    Mg        0.00123
    Na        0.04089
    S         0.005
    Si        0.00285
    -water    1 # kg

EXCHANGE_MASTER_SPECIES
X X-
EXCHANGE_SPECIES                                 #adding different exchange coefficients
X- = X-; -log_k 0
Cu+2 + 2X- = CuX2; -log_k 1.0; -gamma 6.0 0.0
Zn+2 + 2X- = ZnX2; -log_k 0.65;-gamma 5.0 0.0

EXCHANGE 1-10
 X  0.494               #[mol/kgw], CEC=29 [g/L soil] data sheet TP220
 -equil 1               #calculate iX in equilibrium with SOLUTION 1
END

SOLUTION 0
    temp      20
    pH        3.27
    pe        6.063
    redox     pe
    units     mmol/kgw
    density   1
    Cu        0.01542
    S         0.651412 charge
    Zn        12.5804
    -water    1 # kg
END

PRINT
 #-reset false               #to see the output in more detail, exchange composition
 -status false

TRANSPORT
 -cells  10
 -lengths  0.0057          #0.057 [m] column
 -dispersivities  2e-3     #2mm
 -shifts  300                  #30 BV
 -time_step 60               #v=6BV/h
 -flow_direction forward
 -boundary_conditions flux flux
 -correct_disp true
 -punch_cells  10
 -punch_frequency 10            #to select shifts for which results will be written to the selected output file
 -print_cells 10
 -print_frequency 10

 USER_GRAPH 1
 -headings t Zn Cu
 -axis_titles  "Pore volumes" "Zn, mmol/kgw" "Cu, mmol/kgw"
 -chart_title "Outflow Solution Composition - Column B TP220"
 -initial_solutions  false
 -connect_simulations false
 -plot_concentration_vs t
 -start
  5 PUT(GET(0) + 1, 0)
 10 graph_x GET(0)
 20 graph_y tot("Zn")*1e3
 30 graph_sy tot("Cu")*1e3
 -end

 USER_GRAPH 2
 -headings t ZnX2 CuX2
 -plot_concentration_vs t
 -axis_titles  "Pore Volumes" "ZnX2, mmol/kgw" "CuX2, mmol/kgw"
 -chart_title "Exchanger Composition - Column B TP220"
 -connect_simulations false
-start
 10 graph_x GET(0)
 20 graph_y  MOL("ZnX2")*1e3
 30 graph_sy MOL("CuX2")*1e3
 -end
END

TRANSPORT
 -shifts  200                  #20 BV
 -time_step 360               #v=1BV/h
END

 TRANSPORT
 -shifts  500                  #50 BV
 -time_step 60               #v=6BV/h
END

[Angie]

(1) It works perfectly! thanks!
(2) I am still obtaining GRAPH 2 (Exchanger composition) only. I am using version 2.18.00, is it possible that this is the problem?

[dlparkhurst]

Get version 3. https://www.usgs.gov/software/phreeqc-version-3