problem with strong and weak sites

[somonka1]

Hi to everyone!
I try to model Europium sorption isotherm on Goethite with the code below.

Code: [Select]

SURFACE_MASTER_SPECIES
    Get_s         Get_sOH
    Get_w         Get_wOH

EQUILIBRIUM_PHASES 1
    Goethite  0 10

SURFACE_SPECIES
#strong_sites
#experimental_res
H+ + Get_sOH = Get_sOH2+
    log_k     6.7
Get_sOH = Get_sOH
    log_k     0
Get_sOH = Get_sO- + H+
    log_k     -10.7
#to_calculate
Eu+3 + Get_sOH = Get_sOEu+2 + H+
    log_k     4.087
#weak_sites
#experimental_res
H+ + Get_wOH = Get_wOH2+
    log_k     6.7
Get_wOH = Get_wOH
    log_k     0
Get_wOH = Get_wO- + H+
    log_k     -10.7
#to_calculate
Eu+3 + Get_wOH = Get_wOEu+2 + H+
    log_k     -50

SOLUTION 1
    temp      25
    pH        7
    units     mol/l
    density   1
    Cl        0.1
    Na        0.1 charge
    -water    1 # kg

SOLUTION_MASTER_SPECIES
    Eu            Eu+3            0     152             152

SOLUTION_SPECIES
#from_HYDRA
Eu+3 = Eu+3
    log_k     0
Eu+3 + H2O = Eu(OH)+2 + H+
    log_k     -7.9
Eu+3 + 2H2O = Eu(OH)2+ + 2H+
    log_k     -16.38
Eu+3 + 3H2O = Eu(OH)3 + 3H+
    log_k     -25.42

END

USE solution 1
SURFACE 1
    -equilibrate with solution 1
    -sites DENSITY
    Get_s      0.01      105       0.3
    Get_w      2.31
SAVE surface 2
END

SOLUTION 2
    temp      25
    pH        7
    units     mol/l
    density   1
    Eu        9.8e-015
    Cl        0.1
    Na        0.1 charge

SELECTED_OUTPUT 1
    -file                 dump.sel
    -reset                false

USER_PUNCH
-start
  10 FOR i = 1 to 15 STEP 1
  20 READ Conc_Eu
  30 DATA 2.50E-08, 1.00E-07, 5.00E-07, 2.50E-06, 5.00E-06, 1.00E-05, 7.50E-06, 2.85E-10, 1.48E-09, 2.90E-09, 5.76E-09, 5.47E-10, 2.05E-09, 4.05E-09, 4.77E-11
  40 a$ = EOL$ +  "SOLUTION " + STR$(i+2) + EOL$
  50 a$ = a$ + " temp 25" + EOL$
  60 a$ = a$ + " pH 4.8" + EOL$
  70 a$ = a$ + " units mol/l" + EOL$
  80 a$ = a$ + " density 1" + EOL$
  90 a$ = a$ + "   Eu " + STR$(Conc_Eu)+ EOL$
  100 a$ = a$ + "   Cl 0.1 charge" + EOL$
  110 a$ = a$ + "   Na  0.1" + EOL$
  120 a$ = a$ + "    -water    1 " + CHR$(35) + " kg" + EOL$
  130 a$ = a$ + 'PRINT ' + CHR$(59) + ' -reset false' + CHR$(59) + ' -echo_input true' + CHR$(59) + ' -selected_output false' + EOL$
  140 a$ = a$ + "END" + EOL$
  150 a$ = a$ + "PHASES" + EOL$                                 # pH = constant" + EOL$
  160 a$ = a$ + "Fix_H+" + EOL$
  170 a$ = a$ + "H+ = H+" + EOL$
  180 a$ = a$ + "log_k  0.0" + EOL$
  190 a$ = a$ + "USE solution " + STR$(i+2) + EOL$
  200 a$ = a$ + "USE surface 2" + EOL$
  210 a$ = a$ + "EQUILIBRIUM_PHASES 1" + EOL$ 
  220 a$ = a$ + "Fix_H+ -4.8  NaOH 10.0" + EOL$
  230 a$ = a$ + 'PRINT ' + CHR$(59) + ' -reset false' + CHR$(59) + ' -echo_input true' + CHR$(59) + ' -selected_output true' + EOL$
  240 a$ = a$ + "END" + EOL$
  250 PUNCH a$
  260 NEXT i
-end
END

SELECTED_OUTPUT
    -file                 Eu+3_on_Get.txt
    -reset                false
    -ph                   true

USER_PUNCH 1
    -headings Eu_sorbed Get_sOEu+2_sorbed Get_wOEu+2_sorbed 
    -start
  10 punch (MOL("Get_sOEu+2") + MOL("Get_wOEu+2")) / (TOT("Eu") + MOL("Get_sOEu+2") + MOL("Get_wOEu+2")) * 100
  20 punch (MOL("Get_sOEu+2")) / (TOT("Eu") + MOL("Get_sOEu+2") + MOL("Get_wOEu+2")) * 100
  30 punch (MOL("Get_wOEu+2")) / (TOT("Eu") + MOL("Get_sOEu+2") + MOL("Get_wOEu+2")) * 100
    -end

INCLUDE$ dump.sel

END
This code must simulate the saturation of the surface, while we add more Europium, using the small concentration of strong sites.

Code: [Select]

SURFACE 1
    -equilibrate with solution 1
    -sites DENSITY
    Get_s      0.01      105       0.3
    Get_w      2.31
But with this option there is no saturation, sorption procent remains constatn and very low.
But if I exclude weak sites in such way

Code: [Select]

SURFACE 1
    -equilibrate with solution 1
    -sites DENSITY
    Get_s      0.01      105       0.3
#    Get_w      2.31
the code start to work OK.
Also, the constant for weak sites manually made very low

Code: [Select]

#to_calculate
Eu+3 + Get_wOH = Get_wOEu+2 + H+
    log_k     -50
Is there any way to fix this problem with weak sites?
Thanks a lot!

[dlparkhurst]

Could you make a simple input file that demonstrates your issue? One that   has only one or two simulations and does not include a file?

[somonka1]

The problem was solved.
The answer was that constant of sorption in the presence and in the the absence of weak sites differ noticeably, due to the different total charge of the surface. So, to model surface saturation in the presence of weak sites I correct the constant of Eu sorption.
Thank you for help!