mirror of
https://git.gfz-potsdam.de/naaice/iphreeqc.git
synced 2025-12-16 16:44:49 +01:00
Updating to Tony's latest
git-svn-id: svn://136.177.114.72/svn_GW/phreeqc3/trunk@7069 1feff8c3-07ed-0310-ac33-dd36852eb9cd
This commit is contained in:
David L Parkhurst
parent
2671345a70
commit
244dcf43bf
176
ex15a
Normal file
176
ex15a
Normal file
@ -0,0 +1,176 @@
|
||||
DATABASE ex15.dat
|
||||
TITLE Example 15.--1D Transport: Kinetic Biodegradation, Cell Growth, and Sorption
|
||||
***********
|
||||
PLEASE NOTE: This problem requires database file ex15.dat!!
|
||||
***********
|
||||
PRINT
|
||||
-reset false
|
||||
-echo_input true
|
||||
-status false
|
||||
SOLUTION 0 Pulse solution with NTA and cobalt
|
||||
units umol/L
|
||||
pH 6
|
||||
C .49
|
||||
O(0) 62.5
|
||||
Nta 5.23
|
||||
Co 5.23
|
||||
Na 1000
|
||||
Cl 1000
|
||||
SOLUTION 1-10 Background solution initially filling column
|
||||
units umol/L
|
||||
pH 6
|
||||
C .49
|
||||
O(0) 62.5
|
||||
Na 1000
|
||||
Cl 1000
|
||||
COPY solution 0 100 # for use later on, and in
|
||||
COPY solution 1 101 # 20 cells model
|
||||
END
|
||||
RATES Rate expressions for the four kinetic reactions
|
||||
#
|
||||
HNTA-2
|
||||
-start
|
||||
10 Ks = 7.64e-7
|
||||
20 Ka = 6.25e-6
|
||||
30 qm = 1.407e-3/3600
|
||||
40 f1 = MOL("HNta-2")/(Ks + MOL("HNta-2"))
|
||||
50 f2 = MOL("O2")/(Ka + MOL("O2"))
|
||||
60 rate = -qm * KIN("Biomass") * f1 * f2
|
||||
70 moles = rate * TIME
|
||||
80 PUT(rate, 1) # save the rate for use in Biomass rate calculation
|
||||
90 SAVE moles
|
||||
-end
|
||||
#
|
||||
Biomass
|
||||
-start
|
||||
10 Y = 65.14
|
||||
20 b = 0.00208/3600
|
||||
30 rate = GET(1) # uses rate calculated in HTNA-2 rate calculation
|
||||
40 rate = -Y*rate -b*M
|
||||
50 moles = -rate * TIME
|
||||
60 if (M + moles) < 0 then moles = -M
|
||||
70 SAVE moles
|
||||
-end
|
||||
#
|
||||
Co_sorption
|
||||
-start
|
||||
10 km = 1/3600
|
||||
20 kd = 5.07e-3
|
||||
30 solids = 3.75e3
|
||||
40 rate = -km*(MOL("Co+2") - (M/solids)/kd)
|
||||
50 moles = rate * TIME
|
||||
60 if (M - moles) < 0 then moles = M
|
||||
70 SAVE moles
|
||||
-end
|
||||
#
|
||||
CoNta_sorption
|
||||
-start
|
||||
10 km = 1/3600
|
||||
20 kd = 5.33e-4
|
||||
30 solids = 3.75e3
|
||||
40 rate = -km*(MOL("CoNta-") - (M/solids)/kd)
|
||||
50 moles = rate * TIME
|
||||
60 if (M - moles) < 0 then moles = M
|
||||
70 SAVE moles
|
||||
-end
|
||||
KINETICS 1-10 Four kinetic reactions for all cells
|
||||
HNTA-2
|
||||
-formula C -3.12 H -1.968 O -4.848 N -0.424 Nta 1.
|
||||
Biomass
|
||||
-formula H 0.0
|
||||
-m 1.36e-4
|
||||
Co_sorption
|
||||
-formula CoCl2
|
||||
-m 0.0
|
||||
-tol 1e-11
|
||||
CoNta_sorption
|
||||
-formula NaCoNta
|
||||
-m 0.0
|
||||
-tol 1e-11
|
||||
COPY kinetics 1 101 # to use with 20 cells
|
||||
END
|
||||
SELECTED_OUTPUT
|
||||
-file ex15.sel
|
||||
-mol Nta-3 CoNta- HNta-2 Co+2
|
||||
USER_PUNCH
|
||||
-headings hours Co_sorb CoNta_sorb Biomass
|
||||
-start
|
||||
10 punch TOTAL_TIME/3600 + 3600/2/3600
|
||||
20 punch KIN("Co_sorption")/3.75e3
|
||||
30 punch KIN("CoNta_sorption")/3.75e3
|
||||
40 punch KIN("Biomass")
|
||||
-end
|
||||
USER_GRAPH
|
||||
-headings 10_cells: Co+2 CoNTA- HNTA-2 pH
|
||||
-chart_title "Example 15"
|
||||
-axis_titles "Time / hours" "umol / kgw" "pH"
|
||||
-axis_scale x_axis 0 75
|
||||
-axis_scale y_axis 0 4
|
||||
-axis_scale secondary_y_axis 5.799 6.8 0.2 0.1
|
||||
-plot_concentration_vs t
|
||||
-start
|
||||
10 x = TOTAL_TIME/3600 + 3600/2/3600
|
||||
20 PLOT_XY -1, -1, line_width = 0, symbol_size = 0
|
||||
30 PLOT_XY x, MOL("Co+2") * 1e6, color = Red, line_width = 0, symbol_size = 4
|
||||
40 PLOT_XY x, MOL("CoNta-") * 1e6, color = Green, line_width = 0, symbol_size = 4
|
||||
50 PLOT_XY x, MOL("HNta-2") * 1e6, color = Blue, line_width = 0, symbol_size = 4
|
||||
60 PLOT_XY x, -LA("H+"), y-axis = 2, color = Magenta, line_width = 0, symbol_size = 4
|
||||
-end
|
||||
TRANSPORT First 20 hours have NTA and cobalt in infilling solution
|
||||
-cells 10
|
||||
-lengths 1
|
||||
-shifts 20
|
||||
-time_step 3600
|
||||
-flow_direction forward
|
||||
-boundary_conditions flux flux
|
||||
-dispersivities .05
|
||||
-correct_disp true
|
||||
-diffusion_coefficient 0.0e-9
|
||||
-punch_cells 10
|
||||
-punch_frequency 1
|
||||
-print_cells 10
|
||||
-print_frequency 5
|
||||
COPY solution 101 0 # initial column solution becomes influent
|
||||
END
|
||||
TRANSPORT Last 55 hours with background infilling solution
|
||||
-shifts 55
|
||||
COPY cell 100 0 # for the 20 cell model...
|
||||
COPY cell 101 1-20
|
||||
END
|
||||
USER_PUNCH
|
||||
-start
|
||||
10 punch TOTAL_TIME/3600 + 3600/4/3600
|
||||
20 punch KIN("Co_sorption")/3.75e3
|
||||
30 punch KIN("CoNta_sorption")/3.75e3
|
||||
40 punch KIN("Biomass")
|
||||
-end
|
||||
USER_GRAPH
|
||||
-headings 20_cells: Co+2 CoNTA- HNTA-2 pH
|
||||
-start
|
||||
10 x = TOTAL_TIME/3600 + 3600/4/3600
|
||||
20 PLOT_XY -1, -1, line_width = 0, symbol_size = 0
|
||||
30 PLOT_XY x, MOL("Co+2") * 1e6, color = Red, symbol_size = 0
|
||||
40 PLOT_XY x, MOL("CoNta-") * 1e6, color = Green, symbol_size = 0
|
||||
50 PLOT_XY x, MOL("HNta-2") * 1e6, color = Blue, symbol_size = 0
|
||||
60 PLOT_XY x, -LA("H+"), y-axis = 2, color = Magenta, symbol_size = 0
|
||||
-end
|
||||
TRANSPORT First 20 hours have NTA and cobalt in infilling solution
|
||||
-cells 20
|
||||
-lengths 0.5
|
||||
-shifts 40
|
||||
-initial_time 0
|
||||
-time_step 1800
|
||||
-flow_direction forward
|
||||
-boundary_conditions flux flux
|
||||
-dispersivities .05
|
||||
-correct_disp true
|
||||
-diffusion_coefficient 0.0e-9
|
||||
-punch_cells 20
|
||||
-punch_frequency 2
|
||||
-print_cells 20
|
||||
-print_frequency 10
|
||||
COPY cell 101 0
|
||||
END
|
||||
TRANSPORT Last 55 hours with background infilling solution
|
||||
-shifts 110
|
||||
END
|
||||
176
ex15b
Normal file
176
ex15b
Normal file
@ -0,0 +1,176 @@
|
||||
DATABASE ex15.dat
|
||||
TITLE Example 15.--1D Transport: Kinetic Biodegradation, Cell Growth, and Sorption
|
||||
***********
|
||||
PLEASE NOTE: This problem requires database file ex15.dat!!
|
||||
***********
|
||||
PRINT
|
||||
-reset false
|
||||
-echo_input true
|
||||
-status false
|
||||
SOLUTION 0 Pulse solution with NTA and cobalt
|
||||
units umol/L
|
||||
pH 6
|
||||
C .49
|
||||
O(0) 62.5
|
||||
Nta 5.23
|
||||
Co 5.23
|
||||
Na 1000
|
||||
Cl 1000
|
||||
SOLUTION 1-10 Background solution initially filling column
|
||||
units umol/L
|
||||
pH 6
|
||||
C .49
|
||||
O(0) 62.5
|
||||
Na 1000
|
||||
Cl 1000
|
||||
COPY solution 0 100 # to use with 20 cells
|
||||
COPY solution 1 101
|
||||
END
|
||||
RATES Rate expressions for the four kinetic reactions
|
||||
#
|
||||
HNTA-2
|
||||
-start
|
||||
10 Ks = 7.64e-7
|
||||
20 Ka = 6.25e-6
|
||||
30 qm = 1.407e-3/3600
|
||||
40 f1 = MOL("HNta-2")/(Ks + MOL("HNta-2"))
|
||||
50 f2 = MOL("O2")/(Ka + MOL("O2"))
|
||||
60 rate = -qm * KIN("Biomass") * f1 * f2
|
||||
70 moles = rate * TIME
|
||||
80 PUT(rate, 1) # save the rate for use in Biomass rate calculation
|
||||
90 SAVE moles
|
||||
-end
|
||||
#
|
||||
Biomass
|
||||
-start
|
||||
10 Y = 65.14
|
||||
20 b = 0.00208/3600
|
||||
30 rate = GET(1) # uses rate calculated in HTNA-2 rate calculation
|
||||
40 rate = -Y*rate -b*M
|
||||
50 moles = -rate * TIME
|
||||
60 if (M + moles) < 0 then moles = -M
|
||||
70 SAVE moles
|
||||
-end
|
||||
#
|
||||
Co_sorption
|
||||
-start
|
||||
10 km = 1/3600
|
||||
20 kd = 5.07e-3
|
||||
30 solids = 3.75e3
|
||||
40 rate = -km*(MOL("Co+2") - (M/solids)/kd)
|
||||
50 moles = rate * TIME
|
||||
60 if (M - moles) < 0 then moles = M
|
||||
70 SAVE moles
|
||||
-end
|
||||
#
|
||||
CoNta_sorption
|
||||
-start
|
||||
10 km = 1/3600
|
||||
20 kd = 5.33e-4
|
||||
30 solids = 3.75e3
|
||||
40 rate = -km*(MOL("CoNta-") - (M/solids)/kd)
|
||||
50 moles = rate * TIME
|
||||
60 if (M - moles) < 0 then moles = M
|
||||
70 SAVE moles
|
||||
-end
|
||||
KINETICS 1-10 Four kinetic reactions for all cells
|
||||
HNTA-2
|
||||
-formula C -3.12 H -1.968 O -4.848 N -0.424 Nta 1.
|
||||
Biomass
|
||||
-formula H 0.0
|
||||
-m 1.36e-4
|
||||
Co_sorption
|
||||
-formula CoCl2
|
||||
-m 0.0
|
||||
-tol 1e-11
|
||||
CoNta_sorption
|
||||
-formula NaCoNta
|
||||
-m 0.0
|
||||
-tol 1e-11
|
||||
# -cvode true; -cvode_order 3 # uncomment with 1000 times higher sorption rates
|
||||
COPY kinetics 1 101 # to use with 20 cells
|
||||
END
|
||||
SELECTED_OUTPUT
|
||||
-file ex15.sel
|
||||
-mol Nta-3 CoNta- HNta-2 Co+2
|
||||
USER_PUNCH
|
||||
-headings hours Co_sorb CoNta_sorb Biomass
|
||||
-start
|
||||
10 punch TOTAL_TIME/3600 + 3600/2/3600
|
||||
20 punch KIN("Co_sorption")/3.75e3
|
||||
30 punch KIN("CoNta_sorption")/3.75e3
|
||||
40 punch KIN("Biomass")
|
||||
-end
|
||||
TRANSPORT First 20 hours have NTA and cobalt in infilling solution
|
||||
-cells 10
|
||||
-lengths 1
|
||||
-shifts 20
|
||||
-time_step 3600
|
||||
-flow_direction forward
|
||||
-boundary_conditions flux flux
|
||||
-dispersivities .05
|
||||
-correct_disp true
|
||||
-diffusion_coefficient 0.0e-9
|
||||
-punch_cells 10
|
||||
-punch_frequency 1
|
||||
-print_cells 10
|
||||
-print_frequency 5
|
||||
-warnings false
|
||||
USER_GRAPH
|
||||
-headings 10_cells: Co+2 CoNTA- Biomass
|
||||
-chart_title "Example 15, Sorbed Species"
|
||||
-axis_titles "Time / hours" "nmol / kgw" "Biomass / (mg/L)"
|
||||
-axis_scale x_axis 0 75
|
||||
-axis_scale y_axis 0 2
|
||||
-axis_scale secondary_y_axis 0 0.4
|
||||
-plot_concentration_vs t
|
||||
-start
|
||||
10 x = TOTAL_TIME/3600 + 3600/2/3600
|
||||
20 PLOT_XY -1, -1, line_width = 0, symbol_size = 0
|
||||
30 PLOT_XY x, KIN("Co_sorption") / 3.75e3 * 1e9, color = Red, line_width = 0, symbol_size = 4
|
||||
40 PLOT_XY x, KIN("CoNta_sorption") / 3.75e3 * 1e9, color = Green, line_width = 0, symbol_size = 4
|
||||
50 PLOT_XY x, KIN("Biomass") * 1e3, y-axis = 2, color = Magenta, line_width = 0, symbol_size = 4
|
||||
-end
|
||||
COPY solution 101 0
|
||||
END
|
||||
TRANSPORT Last 55 hours with background infilling solution
|
||||
-shifts 55
|
||||
COPY cell 100 0
|
||||
COPY cell 101 1-20
|
||||
END
|
||||
USER_PUNCH
|
||||
-start
|
||||
10 punch TOTAL_TIME/3600 + 3600/4/3600
|
||||
20 punch KIN("Co_sorption")/3.75e3
|
||||
30 punch KIN("CoNta_sorption")/3.75e3
|
||||
40 punch KIN("Biomass")
|
||||
-end
|
||||
TRANSPORT First 20 hours have NTA and cobalt in infilling solution
|
||||
-cells 20
|
||||
-lengths 0.5
|
||||
-shifts 40
|
||||
-initial_time 0
|
||||
-time_step 1800
|
||||
-flow_direction forward
|
||||
-boundary_conditions flux flux
|
||||
-dispersivities .05
|
||||
-correct_disp true
|
||||
-diffusion_coefficient 0.0e-9
|
||||
-punch_cells 20
|
||||
-punch_frequency 2
|
||||
-print_cells 20
|
||||
-print_frequency 10
|
||||
USER_GRAPH
|
||||
-headings 20_cells: Co+2 CoNTA- Biomass
|
||||
-start
|
||||
10 x = TOTAL_TIME/3600 + 3600/4/3600
|
||||
20 PLOT_XY -1, -1, line_width = 0, symbol_size = 0
|
||||
30 PLOT_XY x, KIN("Co_sorption") / 3.75e3 * 1e9, color = Red, symbol_size = 0
|
||||
40 PLOT_XY x, KIN("CoNta_sorption") / 3.75e3 * 1e9, color = Green, symbol_size = 0
|
||||
60 PLOT_XY x, KIN("Biomass") * 1e3, y-axis = 2, color = Magenta, symbol_size = 0
|
||||
-end
|
||||
COPY cell 101 0
|
||||
END
|
||||
TRANSPORT Last 55 hours with background infilling solution
|
||||
-shifts 110
|
||||
END
|
||||
349
ex21
Normal file
349
ex21
Normal file
@ -0,0 +1,349 @@
|
||||
TITLE Diffusion through Opalinus Clay in a radial diffusion cell, Appelo, Van Loon and Wersin,
|
||||
2010, GCA 74, 1201
|
||||
SOLUTION_MASTER_SPECIES
|
||||
# element species alk gfw_formula element_gfw
|
||||
Hto Hto 0.0 20 20
|
||||
Na_tr Na_tr+ 0.0 22 22
|
||||
Cl_tr Cl_tr- 0.0 36 36
|
||||
Cs Cs+ 0.0 132.905 132.905
|
||||
SOLUTION_SPECIES
|
||||
Hto = Hto; log_k 0; -gamma 1e6 0; -dw 2.236e-9
|
||||
Na_tr+ = Na_tr+; log_k 0; -gamma 4.0 0.075; -dw 1.33e-9; -erm_ddl 1.23
|
||||
Cl_tr- = Cl_tr-; log_k 0; -gamma 3.5 0.015; -dw 1.31e-9 # dw = dw(water) / 1.55 = 2.03e-9 / 1.55
|
||||
Cs+ = Cs+; log_k 0; -gamma 3.5 0.015; -dw 2.07e-9; -erm_ddl 1.23
|
||||
SURFACE_MASTER_SPECIES
|
||||
Su_fes Su_fes- # Frayed Edge Sites
|
||||
Su_ii Su_ii- # Type II sites of intermediate strength
|
||||
Su_ Su_- # Double layer, planar sites are modeled with EXCHANGE
|
||||
SURFACE_SPECIES
|
||||
Su_fes- = Su_fes-; log_k 0
|
||||
Na+ + Su_fes- = NaSu_fes; log_k 10
|
||||
Na_tr+ + Su_fes- = Na_trSu_fes; log_k 10
|
||||
K+ + Su_fes- = KSu_fes; log_k 12.4
|
||||
Cs+ + Su_fes- = CsSu_fes; log_k 17.14
|
||||
|
||||
Su_ii- = Su_ii-; log_k 0
|
||||
Na+ + Su_ii- = NaSu_ii; log_k 10
|
||||
Na_tr+ + Su_ii- = Na_trSu_ii; log_k 10
|
||||
K+ + Su_ii- = KSu_ii; log_k 12.1
|
||||
Cs+ + Su_ii- = CsSu_ii; log_k 14.6
|
||||
|
||||
Su_- = Su_-; log_k 0
|
||||
|
||||
EXCHANGE_SPECIES
|
||||
Na_tr+ + X- = Na_trX; log_k 0.0; -gamma 4.0 0.075
|
||||
Cs+ + X- = CsX; log_k 2.04; -gamma 3.5 0.015
|
||||
|
||||
SOLUTION 0-2 column with only cell 1, two boundary solutions 0 and 2.
|
||||
Na 1; Cl 1
|
||||
END
|
||||
|
||||
KNOBS; -diagonal_scale true #; -tolerance 1e-17 # if necessary, try lower tolereance
|
||||
|
||||
SOLUTION 3 tracer solution
|
||||
pH 7.6; pe 14 O2(g) -1.0; temp 23
|
||||
Na 240; K 1.61; Mg 16.9; Ca 25.8; Sr 0.505
|
||||
Cl 300; S(6) 14.1; Fe(2) 0.0; Alkalinity 0.476
|
||||
# uncomment tracer concentrations and kg water 1 by 1...
|
||||
Hto 1.14e-6; -water 0.2
|
||||
# Cl_tr 2.505e-2; -water 0.502
|
||||
# Cs 1; Na_tr 1.87e-7; -water 1.02
|
||||
SELECTED_OUTPUT
|
||||
-file radial; -reset false
|
||||
USER_PUNCH
|
||||
# Define symbols and pi...
|
||||
1 nl$ = EOL$ # newline
|
||||
2 x$ = CHR$(35) # cross '#'
|
||||
3 sc$ = CHR$(59) # semicolon ';'
|
||||
4 pi = 2 * ARCTAN(1e10) # 3.14159...
|
||||
|
||||
# Define experimental parameters...
|
||||
10 height = 0.052 # length of the clay cylinder / m
|
||||
20 r_int = 6.58e-3 # inner radius of clay cylinder / m
|
||||
30 r_ext = 25.4e-3 # outer radius
|
||||
40 thickn_filter1 = 1.8e-3 # tracer-in filter thickness / m
|
||||
50 thickn_filter2 = 1.6e-3 # tracer-out filter thickness / m
|
||||
60 por_filter1 = 0.418 # porosity
|
||||
70 por_filter2 = 0.367
|
||||
80 G_filter1 = 4.18 # geometrical factor. (for filters, G = por / 10)
|
||||
90 G_filter2 = 3.67
|
||||
100 V_end = 0.2 # volume of the tracer-out solution / L
|
||||
110 thickn_clay = r_ext - r_int # clay thickness / m
|
||||
120 por_clay = 0.159
|
||||
130 rho_b_eps = 2.7 * (1 - por_clay) / por_clay # clay bulk density / porosity / (kg/L)
|
||||
140 CEC = 0.12 * rho_b_eps # CEC / (eq/L porewater)
|
||||
150 A_por = 37e3 * rho_b_eps # pore surface area / (m2/L porewater)
|
||||
|
||||
160 DIM tracer$(4), exp_time(4), scale_y1$(4), scale_y2$(4), profile_y1$(4), profile_y2$(4)
|
||||
170 DATA 'Hto', 'Cl_tr', 'Na_tr', 'Cs'
|
||||
180 READ tracer$(1), tracer$(2), tracer$(3), tracer$(4)
|
||||
# experimental times (seconds) for HTO, 36Cl, 22Na and Cs, respectively,
|
||||
# in order of increasing times...
|
||||
200 DATA 86400 * 20, 86400 * 40, 86400 * 45, 86400 * 1000
|
||||
210 READ exp_time(1), exp_time(2), exp_time(3), exp_time(4)
|
||||
# scale y1-axis (flux) (not used)...
|
||||
230 DATA '1', '1', '1', '1'
|
||||
240 READ scale_y1$(1), scale_y1$(2), scale_y1$(3), scale_y1$(4)
|
||||
# scale y2-axis (mass) (not used)...
|
||||
260 DATA '1', '1', '1', '1'
|
||||
270 READ scale_y2$(1), scale_y2$(2), scale_y2$(3), scale_y2$(4)
|
||||
# scale max of the profile y axes...
|
||||
280 DATA '0 1.2e-9', '0 2.5e-5', '0 2e-10', '0 auto'
|
||||
290 READ profile_y1$(1), profile_y1$(2), profile_y1$(3), profile_y1$(4)
|
||||
300 DATA '0 1.2e-9', '0 2.5e-5', '0 6e-10', '0 auto'
|
||||
310 READ profile_y2$(1), profile_y2$(2), profile_y2$(3), profile_y2$(4)
|
||||
|
||||
# Define model parameters...
|
||||
350 Dw = 2.5e-9 # default tracer diffusion coefficient / (m2/s)
|
||||
360 nfilt1 = 1 # number of cells in filter 1
|
||||
370 nfilt2 = 1 # number of cells in filter 2
|
||||
380 nclay = 11 # number of clay cells
|
||||
390 f_free = 0.117 # fraction of free pore water (0.01 - 1)
|
||||
400 f_DL_charge = 0.45 # fraction of CEC charge in electrical double layer
|
||||
410 tort_n = -0.99 # exponent in Archie's law, -1.045 without filters
|
||||
420 G_clay = por_clay^tort_n # geometrical factor
|
||||
430 interlayer_D$ = 'false' # 'true' or 'false' for interlayer diffusion
|
||||
440 G_IL = 700 # geometrical factor for clay interlayers
|
||||
450 punch_time = 60 * 60 * 6 # punch time / seconds
|
||||
460 profile$ = 'true' # 'true' or 'false' for c/x profile visualization
|
||||
470 IF nfilt1 = 0 THEN thickn_filter1 = 0
|
||||
480 IF nfilt2 = 0 THEN thickn_filter2 = 0
|
||||
|
||||
# See which tracer is present...
|
||||
490 IF tot("Hto") > 1e-10 THEN tracer = 1 ELSE \
|
||||
IF tot("Cl_tr") > 1e-10 THEN tracer = 2 ELSE tracer = 3
|
||||
|
||||
# Define clay pore water composition...
|
||||
520 sol$ = nl$ + ' pH 7.6' + sc$ +' pe 14 O2(g) -1.0' + sc$ +' temp 23'
|
||||
530 sol$ = sol$ + nl$ + ' Na 240' + sc$ +' K 1.61' + sc$ +' Mg 16.9' + sc$ +' Ca 25.8' + sc$ +' Sr 0.505'
|
||||
540 sol$ = sol$ + nl$ + ' Cl 300' + sc$ +' S(6) 14.1' + sc$ +' Fe(2) 0.0' + sc$ +' Alkalinity 0.476'
|
||||
|
||||
# Define phases in which the tracers precipitate...
|
||||
550 tracer_phases$ = nl$ + 'PHASES '
|
||||
560 tracer_phases$ = tracer_phases$ + nl$ + ' A_Hto' + nl$ + ' Hto = Hto' + sc$ +' log_k -15'
|
||||
570 tracer_phases$ = tracer_phases$ + nl$ + ' A_Na_tr' + nl$ + ' Na_trCl = Na_tr+ + Cl-' + sc$ +' log_k -14'
|
||||
580 tracer_phases$ = tracer_phases$ + nl$ + ' A_Cl_tr' + nl$ + ' NaCl_tr = Na+ + Cl_tr-' + sc$ +' log_k -14'
|
||||
590 tracer_phases$ = tracer_phases$ + nl$ + ' A_Cs' + nl$ + ' CsCl = Cs+ + Cl-' + sc$ + ' log_k -13'
|
||||
600 DIM tracer_equi$(4)
|
||||
610 FOR i = 1 TO 4
|
||||
620 tracer_equi$(i) = nl$ + 'A_' + tracer$(i) + ' 0 0'
|
||||
630 NEXT i
|
||||
|
||||
# Write solutions for the cells...
|
||||
650 punch nl$ + 'PRINT ' + sc$ + ' -reset false' + sc$ + ' -echo_input true' + sc$ + ' -user_print true'
|
||||
660 IF nfilt1 = 0 THEN GOTO 800
|
||||
670 punch nl$ + x$ + ' filter cells at tracer-in side...'
|
||||
680 r1 = r_int - thickn_filter1
|
||||
690 xf1 = thickn_filter1 / nfilt1
|
||||
700 FOR i = 1 TO nfilt1
|
||||
710 num$ = TRIM(STR$(i + 3)) + sc$
|
||||
720 V_water = 1e3 * height * por_filter1 * pi * (SQR(r1 + xf1) - SQR(r1))
|
||||
730 punch nl$ + 'SOLUTION ' + num$ + ' -water ' + STR$(V_water)
|
||||
740 punch sol$ + nl$
|
||||
750 r1 = r1 + xf1
|
||||
760 NEXT i
|
||||
|
||||
800 punch nl$ + nl$ + x$ + ' cells in Opalinus Clay...'
|
||||
810 r1 = r_int
|
||||
820 x = thickn_clay / nclay
|
||||
830 FOR i = 1 TO nclay
|
||||
840 num$ = TRIM(STR$(i + 3 + nfilt1)) + sc$
|
||||
850 V_water = 1e3 * height * por_clay * pi * (SQR(r1 + x) - SQR(r1))
|
||||
860 punch nl$ + 'SOLUTION ' + num$ + ' -water ' + STR$(V_water * f_free)
|
||||
870 punch sol$
|
||||
880 IF f_free = 1 and tracer = 1 THEN GOTO 960
|
||||
890 punch nl$ + 'SURFACE ' + num$ + ' -equil ' + num$
|
||||
900 punch nl$ + ' Su_ ' + TRIM(STR$(f_DL_charge * CEC * V_water)) + STR$(A_por) + ' ' + STR$(V_water)
|
||||
910 punch nl$ + ' Su_ii ' + TRIM(STR$(7.88e-4 * rho_b_eps * V_water))
|
||||
920 punch nl$ + ' Su_fes ' + TRIM(STR$(7.4e-5 * rho_b_eps * V_water))
|
||||
930 IF f_free < 1 THEN punch nl$ + ' -Donnan ' + TRIM(STR$((1 - f_free) * 1e-3 / A_por))
|
||||
940 punch nl$ + 'EXCHANGE ' + num$ + ' -equil ' + num$
|
||||
950 punch nl$ + ' X ' + TRIM(STR$((1 - f_DL_charge) * CEC * V_water)) + nl$
|
||||
960 r1 = r1 + x
|
||||
970 NEXT i
|
||||
|
||||
1000 IF nfilt2 = 0 THEN GOTO 1200
|
||||
1010 punch nl$ + nl$ + x$ + ' tracer-out filter cells...'
|
||||
1020 r1 = r_ext
|
||||
1030 xf2 = thickn_filter2 / nfilt2
|
||||
1040 FOR i = 1 TO nfilt2
|
||||
1050 num$ = TRIM(STR$(i + 3 + nfilt1 + nclay)) + sc$
|
||||
1060 V_water = 1e3 * height * por_filter2 * pi * (SQR(r1 + xf2) - SQR(r1))
|
||||
1070 punch nl$ + 'SOLUTION ' + num$ + ' -water ' + STR$(V_water)
|
||||
1080 punch sol$ + nl$
|
||||
1090 r1 = r1 + xf2
|
||||
1100 NEXT i
|
||||
|
||||
1200 punch nl$ + x$ + ' outside solution...'
|
||||
1210 num$ = TRIM(STR$(4 + nfilt1 + nclay + nfilt2)) + sc$
|
||||
1220 punch nl$ + 'SOLUTION ' + num$ + ' -water ' + STR$(V_end)
|
||||
1230 punch sol$
|
||||
1240 punch nl$ + 'END'
|
||||
|
||||
# Write phases in which the tracers precipitate...
|
||||
1300 punch nl$ + tracer_phases$
|
||||
1310 punch nl$ + 'EQUILIBRIUM_PHASES ' + num$ + tracer_equi$(tracer)
|
||||
1312 If tracer = 3 THEN punch nl$ + tracer_equi$(tracer + 1)
|
||||
1320 punch nl$ + 'END'
|
||||
|
||||
# Define mixing factors for the diffusive flux between cells 1 and 2:
|
||||
# J_12 = -2 * Dw / (x_1 / g_1 + x_2 / g_2) * (c_2 - c_1)
|
||||
# Multiply with dt * A / (V = 1e-3 m3). (Actual volumes are given with SOLUTION; -water)
|
||||
# Use harmonic mean: g_1 = por_1 / G_1, g_2 = por_2 / G_2, x_1 = Delta(x_1), etc.
|
||||
1400 IF nfilt1 > 0 THEN gf1 = por_filter1 / G_filter1
|
||||
1410 IF nfilt2 > 0 THEN gf2 = por_filter2 / G_filter2
|
||||
1420 g = por_clay / G_clay
|
||||
# Find max time step = 0.5 * V_water * dx * G_factor / (Dw * por * A * fbc)
|
||||
# V_water = por * pi * height * ((r + dr)^2 - r^2)
|
||||
# A = por * pi * height * r * 2
|
||||
# At the inlet of the tracers, fbc = 2...
|
||||
1500 IF nfilt1 = 0 THEN GOTO 1530
|
||||
1510 r1 = r_int - thickn_filter1
|
||||
1520 ff = (SQR(r1 + xf1) - SQR(r1)) * xf1 * G_filter1 / (r1 * 2) / 2
|
||||
1530 ff1 = (SQR(r_int + x) - SQR(r_int)) * x * G_clay / (r_int * 2) / 2
|
||||
# Perhaps the clay has very small cells...
|
||||
1540 IF nfilt1 = 0 THEN ff = ff1 ELSE IF ff1 * 2 < ff THEN ff = ff1 * 2
|
||||
# Or at the filter1-clay transition, fbc = 1...
|
||||
1550 IF nfilt1 > 0 THEN ff1 = (SQR(r_int + x) - SQR(r_int)) * (xf1 / gf1 + x / g) / (2 * r_int * 2)
|
||||
1560 IF nfilt1 > 0 AND ff1 < ff THEN ff = ff1
|
||||
# Perhaps filter2 has very small cells...
|
||||
1570 IF nfilt2 > 0 THEN ff1 = (SQR(r_ext + xf2) - SQR(r_ext)) * xf2 * G_filter2 / (r_ext * 2)
|
||||
1580 IF nfilt2 > 0 AND ff1 < ff THEN ff = ff1
|
||||
1590 dt_max = 0.5 * ff / Dw
|
||||
# Check with punch times, set shifts...
|
||||
1610 IF punch_time < dt_max THEN dt = punch_time ELSE dt = dt_max
|
||||
1620 punch_fr = 1
|
||||
1630 IF dt < punch_time THEN punch_fr = ceil(punch_time / dt)
|
||||
1640 dt = punch_time / punch_fr
|
||||
1650 shifts = ceil(exp_time(tracer) / dt)
|
||||
# Write mixing factors...
|
||||
1700 punch nl$ + nl$ + x$ + ' mixing factors...'
|
||||
1710 r1 = r_int
|
||||
1720 IF nfilt1 > 0 THEN r1 = r_int - thickn_filter1
|
||||
1730 A = height * 2 * pi
|
||||
1740 FOR i = 0 TO nfilt1 + nclay + nfilt2
|
||||
1750 IF i = 0 OR i = nfilt1 + nclay + nfilt2 THEN fbc = 2 ELSE fbc = 1
|
||||
1760 IF i > nfilt1 OR nfilt1 = 0 THEN GOTO 1810
|
||||
1770 IF i < nfilt1 THEN mixf = Dw * fbc / (xf1 / gf1) * dt * A * r1 / 1e-3
|
||||
1780 IF i = nfilt1 THEN mixf = 2 * Dw / (xf1 / gf1 + x / g) * dt * A * r1 / 1e-3
|
||||
1790 IF i < nfilt1 THEN r1 = r1 + xf1 ELSE r1 = r1 + x
|
||||
1800 GOTO 1880
|
||||
1810 IF i > nfilt1 + nclay THEN GOTO 1860
|
||||
1820 mixf = Dw * fbc / (x / g) * dt * A * r1 / 1e-3
|
||||
1830 IF i = nfilt1 + nclay AND nfilt2 > 0 THEN mixf = 2 * Dw / (xf2 / gf2 + x / g) * dt * A * r1 / 1e-3
|
||||
1840 IF i < nfilt1 + nclay THEN r1 = r1 + x ELSE r1 = r1 + xf2
|
||||
1850 GOTO 1880
|
||||
1860 mixf = Dw * fbc / (xf2 / gf2) * dt * A * r1 / 1e-3
|
||||
1870 r1 = r1 + xf2
|
||||
1880 punch nl$ + 'MIX ' + TRIM(STR$(i + 3)) + sc$ + STR$(i + 4) + STR$(mixf)
|
||||
1890 NEXT i
|
||||
1900 punch nl$ + 'END'
|
||||
|
||||
# Write TRANSPORT...
|
||||
2000 punch nl$ + 'TRANSPORT'
|
||||
2010 stag = 2 + nfilt1 + nclay + nfilt2
|
||||
2020 punch nl$ + ' -warnings true'
|
||||
2030 punch nl$ + ' -shifts ' + TRIM(STR$(shifts))
|
||||
2040 punch nl$ + ' -flow diff' + sc$ + ' -cells 1' + sc$ + ' -bcon 1 2' + sc$ + ' -stag ' + TRIM(STR$(stag))
|
||||
2050 punch nl$ + ' -time ' + STR$(dt)
|
||||
2060 punch nl$ + ' -multi_D true ' + STR$(Dw) + STR$(por_clay) + ' 0.0 ' + TRIM(STR$(-tort_n))
|
||||
2070 punch nl$ + ' -interlayer_D ' + interlayer_D$ + ' 0.001 0.0 ' + TRIM(STR$(G_IL))
|
||||
2080 punch nl$ + ' -punch_fr ' + TRIM(STR$(punch_fr)) + sc$ + ' -punch_c ' + TRIM(STR$(2 + stag))
|
||||
|
||||
# Write USER_GRAPH...
|
||||
2180 FOR i = 0 to 1
|
||||
2190 punch nl$ + 'USER_GRAPH ' + TRIM(STR$(tracer + i)) + ' Example 21' + nl$
|
||||
2200 punch nl$ + ' -chart_title " ' + tracer$(tracer + i) + ' Diffusion to Outer Cell"'
|
||||
2210 punch nl$ + ' -plot_tsv_file ex21_' + tracer$(tracer + i) + '_rad.tsv'
|
||||
2220 punch nl$ + ' -axis_scale x_axis 0 ' + TRIM(STR$(exp_time(tracer + i) / (3600 * 24)))
|
||||
2230 punch nl$ + ' -axis_titles "TIME, IN DAYS" "FLUX, IN MOL/M2/S" "ACCUMULATED MASS, IN MOL"'
|
||||
2240 punch nl$ + ' -plot_concentration_vs time'
|
||||
2250 punch nl$ + ' 10 days = total_time / (3600 * 24)'
|
||||
2260 punch nl$ + ' 20 a = equi("A_' + tracer$(tracer + i) + '")'
|
||||
2270 punch nl$ + ' 30 IF get(1) = 0 AND total_time > 0 THEN put(total_time, 1)'
|
||||
2280 punch nl$ + ' 40 dt = get(1)'
|
||||
2290 A = 2 * pi * r_ext * height
|
||||
2300 i$ = TRIM(STR$(2 + i))
|
||||
2310 punch nl$ + ' 50 plot_xy days - dt / (2 * 3600 * 24), (a - get(' + i$ + ')) / dt /' + STR$(A) + \
|
||||
', color = Green, symbol = None'
|
||||
2320 punch nl$ + ' 60 put(a, ' + i$ + ')'
|
||||
2330 punch nl$ + ' 70 plot_xy days, equi("A_' + tracer$(tracer + i) + \
|
||||
'"), y_axis = 2, color = Red, symbol = None'
|
||||
2340 IF tracer < 3 THEN GOTO 2360
|
||||
2350 NEXT i
|
||||
2360 punch nl$ + 'END'
|
||||
|
||||
2400 IF profile$ = 'true' THEN GOSUB 3000
|
||||
2410 IF tracer < 3 THEN END # finished for Hto and Cl
|
||||
|
||||
# Continue with Cs...
|
||||
2420 IF profile$ = 'false' THEN punch nl$ + 'USER_GRAPH ' + TRIM(STR$(tracer)) + sc$ + ' -detach' ELSE \
|
||||
punch nl$ + 'USER_GRAPH ' + TRIM(STR$(tracer + 4)) + sc$ + ' -detach'
|
||||
2440 tracer = tracer + 1
|
||||
2450 punch nl$ + 'TRANSPORT'
|
||||
2460 shifts = ceil((exp_time(tracer) - exp_time(tracer - 1))/ dt)
|
||||
2480 punch nl$ + ' -shifts ' + TRIM(STR$(shifts))
|
||||
2490 punch nl$ + ' -punch_fr ' + TRIM(STR$(punch_fr)) + sc$ + ' -punch_c ' + TRIM(STR$(2 + stag))
|
||||
2500 punch nl$ + 'END'
|
||||
2510 IF profile$ = 'true' THEN GOSUB 3000
|
||||
2520 END # finished...
|
||||
|
||||
# Write TRANSPORT and USER_GRAPH for concentration profile...
|
||||
3000 punch nl$ + 'TRANSPORT'
|
||||
3010 punch nl$ + ' -shifts 0'
|
||||
3020 punch nl$ + ' -punch_fr 2' + sc$ + ' -punch_c 3-' + TRIM(STR$(2 + stag))
|
||||
# Write USER_GRAPH...
|
||||
3030 punch nl$ + 'USER_GRAPH ' + TRIM(STR$(tracer)) + sc$ + ' -detach'
|
||||
3040 punch nl$ + 'USER_GRAPH ' + TRIM(STR$(tracer + 4)) + ' Example 21' + nl$
|
||||
3050 punch nl$ + ' -chart_title "' + tracer$(tracer) + ' Concentration Profile: Filter1 | Clay | Filter2"'
|
||||
3060 REM punch nl$ + ' -plot_tsv_file + tracer$(tracer) + '_prof.tsv'
|
||||
3070 punch nl$ + ' -axis_scale x_axis 0 ' + \
|
||||
TRIM(STR$((thickn_filter1 + thickn_clay + thickn_filter2) * 1e3))
|
||||
3080 punch nl$ + ' -axis_scale y_axis ' + profile_y1$(tracer)
|
||||
3090 punch nl$ + ' -axis_scale sy_axis ' + profile_y2$(tracer)
|
||||
3100 punch nl$ + ' -axis_titles ' + '"DISTANCE, IN MILLIMETERS" "FREE PORE-WATER MOLALITY" "TOTAL MOLALITY"'
|
||||
3110 punch nl$ + ' -headings ' + tracer$(tracer) + '_free ' + tracer$(tracer) + '_tot'
|
||||
3120 punch nl$ + ' -plot_concentration_vs x'
|
||||
3130 punch nl$ + ' -initial_solutions true'
|
||||
3140 punch nl$ + ' 10 IF cell_no = 3 THEN xval = 0 ELSE xval = get(14)'
|
||||
3150 punch nl$ + ' 20 IF (' + TRIM(STR$(nfilt1)) + ' = 0 OR cell_no > ' + \
|
||||
TRIM(STR$(nfilt1 + 3)) + ') THEN GOTO 60'
|
||||
3160 punch nl$ + ' 30 IF (cell_no = 4) THEN xval = xval + 0.5 * ' + TRIM(STR$(xf1))
|
||||
3170 punch nl$ + ' 40 IF (cell_no > 4 AND cell_no < ' + \
|
||||
TRIM(STR$(nfilt1 + 4)) + ') THEN xval = xval + ' + TRIM(STR$(xf1))
|
||||
3180 punch nl$ + ' 50 GOTO 200'
|
||||
3190 punch nl$ + ' 60 IF (cell_no = ' + TRIM(STR$(4 + nfilt1)) + ') THEN \
|
||||
xval = xval + 0.5 * ' + TRIM(STR$(xf1)) + ' + 0.5 * ' + TRIM(STR$(x))
|
||||
3200 punch nl$ + ' 70 IF (cell_no > ' + TRIM(STR$(4 + nfilt1)) + ' AND cell_no < ' + \
|
||||
TRIM(STR$(4 + nfilt1 + nclay)) + ') THEN xval = xval + ' + TRIM(STR$(x)) + ' ELSE GOTO 90'
|
||||
3210 punch nl$ + ' 80 GOTO 200'
|
||||
3220 punch nl$ + ' 90 IF (cell_no = ' + TRIM(STR$(4 + nfilt1 + nclay)) + ') THEN \
|
||||
xval = xval + 0.5 * ' + TRIM(STR$(x)) + ' + 0.5 * ' + TRIM(STR$(xf2))
|
||||
3230 punch nl$ + ' 100 IF (cell_no > ' + TRIM(STR$(4 + nfilt1 + nclay)) + ' AND cell_no <= ' + \
|
||||
TRIM(STR$(3 + nfilt1 + nclay + nfilt2)) + ') THEN xval = xval + ' + TRIM(STR$(xf2))
|
||||
3240 punch nl$ + ' 110 IF (cell_no = ' + TRIM(STR$(4 + nfilt1 + nclay + nfilt2)) + ') THEN \
|
||||
xval = xval + 0.5 * ' + TRIM(STR$(xf2))
|
||||
3250 punch nl$ + ' 200 y1 = TOT("' + tracer$(tracer) + '")'
|
||||
3260 punch nl$ + ' 210 plot_xy xval * 1e3, y1, color = Blue, symbol = Plus'
|
||||
3270 punch nl$ + ' 220 IF cell_no = 3 THEN put(y1, 15)'
|
||||
3280 punch nl$ + ' 230 IF (cell_no < ' + TRIM(STR$(4 + nfilt1)) + ' OR cell_no > ' + \
|
||||
TRIM(STR$(3 + nfilt1 + nclay)) + ') THEN GOTO 400'
|
||||
3290 punch nl$ + ' 240 y2 = SYS("' + tracer$(tracer) + '") / (tot("water") + edl("water"))'
|
||||
# Remove REM if total conc's per kg solid must be plotted (and adapt axis_titles)...
|
||||
3310 punch nl$ + ' 250 REM y2 = y2 / ' + TRIM(STR$(rho_b_eps)) + x$ + ' conc / kg solid'
|
||||
3320 punch nl$ + ' 260 plot_xy xval * 1e3, y2, symbol = Circle, y_axis = 2'
|
||||
3330 punch nl$ + ' 270 IF (cell_no > ' + TRIM(STR$(5 + nfilt1)) + ') THEN GOTO 400'
|
||||
3340 punch nl$ + ' 280 IF ' + TRIM(STR$(nfilt1)) + ' THEN plot_xy ' + \
|
||||
TRIM(STR$(thickn_filter1 * 1e3)) + ', get(15), color = Black, symbol = None'
|
||||
3350 punch nl$ + ' 290 IF ' + TRIM(STR$(nfilt2)) + ' THEN plot_xy ' + \
|
||||
TRIM(STR$((thickn_filter1 + thickn_clay) * 1e3)) + \
|
||||
', get(15), color = Black, symbol = None'
|
||||
3360 punch nl$ + ' 300 put(0, 15)'
|
||||
3370 punch nl$ + ' 400 put(xval, 14)'
|
||||
3380 punch nl$ + 'END'
|
||||
3390 RETURN
|
||||
END
|
||||
PRINT
|
||||
-selected_out false
|
||||
INCLUDE$ radial
|
||||
END
|
||||
971
ex21.out
Normal file
971
ex21.out
Normal file
@ -0,0 +1,971 @@
|
||||
Input file: ../examples/ex21_radial
|
||||
Output file: ex21_radial.out
|
||||
Database file: ../database/phreeqc.dat
|
||||
|
||||
------------------
|
||||
Reading data base.
|
||||
------------------
|
||||
|
||||
SOLUTION_MASTER_SPECIES
|
||||
SOLUTION_SPECIES
|
||||
PHASES
|
||||
EXCHANGE_MASTER_SPECIES
|
||||
EXCHANGE_SPECIES
|
||||
SURFACE_MASTER_SPECIES
|
||||
SURFACE_SPECIES
|
||||
RATES
|
||||
END
|
||||
------------------------------------
|
||||
Reading input data for simulation 1.
|
||||
------------------------------------
|
||||
|
||||
TITLE Diffusion through Opalinus Clay in a radial diffusion cell, Appelo, Van Loon and Wersin,
|
||||
2010, GCA 74, 1201
|
||||
SOLUTION_MASTER_SPECIES
|
||||
Hto Hto 0.0 20 20
|
||||
Na_tr Na_tr+ 0.0 22 22
|
||||
Cl_tr Cl_tr- 0.0 36 36
|
||||
Cs Cs+ 0.0 132.905 132.905
|
||||
SOLUTION_SPECIES
|
||||
Hto = Hto
|
||||
log_k 0
|
||||
gamma 1e6 0
|
||||
dw 2.236e-9
|
||||
Na_tr+ = Na_tr+
|
||||
log_k 0
|
||||
gamma 4.0 0.075
|
||||
dw 1.33e-9
|
||||
erm_ddl 1.23
|
||||
Cl_tr- = Cl_tr-
|
||||
log_k 0
|
||||
gamma 3.5 0.015
|
||||
dw 1.31e-9 # dw = dw(water) / 1.55 = 2.03e-9 / 1.55
|
||||
Cs+ = Cs+
|
||||
log_k 0
|
||||
gamma 3.5 0.015
|
||||
dw 2.07e-9
|
||||
erm_ddl 1.23
|
||||
SURFACE_MASTER_SPECIES
|
||||
Su_fes Su_fes- # Frayed Edge Sites
|
||||
Su_ii Su_ii- # Type II sites of intermediate strength
|
||||
Su_ Su_- # Double layer, planar sites are modeled with EXCHANGE
|
||||
SURFACE_SPECIES
|
||||
Su_fes- = Su_fes-
|
||||
log_k 0
|
||||
Na+ + Su_fes- = NaSu_fes
|
||||
log_k 10
|
||||
Na_tr+ + Su_fes- = Na_trSu_fes
|
||||
log_k 10
|
||||
K+ + Su_fes- = KSu_fes
|
||||
log_k 12.4
|
||||
Cs+ + Su_fes- = CsSu_fes
|
||||
log_k 17.14
|
||||
Su_ii- = Su_ii-
|
||||
log_k 0
|
||||
Na+ + Su_ii- = NaSu_ii
|
||||
log_k 10
|
||||
Na_tr+ + Su_ii- = Na_trSu_ii
|
||||
log_k 10
|
||||
K+ + Su_ii- = KSu_ii
|
||||
log_k 12.1
|
||||
Cs+ + Su_ii- = CsSu_ii
|
||||
log_k 14.6
|
||||
Su_- = Su_-
|
||||
log_k 0
|
||||
EXCHANGE_SPECIES
|
||||
Na_tr+ + X- = Na_trX
|
||||
log_k 0.0
|
||||
gamma 4.0 0.075
|
||||
Cs+ + X- = CsX
|
||||
log_k 2.04
|
||||
gamma 3.5 0.015
|
||||
SOLUTION 0-2 column with only cell 1, two boundary solutions 0 and 2.
|
||||
Na 1
|
||||
Cl 1
|
||||
END
|
||||
-----
|
||||
TITLE
|
||||
-----
|
||||
|
||||
Diffusion through Opalinus Clay in a radial diffusion cell, Appelo, Van Loon and Wersin,
|
||||
2010, GCA 74, 1201
|
||||
|
||||
-------------------------------------------
|
||||
Beginning of initial solution calculations.
|
||||
-------------------------------------------
|
||||
|
||||
Initial solution 0. column with only cell 1, two boundary solutions 0 and 2.
|
||||
|
||||
-----------------------------Solution composition------------------------------
|
||||
|
||||
Elements Molality Moles
|
||||
|
||||
Cl 1.000e-03 1.000e-03
|
||||
Na 1.000e-03 1.000e-03
|
||||
|
||||
----------------------------Description of solution----------------------------
|
||||
|
||||
pH = 7.000
|
||||
pe = 4.000
|
||||
Specific Conductance (uS/cm, 25 oC) = 123
|
||||
Density (g/cm3) = 0.99708
|
||||
Activity of water = 1.000
|
||||
Ionic strength = 1.000e-03
|
||||
Mass of water (kg) = 1.000e+00
|
||||
Total alkalinity (eq/kg) = 4.191e-10
|
||||
Total carbon (mol/kg) = 0.000e+00
|
||||
Total CO2 (mol/kg) = 0.000e+00
|
||||
Temperature (deg C) = 25.00
|
||||
Electrical balance (eq) = -4.191e-10
|
||||
Percent error, 100*(Cat-|An|)/(Cat+|An|) = -0.00
|
||||
Iterations = 3
|
||||
Total H = 1.110124e+02
|
||||
Total O = 5.550622e+01
|
||||
|
||||
----------------------------Distribution of species----------------------------
|
||||
|
||||
Log Log Log mole V
|
||||
Species Molality Activity Molality Activity Gamma cm3/mol
|
||||
|
||||
OH- 1.039e-07 1.002e-07 -6.983 -6.999 -0.016 -4.08
|
||||
H+ 1.035e-07 1.000e-07 -6.985 -7.000 -0.015 0.00
|
||||
H2O 5.551e+01 1.000e+00 1.744 -0.000 0.000 18.07
|
||||
Cl 1.000e-03
|
||||
Cl- 1.000e-03 9.648e-04 -3.000 -3.016 -0.016 18.06
|
||||
H(0) 1.416e-25
|
||||
H2 7.078e-26 7.079e-26 -25.150 -25.150 0.000 28.61
|
||||
Na 1.000e-03
|
||||
Na+ 1.000e-03 9.651e-04 -3.000 -3.015 -0.015 -1.20
|
||||
NaOH 9.670e-21 9.672e-21 -20.015 -20.014 0.000 (0)
|
||||
O(0) 0.000e+00
|
||||
O2 0.000e+00 0.000e+00 -42.080 -42.080 0.000 30.40
|
||||
|
||||
------------------------------Saturation indices-------------------------------
|
||||
|
||||
Phase SI log IAP log K(298 K, 1 atm)
|
||||
|
||||
H2(g) -22.05 -25.15 -3.10 H2
|
||||
H2O(g) -1.50 -0.00 1.50 H2O
|
||||
Halite -7.61 -6.03 1.58 NaCl
|
||||
O2(g) -39.19 -42.08 -2.89 O2
|
||||
|
||||
|
||||
------------------
|
||||
End of simulation.
|
||||
------------------
|
||||
|
||||
------------------------------------
|
||||
Reading input data for simulation 2.
|
||||
------------------------------------
|
||||
|
||||
KNOBS
|
||||
diagonal_scale true #; -tolerance 1e-17 # if necessary, try lower tolereance
|
||||
SOLUTION 3 tracer solution
|
||||
pH 7.6
|
||||
pe 14 O2(g) -1.0
|
||||
temp 23
|
||||
Na 240
|
||||
K 1.61
|
||||
Mg 16.9
|
||||
Ca 25.8
|
||||
Sr 0.505
|
||||
Cl 300
|
||||
S(6) 14.1
|
||||
Fe(2) 0.0
|
||||
Alkalinity 0.476
|
||||
Hto 1.14e-6
|
||||
water 0.2
|
||||
SELECTED_OUTPUT
|
||||
file radial
|
||||
reset false
|
||||
USER_PUNCH
|
||||
1 nl$ = EOL$ # newline
|
||||
2 x$ = CHR$(35) # cross '#'
|
||||
3 sc$ = CHR$(59) # semicolon ';'
|
||||
4 pi = 2 * ARCTAN(1e10) # 3.14159...
|
||||
10 height = 0.052 # length of the clay cylinder / m
|
||||
20 r_int = 6.58e-3 # inner radius of clay cylinder / m
|
||||
30 r_ext = 25.4e-3 # outer radius
|
||||
40 thickn_filter1 = 1.8e-3 # tracer-in filter thickness / m
|
||||
50 thickn_filter2 = 1.6e-3 # tracer-out filter thickness / m
|
||||
60 por_filter1 = 0.418 # porosity
|
||||
70 por_filter2 = 0.367
|
||||
80 G_filter1 = 4.18 # geometrical factor. (for filters, G = por / 10)
|
||||
90 G_filter2 = 3.67
|
||||
100 V_end = 0.2 # volume of the tracer-out solution / L
|
||||
110 thickn_clay = r_ext - r_int # clay thickness / m
|
||||
120 por_clay = 0.159
|
||||
130 rho_b_eps = 2.7 * (1 - por_clay) / por_clay # clay bulk density / porosity / (kg/L)
|
||||
140 CEC = 0.12 * rho_b_eps # CEC / (eq/L porewater)
|
||||
150 A_por = 37e3 * rho_b_eps # pore surface area / (m2/L porewater)
|
||||
160 DIM tracer$(4), exp_time(4), scale_y1$(4), scale_y2$(4), profile_y1$(4), profile_y2$(4)
|
||||
170 DATA 'Hto', 'Cl_tr', 'Na_tr', 'Cs'
|
||||
180 READ tracer$(1), tracer$(2), tracer$(3), tracer$(4)
|
||||
200 DATA 86400 * 20, 86400 * 40, 86400 * 45, 86400 * 1000
|
||||
210 READ exp_time(1), exp_time(2), exp_time(3), exp_time(4)
|
||||
230 DATA '1', '1', '1', '1'
|
||||
240 READ scale_y1$(1), scale_y1$(2), scale_y1$(3), scale_y1$(4)
|
||||
260 DATA '1', '1', '1', '1'
|
||||
270 READ scale_y2$(1), scale_y2$(2), scale_y2$(3), scale_y2$(4)
|
||||
280 DATA '0 1.2e-9', '0 2.5e-5', '0 2e-10', '0 auto'
|
||||
290 READ profile_y1$(1), profile_y1$(2), profile_y1$(3), profile_y1$(4)
|
||||
300 DATA '0 1.2e-9', '0 2.5e-5', '0 6e-10', '0 auto'
|
||||
310 READ profile_y2$(1), profile_y2$(2), profile_y2$(3), profile_y2$(4)
|
||||
350 Dw = 2.5e-9 # default tracer diffusion coefficient / (m2/s)
|
||||
360 nfilt1 = 1 # number of cells in filter 1
|
||||
370 nfilt2 = 1 # number of cells in filter 2
|
||||
380 nclay = 11 # number of clay cells
|
||||
390 f_free = 0.117 # fraction of free pore water (0.01 - 1)
|
||||
400 f_DL_charge = 0.45 # fraction of CEC charge in electrical double layer
|
||||
410 tort_n = -0.99 # exponent in Archie's law, -1.045 without filters
|
||||
420 G_clay = por_clay^tort_n # geometrical factor
|
||||
430 interlayer_D$ = 'false' # 'true' or 'false' for interlayer diffusion
|
||||
440 G_IL = 700 # geometrical factor for clay interlayers
|
||||
450 punch_time = 60 * 60 * 6 # punch time / seconds
|
||||
460 profile$ = 'true' # 'true' or 'false' for c/x profile visualization
|
||||
470 IF nfilt1 = 0 THEN thickn_filter1 = 0
|
||||
480 IF nfilt2 = 0 THEN thickn_filter2 = 0
|
||||
490 IF tot("Hto") > 1e-10 THEN tracer = 1 ELSE IF tot("Cl_tr") > 1e-10 THEN tracer = 2 ELSE tracer = 3
|
||||
520 sol$ = nl$ + ' pH 7.6' + sc$ +' pe 14 O2(g) -1.0' + sc$ +' temp 23'
|
||||
530 sol$ = sol$ + nl$ + ' Na 240' + sc$ +' K 1.61' + sc$ +' Mg 16.9' + sc$ +' Ca 25.8' + sc$ +' Sr 0.505'
|
||||
540 sol$ = sol$ + nl$ + ' Cl 300' + sc$ +' S(6) 14.1' + sc$ +' Fe(2) 0.0' + sc$ +' Alkalinity 0.476'
|
||||
550 tracer_phases$ = nl$ + 'PHASES '
|
||||
560 tracer_phases$ = tracer_phases$ + nl$ + ' A_Hto' + nl$ + ' Hto = Hto' + sc$ +' log_k -15'
|
||||
570 tracer_phases$ = tracer_phases$ + nl$ + ' A_Na_tr' + nl$ + ' Na_trCl = Na_tr+ + Cl-' + sc$ +' log_k -14'
|
||||
580 tracer_phases$ = tracer_phases$ + nl$ + ' A_Cl_tr' + nl$ + ' NaCl_tr = Na+ + Cl_tr-' + sc$ +' log_k -14'
|
||||
590 tracer_phases$ = tracer_phases$ + nl$ + ' A_Cs' + nl$ + ' CsCl = Cs+ + Cl-' + sc$ + ' log_k -13'
|
||||
600 DIM tracer_equi$(4)
|
||||
610 FOR i = 1 TO 4
|
||||
620 tracer_equi$(i) = nl$ + 'A_' + tracer$(i) + ' 0 0'
|
||||
630 NEXT i
|
||||
650 punch nl$ + 'PRINT ' + sc$ + ' -reset false' + sc$ + ' -echo_input true' + sc$ + ' -user_print true'
|
||||
660 IF nfilt1 = 0 THEN GOTO 800
|
||||
670 punch nl$ + x$ + ' filter cells at tracer-in side...'
|
||||
680 r1 = r_int - thickn_filter1
|
||||
690 xf1 = thickn_filter1 / nfilt1
|
||||
700 FOR i = 1 TO nfilt1
|
||||
710 num$ = TRIM(STR$(i + 3)) + sc$
|
||||
720 V_water = 1e3 * height * por_filter1 * pi * (SQR(r1 + xf1) - SQR(r1))
|
||||
730 punch nl$ + 'SOLUTION ' + num$ + ' -water ' + STR$(V_water)
|
||||
740 punch sol$ + nl$
|
||||
750 r1 = r1 + xf1
|
||||
760 NEXT i
|
||||
800 punch nl$ + nl$ + x$ + ' cells in Opalinus Clay...'
|
||||
810 r1 = r_int
|
||||
820 x = thickn_clay / nclay
|
||||
830 FOR i = 1 TO nclay
|
||||
840 num$ = TRIM(STR$(i + 3 + nfilt1)) + sc$
|
||||
850 V_water = 1e3 * height * por_clay * pi * (SQR(r1 + x) - SQR(r1))
|
||||
860 punch nl$ + 'SOLUTION ' + num$ + ' -water ' + STR$(V_water * f_free)
|
||||
870 punch sol$
|
||||
880 IF f_free = 1 and tracer = 1 THEN GOTO 960
|
||||
890 punch nl$ + 'SURFACE ' + num$ + ' -equil ' + num$
|
||||
900 punch nl$ + ' Su_ ' + TRIM(STR$(f_DL_charge * CEC * V_water)) + STR$(A_por) + ' ' + STR$(V_water)
|
||||
910 punch nl$ + ' Su_ii ' + TRIM(STR$(7.88e-4 * rho_b_eps * V_water))
|
||||
920 punch nl$ + ' Su_fes ' + TRIM(STR$(7.4e-5 * rho_b_eps * V_water))
|
||||
930 IF f_free < 1 THEN punch nl$ + ' -Donnan ' + TRIM(STR$((1 - f_free) * 1e-3 / A_por))
|
||||
940 punch nl$ + 'EXCHANGE ' + num$ + ' -equil ' + num$
|
||||
950 punch nl$ + ' X ' + TRIM(STR$((1 - f_DL_charge) * CEC * V_water)) + nl$
|
||||
960 r1 = r1 + x
|
||||
970 NEXT i
|
||||
1000 IF nfilt2 = 0 THEN GOTO 1200
|
||||
1010 punch nl$ + nl$ + x$ + ' tracer-out filter cells...'
|
||||
1020 r1 = r_ext
|
||||
1030 xf2 = thickn_filter2 / nfilt2
|
||||
1040 FOR i = 1 TO nfilt2
|
||||
1050 num$ = TRIM(STR$(i + 3 + nfilt1 + nclay)) + sc$
|
||||
1060 V_water = 1e3 * height * por_filter2 * pi * (SQR(r1 + xf2) - SQR(r1))
|
||||
1070 punch nl$ + 'SOLUTION ' + num$ + ' -water ' + STR$(V_water)
|
||||
1080 punch sol$ + nl$
|
||||
1090 r1 = r1 + xf2
|
||||
1100 NEXT i
|
||||
1200 punch nl$ + x$ + ' outside solution...'
|
||||
1210 num$ = TRIM(STR$(4 + nfilt1 + nclay + nfilt2)) + sc$
|
||||
1220 punch nl$ + 'SOLUTION ' + num$ + ' -water ' + STR$(V_end)
|
||||
1230 punch sol$
|
||||
1240 punch nl$ + 'END'
|
||||
1300 punch nl$ + tracer_phases$
|
||||
1310 punch nl$ + 'EQUILIBRIUM_PHASES ' + num$ + tracer_equi$(tracer)
|
||||
1312 If tracer = 3 THEN punch nl$ + tracer_equi$(tracer + 1)
|
||||
1320 punch nl$ + 'END'
|
||||
1400 IF nfilt1 > 0 THEN gf1 = por_filter1 / G_filter1
|
||||
1410 IF nfilt2 > 0 THEN gf2 = por_filter2 / G_filter2
|
||||
1420 g = por_clay / G_clay
|
||||
1500 IF nfilt1 = 0 THEN GOTO 1530
|
||||
1510 r1 = r_int - thickn_filter1
|
||||
1520 ff = (SQR(r1 + xf1) - SQR(r1)) * xf1 * G_filter1 / (r1 * 2) / 2
|
||||
1530 ff1 = (SQR(r_int + x) - SQR(r_int)) * x * G_clay / (r_int * 2) / 2
|
||||
1540 IF nfilt1 = 0 THEN ff = ff1 ELSE IF ff1 * 2 < ff THEN ff = ff1 * 2
|
||||
1550 IF nfilt1 > 0 THEN ff1 = (SQR(r_int + x) - SQR(r_int)) * (xf1 / gf1 + x / g) / (2 * r_int * 2)
|
||||
1560 IF nfilt1 > 0 AND ff1 < ff THEN ff = ff1
|
||||
1570 IF nfilt2 > 0 THEN ff1 = (SQR(r_ext + xf2) - SQR(r_ext)) * xf2 * G_filter2 / (r_ext * 2)
|
||||
1580 IF nfilt2 > 0 AND ff1 < ff THEN ff = ff1
|
||||
1590 dt_max = 0.5 * ff / Dw
|
||||
1610 IF punch_time < dt_max THEN dt = punch_time ELSE dt = dt_max
|
||||
1620 punch_fr = 1
|
||||
1630 IF dt < punch_time THEN punch_fr = ceil(punch_time / dt)
|
||||
1640 dt = punch_time / punch_fr
|
||||
1650 shifts = ceil(exp_time(tracer) / dt)
|
||||
1700 punch nl$ + nl$ + x$ + ' mixing factors...'
|
||||
1710 r1 = r_int
|
||||
1720 IF nfilt1 > 0 THEN r1 = r_int - thickn_filter1
|
||||
1730 A = height * 2 * pi
|
||||
1740 FOR i = 0 TO nfilt1 + nclay + nfilt2
|
||||
1750 IF i = 0 OR i = nfilt1 + nclay + nfilt2 THEN fbc = 2 ELSE fbc = 1
|
||||
1760 IF i > nfilt1 OR nfilt1 = 0 THEN GOTO 1810
|
||||
1770 IF i < nfilt1 THEN mixf = Dw * fbc / (xf1 / gf1) * dt * A * r1 / 1e-3
|
||||
1780 IF i = nfilt1 THEN mixf = 2 * Dw / (xf1 / gf1 + x / g) * dt * A * r1 / 1e-3
|
||||
1790 IF i < nfilt1 THEN r1 = r1 + xf1 ELSE r1 = r1 + x
|
||||
1800 GOTO 1880
|
||||
1810 IF i > nfilt1 + nclay THEN GOTO 1860
|
||||
1820 mixf = Dw * fbc / (x / g) * dt * A * r1 / 1e-3
|
||||
1830 IF i = nfilt1 + nclay AND nfilt2 > 0 THEN mixf = 2 * Dw / (xf2 / gf2 + x / g) * dt * A * r1 / 1e-3
|
||||
1840 IF i < nfilt1 + nclay THEN r1 = r1 + x ELSE r1 = r1 + xf2
|
||||
1850 GOTO 1880
|
||||
1860 mixf = Dw * fbc / (xf2 / gf2) * dt * A * r1 / 1e-3
|
||||
1870 r1 = r1 + xf2
|
||||
1880 punch nl$ + 'MIX ' + TRIM(STR$(i + 3)) + sc$ + STR$(i + 4) + STR$(mixf)
|
||||
1890 NEXT i
|
||||
1900 punch nl$ + 'END'
|
||||
2000 punch nl$ + 'TRANSPORT'
|
||||
2010 stag = 2 + nfilt1 + nclay + nfilt2
|
||||
2020 punch nl$ + ' -warnings true'
|
||||
2030 punch nl$ + ' -shifts ' + TRIM(STR$(shifts))
|
||||
2040 punch nl$ + ' -flow diff' + sc$ + ' -cells 1' + sc$ + ' -bcon 1 2' + sc$ + ' -stag ' + TRIM(STR$(stag))
|
||||
2050 punch nl$ + ' -time ' + STR$(dt)
|
||||
2060 punch nl$ + ' -multi_D true ' + STR$(Dw) + STR$(por_clay) + ' 0.0 ' + TRIM(STR$(-tort_n))
|
||||
2070 punch nl$ + ' -interlayer_D ' + interlayer_D$ + ' 0.001 0.0 ' + TRIM(STR$(G_IL))
|
||||
2080 punch nl$ + ' -punch_fr ' + TRIM(STR$(punch_fr)) + sc$ + ' -punch_c ' + TRIM(STR$(2 + stag))
|
||||
2180 FOR i = 0 to 1
|
||||
2190 punch nl$ + 'USER_GRAPH ' + TRIM(STR$(tracer + i)) + ' Example 21' + nl$
|
||||
2200 punch nl$ + ' -chart_title " ' + tracer$(tracer + i) + ' Diffusion to Outer Cell"'
|
||||
2210 punch nl$ + ' -plot_tsv_file ex21_' + tracer$(tracer + i) + '_rad.tsv'
|
||||
2220 punch nl$ + ' -axis_scale x_axis 0 ' + TRIM(STR$(exp_time(tracer + i) / (3600 * 24)))
|
||||
2230 punch nl$ + ' -axis_titles "TIME, IN DAYS" "FLUX, IN MOL/M2/S" "ACCUMULATED MASS, IN MOL"'
|
||||
2240 punch nl$ + ' -plot_concentration_vs time'
|
||||
2250 punch nl$ + ' 10 days = total_time / (3600 * 24)'
|
||||
2260 punch nl$ + ' 20 a = equi("A_' + tracer$(tracer + i) + '")'
|
||||
2270 punch nl$ + ' 30 IF get(1) = 0 AND total_time > 0 THEN put(total_time, 1)'
|
||||
2280 punch nl$ + ' 40 dt = get(1)'
|
||||
2290 A = 2 * pi * r_ext * height
|
||||
2300 i$ = TRIM(STR$(2 + i))
|
||||
2310 punch nl$ + ' 50 plot_xy days - dt / (2 * 3600 * 24), (a - get(' + i$ + ')) / dt /' + STR$(A) + ', color = Green, symbol = None'
|
||||
2320 punch nl$ + ' 60 put(a, ' + i$ + ')'
|
||||
2330 punch nl$ + ' 70 plot_xy days, equi("A_' + tracer$(tracer + i) + '"), y_axis = 2, color = Red, symbol = None'
|
||||
2340 IF tracer < 3 THEN GOTO 2360
|
||||
2350 NEXT i
|
||||
2360 punch nl$ + 'END'
|
||||
2400 IF profile$ = 'true' THEN GOSUB 3000
|
||||
2410 IF tracer < 3 THEN END # finished for Hto and Cl
|
||||
2420 IF profile$ = 'false' THEN punch nl$ + 'USER_GRAPH ' + TRIM(STR$(tracer)) + sc$ + ' -detach' ELSE punch nl$ + 'USER_GRAPH ' + TRIM(STR$(tracer + 4)) + sc$ + ' -detach'
|
||||
2440 tracer = tracer + 1
|
||||
2450 punch nl$ + 'TRANSPORT'
|
||||
2460 shifts = ceil((exp_time(tracer) - exp_time(tracer - 1))/ dt)
|
||||
2480 punch nl$ + ' -shifts ' + TRIM(STR$(shifts))
|
||||
2490 punch nl$ + ' -punch_fr ' + TRIM(STR$(punch_fr)) + sc$ + ' -punch_c ' + TRIM(STR$(2 + stag))
|
||||
2500 punch nl$ + 'END'
|
||||
2510 IF profile$ = 'true' THEN GOSUB 3000
|
||||
2520 END # finished...
|
||||
3000 punch nl$ + 'TRANSPORT'
|
||||
3010 punch nl$ + ' -shifts 0'
|
||||
3020 punch nl$ + ' -punch_fr 2' + sc$ + ' -punch_c 3-' + TRIM(STR$(2 + stag))
|
||||
3030 punch nl$ + 'USER_GRAPH ' + TRIM(STR$(tracer)) + sc$ + ' -detach'
|
||||
3040 punch nl$ + 'USER_GRAPH ' + TRIM(STR$(tracer + 4)) + ' Example 21' + nl$
|
||||
3050 punch nl$ + ' -chart_title "' + tracer$(tracer) + ' Concentration Profile: Filter1 | Clay | Filter2"'
|
||||
3060 REM punch nl$ + ' -plot_tsv_file + tracer$(tracer) + '_prof.tsv'
|
||||
3070 punch nl$ + ' -axis_scale x_axis 0 ' + TRIM(STR$((thickn_filter1 + thickn_clay + thickn_filter2) * 1e3))
|
||||
3080 punch nl$ + ' -axis_scale y_axis ' + profile_y1$(tracer)
|
||||
3090 punch nl$ + ' -axis_scale sy_axis ' + profile_y2$(tracer)
|
||||
3100 punch nl$ + ' -axis_titles ' + '"DISTANCE, IN MILLIMETERS" "FREE PORE-WATER MOLALITY" "TOTAL MOLALITY"'
|
||||
3110 punch nl$ + ' -headings ' + tracer$(tracer) + '_free ' + tracer$(tracer) + '_tot'
|
||||
3120 punch nl$ + ' -plot_concentration_vs x'
|
||||
3130 punch nl$ + ' -initial_solutions true'
|
||||
3140 punch nl$ + ' 10 IF cell_no = 3 THEN xval = 0 ELSE xval = get(14)'
|
||||
3150 punch nl$ + ' 20 IF (' + TRIM(STR$(nfilt1)) + ' = 0 OR cell_no > ' + TRIM(STR$(nfilt1 + 3)) + ') THEN GOTO 60'
|
||||
3160 punch nl$ + ' 30 IF (cell_no = 4) THEN xval = xval + 0.5 * ' + TRIM(STR$(xf1))
|
||||
3170 punch nl$ + ' 40 IF (cell_no > 4 AND cell_no < ' + TRIM(STR$(nfilt1 + 4)) + ') THEN xval = xval + ' + TRIM(STR$(xf1))
|
||||
3180 punch nl$ + ' 50 GOTO 200'
|
||||
3190 punch nl$ + ' 60 IF (cell_no = ' + TRIM(STR$(4 + nfilt1)) + ') THEN xval = xval + 0.5 * ' + TRIM(STR$(xf1)) + ' + 0.5 * ' + TRIM(STR$(x))
|
||||
3200 punch nl$ + ' 70 IF (cell_no > ' + TRIM(STR$(4 + nfilt1)) + ' AND cell_no < ' + TRIM(STR$(4 + nfilt1 + nclay)) + ') THEN xval = xval + ' + TRIM(STR$(x)) + ' ELSE GOTO 90'
|
||||
3210 punch nl$ + ' 80 GOTO 200'
|
||||
3220 punch nl$ + ' 90 IF (cell_no = ' + TRIM(STR$(4 + nfilt1 + nclay)) + ') THEN xval = xval + 0.5 * ' + TRIM(STR$(x)) + ' + 0.5 * ' + TRIM(STR$(xf2))
|
||||
3230 punch nl$ + ' 100 IF (cell_no > ' + TRIM(STR$(4 + nfilt1 + nclay)) + ' AND cell_no <= ' + TRIM(STR$(3 + nfilt1 + nclay + nfilt2)) + ') THEN xval = xval + ' + TRIM(STR$(xf2))
|
||||
3240 punch nl$ + ' 110 IF (cell_no = ' + TRIM(STR$(4 + nfilt1 + nclay + nfilt2)) + ') THEN xval = xval + 0.5 * ' + TRIM(STR$(xf2))
|
||||
3250 punch nl$ + ' 200 y1 = TOT("' + tracer$(tracer) + '")'
|
||||
3260 punch nl$ + ' 210 plot_xy xval * 1e3, y1, color = Blue, symbol = Plus'
|
||||
3270 punch nl$ + ' 220 IF cell_no = 3 THEN put(y1, 15)'
|
||||
3280 punch nl$ + ' 230 IF (cell_no < ' + TRIM(STR$(4 + nfilt1)) + ' OR cell_no > ' + TRIM(STR$(3 + nfilt1 + nclay)) + ') THEN GOTO 400'
|
||||
3290 punch nl$ + ' 240 y2 = SYS("' + tracer$(tracer) + '") / (tot("water") + edl("water"))'
|
||||
3310 punch nl$ + ' 250 REM y2 = y2 / ' + TRIM(STR$(rho_b_eps)) + x$ + ' conc / kg solid'
|
||||
3320 punch nl$ + ' 260 plot_xy xval * 1e3, y2, symbol = Circle, y_axis = 2'
|
||||
3330 punch nl$ + ' 270 IF (cell_no > ' + TRIM(STR$(5 + nfilt1)) + ') THEN GOTO 400'
|
||||
3340 punch nl$ + ' 280 IF ' + TRIM(STR$(nfilt1)) + ' THEN plot_xy ' + TRIM(STR$(thickn_filter1 * 1e3)) + ', get(15), color = Black, symbol = None'
|
||||
3350 punch nl$ + ' 290 IF ' + TRIM(STR$(nfilt2)) + ' THEN plot_xy ' + TRIM(STR$((thickn_filter1 + thickn_clay) * 1e3)) + ', get(15), color = Black, symbol = None'
|
||||
3360 punch nl$ + ' 300 put(0, 15)'
|
||||
3370 punch nl$ + ' 400 put(xval, 14)'
|
||||
3380 punch nl$ + 'END'
|
||||
3390 RETURN
|
||||
END
|
||||
-------------------------------------------
|
||||
Beginning of initial solution calculations.
|
||||
-------------------------------------------
|
||||
|
||||
Initial solution 3. tracer solution
|
||||
|
||||
WARNING: USER_PUNCH: Headings count doesn't match number of calls to PUNCH.
|
||||
|
||||
-----------------------------Solution composition------------------------------
|
||||
|
||||
Elements Molality Moles
|
||||
|
||||
Alkalinity 4.760e-04 9.520e-05
|
||||
Ca 2.580e-02 5.160e-03
|
||||
Cl 3.000e-01 6.000e-02
|
||||
Hto 1.140e-09 2.280e-10
|
||||
K 1.610e-03 3.220e-04
|
||||
Mg 1.690e-02 3.380e-03
|
||||
Na 2.400e-01 4.800e-02
|
||||
S(6) 1.410e-02 2.820e-03
|
||||
Sr 5.050e-04 1.010e-04
|
||||
|
||||
----------------------------Description of solution----------------------------
|
||||
|
||||
pH = 7.600
|
||||
pe = 13.120 Equilibrium with O2(g)
|
||||
Specific Conductance (uS/cm, 23 oC) = 29748
|
||||
Density (g/cm3) = 1.01167
|
||||
Activity of water = 0.990
|
||||
Ionic strength = 3.658e-01
|
||||
Mass of water (kg) = 2.000e-01
|
||||
Total carbon (mol/kg) = 4.786e-04
|
||||
Total CO2 (mol/kg) = 4.786e-04
|
||||
Temperature (deg C) = 23.00
|
||||
Electrical balance (eq) = -1.312e-04
|
||||
Percent error, 100*(Cat-|An|)/(Cat+|An|) = -0.10
|
||||
Iterations = 6
|
||||
Total H = 2.220258e+01
|
||||
Total O = 1.111286e+01
|
||||
|
||||
----------------------------Distribution of species----------------------------
|
||||
|
||||
Log Log Log mole V
|
||||
Species Molality Activity Molality Activity Gamma cm3/mol
|
||||
|
||||
OH- 5.132e-07 3.380e-07 -6.290 -6.471 -0.181 -3.05
|
||||
H+ 3.238e-08 2.512e-08 -7.490 -7.600 -0.110 0.00
|
||||
H2O 5.551e+01 9.899e-01 1.744 -0.004 0.000 18.06
|
||||
C(4) 4.786e-04
|
||||
HCO3- 3.609e-04 2.565e-04 -3.443 -3.591 -0.148 25.80
|
||||
NaHCO3 4.505e-05 4.901e-05 -4.346 -4.310 0.037 28.81
|
||||
CaHCO3+ 2.941e-05 2.090e-05 -4.532 -4.680 -0.148 10.13
|
||||
MgHCO3+ 1.877e-05 1.373e-05 -4.726 -4.862 -0.136 5.95
|
||||
CO2 1.380e-05 1.502e-05 -4.860 -4.823 0.037 30.80
|
||||
CaCO3 4.500e-06 4.895e-06 -5.347 -5.310 0.037 -14.61
|
||||
MgCO3 1.804e-06 1.962e-06 -5.744 -5.707 0.037 -17.09
|
||||
NaCO3- 1.804e-06 1.319e-06 -5.744 -5.880 -0.136 -0.58
|
||||
CO3-2 1.802e-06 4.594e-07 -5.744 -6.338 -0.594 -1.87
|
||||
SrHCO3+ 6.601e-07 4.691e-07 -6.180 -6.329 -0.148 (0)
|
||||
SrCO3 3.254e-08 3.540e-08 -7.488 -7.451 0.037 -14.14
|
||||
Ca 2.580e-02
|
||||
Ca+2 2.376e-02 6.593e-03 -1.624 -2.181 -0.557 -16.69
|
||||
CaSO4 2.002e-03 2.178e-03 -2.699 -2.662 0.037 7.42
|
||||
CaHCO3+ 2.941e-05 2.090e-05 -4.532 -4.680 -0.148 10.13
|
||||
CaCO3 4.500e-06 4.895e-06 -5.347 -5.310 0.037 -14.61
|
||||
CaOH+ 5.895e-08 4.312e-08 -7.230 -7.365 -0.136 (0)
|
||||
CaHSO4+ 4.779e-10 3.496e-10 -9.321 -9.456 -0.136 (0)
|
||||
Cl 3.000e-01
|
||||
Cl- 3.000e-01 2.001e-01 -0.523 -0.699 -0.176 18.52
|
||||
H(0) 0.000e+00
|
||||
H2 0.000e+00 0.000e+00 -44.617 -44.580 0.037 28.61
|
||||
Hto 1.140e-09
|
||||
Hto 1.140e-09 1.140e-09 -8.943 -8.943 0.000 (0)
|
||||
K 1.610e-03
|
||||
K+ 1.591e-03 1.061e-03 -2.798 -2.974 -0.176 9.44
|
||||
KSO4- 1.856e-05 1.358e-05 -4.731 -4.867 -0.136 (0)
|
||||
KOH 1.333e-10 1.450e-10 -9.875 -9.839 0.037 (0)
|
||||
Mg 1.690e-02
|
||||
Mg+2 1.510e-02 4.615e-03 -1.821 -2.336 -0.515 -20.63
|
||||
MgSO4 1.781e-03 1.937e-03 -2.749 -2.713 0.037 5.76
|
||||
MgHCO3+ 1.877e-05 1.373e-05 -4.726 -4.862 -0.136 5.95
|
||||
MgCO3 1.804e-06 1.962e-06 -5.744 -5.707 0.037 -17.09
|
||||
MgOH+ 7.526e-07 5.506e-07 -6.123 -6.259 -0.136 (0)
|
||||
Na 2.400e-01
|
||||
Na+ 2.378e-01 1.707e-01 -0.624 -0.768 -0.144 -0.85
|
||||
NaSO4- 2.178e-03 1.593e-03 -2.662 -2.798 -0.136 29.59
|
||||
NaHCO3 4.505e-05 4.901e-05 -4.346 -4.310 0.037 28.81
|
||||
NaCO3- 1.804e-06 1.319e-06 -5.744 -5.880 -0.136 -0.58
|
||||
NaOH 5.303e-18 5.769e-18 -17.275 -17.239 0.037 (0)
|
||||
O(0) 2.437e-04
|
||||
O2 1.218e-04 1.325e-04 -3.914 -3.878 0.037 30.24
|
||||
S(6) 1.410e-02
|
||||
SO4-2 8.079e-03 1.886e-03 -2.093 -2.724 -0.632 13.42
|
||||
NaSO4- 2.178e-03 1.593e-03 -2.662 -2.798 -0.136 29.59
|
||||
CaSO4 2.002e-03 2.178e-03 -2.699 -2.662 0.037 7.42
|
||||
MgSO4 1.781e-03 1.937e-03 -2.749 -2.713 0.037 5.76
|
||||
SrSO4 4.227e-05 4.598e-05 -4.374 -4.337 0.037 24.16
|
||||
KSO4- 1.856e-05 1.358e-05 -4.731 -4.867 -0.136 (0)
|
||||
HSO4- 6.030e-09 4.411e-09 -8.220 -8.355 -0.136 35.56
|
||||
CaHSO4+ 4.779e-10 3.496e-10 -9.321 -9.456 -0.136 (0)
|
||||
Sr 5.050e-04
|
||||
Sr+2 4.620e-04 1.280e-04 -3.335 -3.893 -0.557 -17.67
|
||||
SrSO4 4.227e-05 4.598e-05 -4.374 -4.337 0.037 24.16
|
||||
SrHCO3+ 6.601e-07 4.691e-07 -6.180 -6.329 -0.148 (0)
|
||||
SrCO3 3.254e-08 3.540e-08 -7.488 -7.451 0.037 -14.14
|
||||
SrOH+ 3.692e-10 2.588e-10 -9.433 -9.587 -0.154 (0)
|
||||
|
||||
------------------------------Saturation indices-------------------------------
|
||||
|
||||
Phase SI log IAP log K(296 K, 1 atm)
|
||||
|
||||
Anhydrite -0.69 -4.91 -4.22 CaSO4
|
||||
Aragonite -0.20 -8.52 -8.32 CaCO3
|
||||
Calcite -0.05 -8.52 -8.47 CaCO3
|
||||
Celestite 0.03 -6.62 -6.65 SrSO4
|
||||
CO2(g) -3.39 -4.82 -1.44 CO2
|
||||
Dolomite -0.15 -17.19 -17.04 CaMg(CO3)2
|
||||
Gypsum -0.33 -4.91 -4.58 CaSO4:2H2O
|
||||
H2(g) -41.48 -44.58 -3.10 H2
|
||||
H2O(g) -1.56 -0.00 1.55 H2O
|
||||
Halite -3.04 -1.47 1.58 NaCl
|
||||
O2(g) -1.00 -3.88 -2.88 O2 Pressure 0.1 atm, phi 1.000.
|
||||
Strontianite -0.96 -10.23 -9.27 SrCO3
|
||||
|
||||
|
||||
------------------
|
||||
End of simulation.
|
||||
------------------
|
||||
|
||||
------------------------------------
|
||||
Reading input data for simulation 3.
|
||||
------------------------------------
|
||||
|
||||
PRINT
|
||||
selected_output false
|
||||
PRINT
|
||||
reset false
|
||||
user_print true
|
||||
SOLUTION 4
|
||||
water 1.3963e-03
|
||||
pH 7.6
|
||||
pe 14 O2(g) -1.0
|
||||
temp 23
|
||||
Na 240
|
||||
K 1.61
|
||||
Mg 16.9
|
||||
Ca 25.8
|
||||
Sr 0.505
|
||||
Cl 300
|
||||
S(6) 14.1
|
||||
Fe(2) 0.0
|
||||
Alkalinity 0.476
|
||||
SOLUTION 5
|
||||
water 7.7322e-05
|
||||
pH 7.6
|
||||
pe 14 O2(g) -1.0
|
||||
temp 23
|
||||
Na 240
|
||||
K 1.61
|
||||
Mg 16.9
|
||||
Ca 25.8
|
||||
Sr 0.505
|
||||
Cl 300
|
||||
S(6) 14.1
|
||||
Fe(2) 0.0
|
||||
Alkalinity 0.476
|
||||
SURFACE 5
|
||||
equilibrate 5
|
||||
Su_ 5.0965e-04 5.2840e+05 6.6087e-04
|
||||
Su_ii 7.4371e-06
|
||||
Su_fes 6.9841e-07
|
||||
donnan 1.6711e-09
|
||||
EXCHANGE 5
|
||||
equilibrate 5
|
||||
X 6.2290e-04
|
||||
SOLUTION 6
|
||||
water 9.5113e-05
|
||||
pH 7.6
|
||||
pe 14 O2(g) -1.0
|
||||
temp 23
|
||||
Na 240
|
||||
K 1.61
|
||||
Mg 16.9
|
||||
Ca 25.8
|
||||
Sr 0.505
|
||||
Cl 300
|
||||
S(6) 14.1
|
||||
Fe(2) 0.0
|
||||
Alkalinity 0.476
|
||||
SURFACE 6
|
||||
equilibrate 6
|
||||
Su_ 6.2692e-04 5.2840e+05 8.1293e-04
|
||||
Su_ii 9.1484e-06
|
||||
Su_fes 8.5911e-07
|
||||
donnan 1.6711e-09
|
||||
EXCHANGE 6
|
||||
equilibrate 6
|
||||
X 7.6624e-04
|
||||
SOLUTION 7
|
||||
water 1.1291e-04
|
||||
pH 7.6
|
||||
pe 14 O2(g) -1.0
|
||||
temp 23
|
||||
Na 240
|
||||
K 1.61
|
||||
Mg 16.9
|
||||
Ca 25.8
|
||||
Sr 0.505
|
||||
Cl 300
|
||||
S(6) 14.1
|
||||
Fe(2) 0.0
|
||||
Alkalinity 0.476
|
||||
SURFACE 7
|
||||
equilibrate 7
|
||||
Su_ 7.4419e-04 5.2840e+05 9.6500e-04
|
||||
Su_ii 1.0860e-05
|
||||
Su_fes 1.0198e-06
|
||||
donnan 1.6711e-09
|
||||
EXCHANGE 7
|
||||
equilibrate 7
|
||||
X 9.0957e-04
|
||||
SOLUTION 8
|
||||
water 1.3070e-04
|
||||
pH 7.6
|
||||
pe 14 O2(g) -1.0
|
||||
temp 23
|
||||
Na 240
|
||||
K 1.61
|
||||
Mg 16.9
|
||||
Ca 25.8
|
||||
Sr 0.505
|
||||
Cl 300
|
||||
S(6) 14.1
|
||||
Fe(2) 0.0
|
||||
Alkalinity 0.476
|
||||
SURFACE 8
|
||||
equilibrate 8
|
||||
Su_ 8.6146e-04 5.2840e+05 1.1171e-03
|
||||
Su_ii 1.2571e-05
|
||||
Su_fes 1.1805e-06
|
||||
donnan 1.6711e-09
|
||||
EXCHANGE 8
|
||||
equilibrate 8
|
||||
X 1.0529e-03
|
||||
SOLUTION 9
|
||||
water 1.4849e-04
|
||||
pH 7.6
|
||||
pe 14 O2(g) -1.0
|
||||
temp 23
|
||||
Na 240
|
||||
K 1.61
|
||||
Mg 16.9
|
||||
Ca 25.8
|
||||
Sr 0.505
|
||||
Cl 300
|
||||
S(6) 14.1
|
||||
Fe(2) 0.0
|
||||
Alkalinity 0.476
|
||||
SURFACE 9
|
||||
equilibrate 9
|
||||
Su_ 9.7873e-04 5.2840e+05 1.2691e-03
|
||||
Su_ii 1.4282e-05
|
||||
Su_fes 1.3412e-06
|
||||
donnan 1.6711e-09
|
||||
EXCHANGE 9
|
||||
equilibrate 9
|
||||
X 1.1962e-03
|
||||
SOLUTION 10
|
||||
water 1.6628e-04
|
||||
pH 7.6
|
||||
pe 14 O2(g) -1.0
|
||||
temp 23
|
||||
Na 240
|
||||
K 1.61
|
||||
Mg 16.9
|
||||
Ca 25.8
|
||||
Sr 0.505
|
||||
Cl 300
|
||||
S(6) 14.1
|
||||
Fe(2) 0.0
|
||||
Alkalinity 0.476
|
||||
SURFACE 10
|
||||
equilibrate 10
|
||||
Su_ 1.0960e-03 5.2840e+05 1.4212e-03
|
||||
Su_ii 1.5994e-05
|
||||
Su_fes 1.5019e-06
|
||||
donnan 1.6711e-09
|
||||
EXCHANGE 10
|
||||
equilibrate 10
|
||||
X 1.3396e-03
|
||||
SOLUTION 11
|
||||
water 1.8407e-04
|
||||
pH 7.6
|
||||
pe 14 O2(g) -1.0
|
||||
temp 23
|
||||
Na 240
|
||||
K 1.61
|
||||
Mg 16.9
|
||||
Ca 25.8
|
||||
Sr 0.505
|
||||
Cl 300
|
||||
S(6) 14.1
|
||||
Fe(2) 0.0
|
||||
Alkalinity 0.476
|
||||
SURFACE 11
|
||||
equilibrate 11
|
||||
Su_ 1.2133e-03 5.2840e+05 1.5733e-03
|
||||
Su_ii 1.7705e-05
|
||||
Su_fes 1.6626e-06
|
||||
donnan 1.6711e-09
|
||||
EXCHANGE 11
|
||||
equilibrate 11
|
||||
X 1.4829e-03
|
||||
SOLUTION 12
|
||||
water 2.0186e-04
|
||||
pH 7.6
|
||||
pe 14 O2(g) -1.0
|
||||
temp 23
|
||||
Na 240
|
||||
K 1.61
|
||||
Mg 16.9
|
||||
Ca 25.8
|
||||
Sr 0.505
|
||||
Cl 300
|
||||
S(6) 14.1
|
||||
Fe(2) 0.0
|
||||
Alkalinity 0.476
|
||||
SURFACE 12
|
||||
equilibrate 12
|
||||
Su_ 1.3305e-03 5.2840e+05 1.7253e-03
|
||||
Su_ii 1.9416e-05
|
||||
Su_fes 1.8233e-06
|
||||
donnan 1.6711e-09
|
||||
EXCHANGE 12
|
||||
equilibrate 12
|
||||
X 1.6262e-03
|
||||
SOLUTION 13
|
||||
water 2.1966e-04
|
||||
pH 7.6
|
||||
pe 14 O2(g) -1.0
|
||||
temp 23
|
||||
Na 240
|
||||
K 1.61
|
||||
Mg 16.9
|
||||
Ca 25.8
|
||||
Sr 0.505
|
||||
Cl 300
|
||||
S(6) 14.1
|
||||
Fe(2) 0.0
|
||||
Alkalinity 0.476
|
||||
SURFACE 13
|
||||
equilibrate 13
|
||||
Su_ 1.4478e-03 5.2840e+05 1.8774e-03
|
||||
Su_ii 2.1127e-05
|
||||
Su_fes 1.9840e-06
|
||||
donnan 1.6711e-09
|
||||
EXCHANGE 13
|
||||
equilibrate 13
|
||||
X 1.7696e-03
|
||||
SOLUTION 14
|
||||
water 2.3745e-04
|
||||
pH 7.6
|
||||
pe 14 O2(g) -1.0
|
||||
temp 23
|
||||
Na 240
|
||||
K 1.61
|
||||
Mg 16.9
|
||||
Ca 25.8
|
||||
Sr 0.505
|
||||
Cl 300
|
||||
S(6) 14.1
|
||||
Fe(2) 0.0
|
||||
Alkalinity 0.476
|
||||
SURFACE 14
|
||||
equilibrate 14
|
||||
Su_ 1.5651e-03 5.2840e+05 2.0295e-03
|
||||
Su_ii 2.2839e-05
|
||||
Su_fes 2.1448e-06
|
||||
donnan 1.6711e-09
|
||||
EXCHANGE 14
|
||||
equilibrate 14
|
||||
X 1.9129e-03
|
||||
SOLUTION 15
|
||||
water 2.5524e-04
|
||||
pH 7.6
|
||||
pe 14 O2(g) -1.0
|
||||
temp 23
|
||||
Na 240
|
||||
K 1.61
|
||||
Mg 16.9
|
||||
Ca 25.8
|
||||
Sr 0.505
|
||||
Cl 300
|
||||
S(6) 14.1
|
||||
Fe(2) 0.0
|
||||
Alkalinity 0.476
|
||||
SURFACE 15
|
||||
equilibrate 15
|
||||
Su_ 1.6824e-03 5.2840e+05 2.1815e-03
|
||||
Su_ii 2.4550e-05
|
||||
Su_fes 2.3055e-06
|
||||
donnan 1.6711e-09
|
||||
EXCHANGE 15
|
||||
equilibrate 15
|
||||
X 2.0562e-03
|
||||
SOLUTION 16
|
||||
water 5.0266e-03
|
||||
pH 7.6
|
||||
pe 14 O2(g) -1.0
|
||||
temp 23
|
||||
Na 240
|
||||
K 1.61
|
||||
Mg 16.9
|
||||
Ca 25.8
|
||||
Sr 0.505
|
||||
Cl 300
|
||||
S(6) 14.1
|
||||
Fe(2) 0.0
|
||||
Alkalinity 0.476
|
||||
SOLUTION 17
|
||||
water 2.0000e-01
|
||||
pH 7.6
|
||||
pe 14 O2(g) -1.0
|
||||
temp 23
|
||||
Na 240
|
||||
K 1.61
|
||||
Mg 16.9
|
||||
Ca 25.8
|
||||
Sr 0.505
|
||||
Cl 300
|
||||
S(6) 14.1
|
||||
Fe(2) 0.0
|
||||
Alkalinity 0.476
|
||||
END
|
||||
PHASES
|
||||
A_Hto
|
||||
Hto = Hto
|
||||
log_k -15
|
||||
A_Na_tr
|
||||
Na_trCl = Na_tr+ + Cl-
|
||||
log_k -14
|
||||
A_Cl_tr
|
||||
NaCl_tr = Na+ + Cl_tr-
|
||||
log_k -14
|
||||
A_Cs
|
||||
CsCl = Cs+ + Cl-
|
||||
log_k -13
|
||||
EQUILIBRIUM_PHASES 17
|
||||
A_Hto 0 0
|
||||
END
|
||||
MIX 3
|
||||
4 6.6932e-04
|
||||
MIX 4
|
||||
5 1.9640e-04
|
||||
MIX 5
|
||||
6 1.5725e-04
|
||||
MIX 6
|
||||
7 1.8971e-04
|
||||
MIX 7
|
||||
8 2.2216e-04
|
||||
MIX 8
|
||||
9 2.5461e-04
|
||||
MIX 9
|
||||
10 2.8706e-04
|
||||
MIX 10
|
||||
11 3.1951e-04
|
||||
MIX 11
|
||||
12 3.5196e-04
|
||||
MIX 12
|
||||
13 3.8441e-04
|
||||
MIX 13
|
||||
14 4.1686e-04
|
||||
MIX 14
|
||||
15 4.4931e-04
|
||||
MIX 15
|
||||
16 7.7653e-04
|
||||
MIX 16
|
||||
17 4.2533e-03
|
||||
END
|
||||
TRANSPORT
|
||||
warnings true
|
||||
shifts 1120
|
||||
flow_direction diff
|
||||
cells 1
|
||||
bcond 1 2
|
||||
stagnant 15
|
||||
timest 1.5429e+03
|
||||
multi_d true 2.5000e-09 1.5900e-01 0.0 9.9000e-01
|
||||
interlayer_d false 0.001 0.0 700
|
||||
punch_frequency 14
|
||||
punch_cells 17
|
||||
USER_GRAPH 1 Example 21
|
||||
WARNING: No cell-lengths were read; length = 1 m assumed.
|
||||
WARNING: No dispersivities were read; disp = 0 assumed.
|
||||
-chart_title " Hto Diffusion to Outer Cell"
|
||||
-plot_tsv_file ex21_Hto_rad.tsv
|
||||
-axis_scale x_axis 0 20
|
||||
-axis_titles "TIME, IN DAYS" "FLUX, IN MOL/M2/S" "ACCUMULATED MASS, IN MOL"
|
||||
-plot_concentration_vs time
|
||||
10 days = total_time / (3600 * 24)
|
||||
20 a = equi("A_Hto")
|
||||
30 IF get(1) = 0 AND total_time > 0 THEN put(total_time, 1)
|
||||
40 dt = get(1)
|
||||
50 plot_xy days - dt / (2 * 3600 * 24), (a - get(2)) / dt / 8.2988e-03, color = Green, symbol = None
|
||||
60 put(a, 2)
|
||||
70 plot_xy days, equi("A_Hto"), y_axis = 2, color = Red, symbol = None
|
||||
END
|
||||
WARNING: Maximum iterations exceeded, 100
|
||||
|
||||
WARNING: Numerical method failed with this set of convergence parameters.
|
||||
|
||||
WARNING: Trying smaller step size, pe step size 10, 5 ...
|
||||
|
||||
WARNING: Maximum iterations exceeded, 200
|
||||
|
||||
WARNING: Numerical method failed with this set of convergence parameters.
|
||||
|
||||
WARNING: Trying reduced tolerance 1e-16 ...
|
||||
|
||||
TRANSPORT
|
||||
shifts 0
|
||||
punch_frequency 2
|
||||
punch_cells 3-17
|
||||
USER_GRAPH 1
|
||||
-detach
|
||||
USER_GRAPH 5 Example 21
|
||||
-chart_title "Hto Concentration Profile: Filter1 | Clay | Filter2"
|
||||
-axis_scale x_axis 0 2.2220e+01
|
||||
-axis_scale y_axis 0 1.2e-9
|
||||
-axis_scale sy_axis 0 1.2e-9
|
||||
-axis_titles "DISTANCE, IN MILLIMETERS" "FREE PORE-WATER MOLALITY" "TOTAL MOLALITY"
|
||||
-headings Hto_free Hto_tot
|
||||
-plot_concentration_vs x
|
||||
-initial_solutions true
|
||||
10 IF cell_no = 3 THEN xval = 0 ELSE xval = get(14)
|
||||
20 IF (1 = 0 OR cell_no > 4) THEN GOTO 60
|
||||
30 IF (cell_no = 4) THEN xval = xval + 0.5 * 1.8000e-03
|
||||
40 IF (cell_no > 4 AND cell_no < 5) THEN xval = xval + 1.8000e-03
|
||||
50 GOTO 200
|
||||
60 IF (cell_no = 5) THEN xval = xval + 0.5 * 1.8000e-03 + 0.5 * 1.7109e-03
|
||||
70 IF (cell_no > 5 AND cell_no < 16) THEN xval = xval + 1.7109e-03 ELSE GOTO 90
|
||||
80 GOTO 200
|
||||
90 IF (cell_no = 16) THEN xval = xval + 0.5 * 1.7109e-03 + 0.5 * 1.6000e-03
|
||||
100 IF (cell_no > 16 AND cell_no <= 16) THEN xval = xval + 1.6000e-03
|
||||
110 IF (cell_no = 17) THEN xval = xval + 0.5 * 1.6000e-03
|
||||
200 y1 = TOT("Hto")
|
||||
210 plot_xy xval * 1e3, y1, color = Blue, symbol = Plus
|
||||
220 IF cell_no = 3 THEN put(y1, 15)
|
||||
230 IF (cell_no < 5 OR cell_no > 15) THEN GOTO 400
|
||||
240 y2 = SYS("Hto") / (tot("water") + edl("water"))
|
||||
250 REM y2 = y2 / 1.4281e+01# conc / kg solid
|
||||
260 plot_xy xval * 1e3, y2, symbol = Circle, y_axis = 2
|
||||
270 IF (cell_no > 6) THEN GOTO 400
|
||||
280 IF 1 THEN plot_xy 1.8000e+00, get(15), color = Black, symbol = None
|
||||
290 IF 1 THEN plot_xy 2.0620e+01, get(15), color = Black, symbol = None
|
||||
300 put(0, 15)
|
||||
400 put(xval, 14)
|
||||
END
|
||||
END
|
||||
-------------------------------
|
||||
End of Run after 32.23 Seconds.
|
||||
-------------------------------
|
||||
|
||||
Loading…
x
Reference in New Issue
Block a user