5
$\begingroup$

I am fitting data on an experiment of potentiometric titration that uses alkali to titrate weak acid in analytical chemistry.

Here is the data:

data = {{0.00, 3.14}, {3.0, 4.18}, {6.0, 4.63}, {9.0, 5.08}, {10.0, 5.27}, {11.0, 5.53}, {11.5, 5.71}, {12.0, 6.00}, {12.2, 6.16}, {12.4, 6.47}, {12.6, 6.83}, {12.75, 8.97}, {13.02, 10.25}, {13.24, 10.64}, {13.40, 10.78}, {13.63, 11.01}, {13.81, 11.14}, {14.00, 11.22}, {14.50, 11.40}, {14.99, 11.52}, {16.02, 11.69}, {17.00, 11.80}, {20.00, 12.02}, {22.02, 12.11}}; ListPlot[data, PlotRange -> All]

enter image description here

The abscissae are volumes of titrant(alkali, $V_\text{OH} \backslash\text{mL}$) added into the titrated acid, and the ordinates are the $\text{pH}$ of the solution. I'd like to make a fitting curve that matches the data in order to find the inflection point(as known as the titration end-point in chemistry) of it.

The theoretical model is like as

$$ 10^{-3\text{pH}}+ \biggl( K_{\mathrm{a}}+ \frac{c_{\text{OH}}V_{\text{OH}}}{V_{\text{H}} + V_{\text{OH}}}\biggr)10^{-2\text{pH}} + \biggl( \frac{c_{\text{OH}}V_{\text{OH}}- c_{\text{H}}V_{\text{H}}}{V_{\text{H}} + V_{\text{OH}}}K_\mathrm{a} - K_\mathrm{w} \biggr)10^{-\text{pH}} - K_\mathrm{w}K_{\mathrm{a}} = 0 $$

Ka = 10^-4.76; VH = 10*10^-3;(* L *) cOH = 0.1;(* mol\[CenterDot]L^-1 *) cH = 0.1;(* mol\[CenterDot]L^-1 *) Kw = 10^-14; ContourPlot[ 10^(-3 pH) + (Ka + (cOH*VOH)/( VH + VOH))*10^(-2 pH) + (((cOH*VOH - cH*VH)/(VH + VOH))*Ka - Kw)*10^(-pH) - Kw*Ka == 0, {VOH, 0, 25*10^-3}, {pH, 0, 14}, PlotPoints -> 100]

enter image description here

At frist, I just tried Interpolation:

ifun = Interpolation[data] Plot[{ifun[x], ifun'[x], ifun''[x]}, {x, 0, 25.00}, Epilog -> Point[points]]

enter image description here

This interpolating function looks OK but since interpolation must traverse all that points, the curve is not that smooth. Hence, its first and second derivatives are very hackly.

The implicit function equation of the model can actually be solved out as the form: $\text{pH}=f(V_{\text{OH}})$ but in a quite complicated expression. In my attempt, I also tried fit data as an implicit function but all fail with either strange errors or wrong outputs. So far, I can only get the inverse function ($V_{\text{OH}} = f(\text{pH})$) fit instead of the original. With the input:

Solve[10^(-3 pH) + (Ka + (cOH*VOH)/( VH + VOH))*10^(-2 pH) + (((cOH*VOH - cH*VH)/(VH + VOH))*Ka - Kw)*10^(-pH) - Kw*Ka == 0, VOH] (* {{VOH -> -(((-1 - 10^pH Ka + 10^(2 pH) cH Ka + 10^(2 pH) Kw + 10^(3 pH) Ka Kw) VH)/((1 + 10^pH Ka) (-1 - 10^pH cOH + 10^(2 pH) Kw)))}} *) model = -(((-1 - 10^pH Ka + 10^(2 pH) cH Ka + 10^(2 pH) Kw + 10^(3 pH) Ka Kw) VH)/((1 + 10^pH Ka) (-1 - 10^pH cOH + 10^(2 pH) Kw))); data = {{3.14, 0.00}, {4.18, 3.0*10^-3}, {4.63, 6.0*10^-3}, {5.08, 9.0*10^-3}, {5.27, 10.0*10^-3}, {5.53, 11.0*10^-3}, {5.71, 11.5*10^-3}, {6.00, 12.0*10^-3}, {6.16, 12.2*10^-3}, {6.47, 12.4*10^-3}, {6.83, 12.6*10^-3}, {8.97, 12.75*10^-3}, {10.25, 13.02*10^-3}, {10.64, 13.24*10^-3}, {10.78, 13.40*10^-3}, {11.01, 13.63*10^-3}, {11.14, 13.81*10^-3}, {11.22, 14.00*10^-3}, {11.40, 14.50*10^-3}, {11.52, 14.99*10^-3}, {11.69, 16.02*10^-3}, {11.80, 17.00*10^-3}, {12.02, 20.00*10^-3}, {12.11, 22.02*10^-3}}; fit = FindFit[data, model, {Ka, {VH, 10.00*10^-3}, {cOH, 0.0782}, cH, {Kw, 1.00*10^-14}}, pH] Show[Plot[Evaluate[model /. fit], {pH, 0.00, 14.00}, PlotStyle -> {Red}], ListPlot[data, PlotRange -> All]]

deservedly I got an unexpected result:

{Ka -> -0.9725140523521845`, VH -> -0.001460202492769591`, cOH -> 0.000023960233606154922`, cH -> 0.0002099762395335789`, Kw -> 7.84017764871689`*^-18}

FindFit::cvmit: Failed to converge to the requested accuracy or precision within 100 iterations.

enter image description here

The parameters $V_\text{H} = 10.00 \backslash\text{mL}$, $c_\text{OH} = 0.0782 \backslash\text{mol·L}^{-1}$ and $K_\text{w} = 1.00 \times 10^{-14} $ are initially known from the experiment while $c_\text{H}$ and $K_\text{a}$ are unknown.

Is there any way to solve this fitting?

I'm new to use Mathematica(English also, I think). If I made any other elementary mistakes, I would be pleased to be pointed out.

the Final Solution:

Clear[Ka, Kw] data = {{0.00, 3.14}, {3.0*10^-3, 4.18}, {6.0*10^-3, 4.63}, {9.0*10^-3, 5.08}, {10.0*10^-3, 5.27}, {11.0*10^-3, 5.53}, {11.5*10^-3, 5.71}, {12.0*10^-3, 6.00}, {12.2*10^-3, 6.16}, {12.4*10^-3, 6.47}, {12.6*10^-3, 6.83}, {12.75*10^-3, 8.97}, {13.02*10^-3, 10.25}, {13.24*10^-3, 10.64}, {13.40*10^-3, 10.78}, {13.63*10^-3, 11.01}, {13.81*10^-3, 11.14}, {14.00*10^-3, 11.22}, {14.50*10^-3, 11.40}, {14.99*10^-3, 11.52}, {16.02*10^-3, 11.69}, {17.00*10^-3, 11.80}, {20.00*10^-3, 12.02}, {22.02*10^-3, 12.11}}; fitfunc[Ka_?NumericQ, cHAc_?NumericQ, VNaOH_?NumericQ] := pH /. FindRoot[ 10^(-3 pH) + (Ka + (cNaOH*VNaOH)/( VHAc + VH2O + VNaOH))*10^(-2 pH) + (((cNaOH*VNaOH - cHAc*VHAc)/( VHAc + VH2O + VNaOH))*Ka - Kw)*10^-pH - Kw*Ka == 0 /. {VHAc -> 1000.*10^-5, cNaOH -> 0.0782, Kw -> 100.*10^-16, VH2O -> 40.*10^-3}, {pH, 1.}] (* add the parameter VH2O to revise the model *) fit = FindFit[data, fitfunc[Ka, cHAc, VNaOH], {{Ka, 10^-5}, {cHAc, 0.1}}, VNaOH] (* {Ka -> 0.0000174085, cHAc -> 0.0996525} *) Show[Plot[ Evaluate[fitfunc[Ka, cHAc, VNaOH] /. fit], {VNaOH, 0, 25*10^-3}, PlotStyle -> {Darker[Green]}, AxesOrigin -> {0, 0}, Frame -> True, BaseStyle -> {FontFamily -> "CMU Serif"}], ListPlot[data, PlotRange -> All]]

enter image description here

Thanks for all your help!

References:

1: jVincent's answer

Improve this question
$\endgroup$
4
  • 1
    $\begingroup$ If you have an idea of the values that your parameters should get you may give Mathematica a hint by suggesting the starting values. You already did this with VH and Kw. This may help even if your guess is well away from the true value. Further, you wrote that several parameters are known from the experiment. Then it is a wrong thing to put them into the list of the fitting parameters. In contrast, fix them to their known values and only fit the unknown variables. The fewer are the number of the fitting parameters - the easier is the fitting procedure. This may help. $\endgroup$ Commented Apr 17, 2021 at 17:15
  • $\begingroup$ "deservedly I got an unexpected result" --why the result is unexpected and what to you want to improve? I think your fit is already quite good. $\endgroup$ Commented Apr 17, 2021 at 17:23
  • 1
    $\begingroup$ @yarchik: If I'm not mistaken, all of the parameters must be positive (they stand for things like concentrations and volumes.) $\endgroup$ Commented Apr 17, 2021 at 18:41
  • 1
    $\begingroup$ Just add MaxIterations->1000 and you get a result that fits the data. That doesn't mean that the coefficients make any sense but you do get a reasonable prediction. Of course, @MichaelSeifert 's warnings should be heeded. $\endgroup$ Commented Apr 18, 2021 at 3:19

1 Answer 1

Reset to default
12
$\begingroup$

I think you've misunderstood what the purpose of providing numbers in the third argument of FindFit is. Your code has

fit = FindFit[data, model, {Ka, {VH, 10.00*10^-3}, {cOH, 0.0782}, cH, {Kw, 1.00*10^-14}}, pH]

But Mathematica interprets this to mean that you want to find the best fit over the variations of all five parameters ($K_a, V_H, c_{OH}, c_H,$ and $K_w$). If you provide numbers in there, as you have, this indicates to Mathematica that your best guess is that (for example) $V_h \approx 10 \times 10^{-3}$. But Mathematica will still feel free to modify these parameters to get a better fit.

If you know for a certainty that $V_h = 10 \times 10^{-3}$, then you should instead set these parameters equal to their known experimental values in the model function, like so:

fit = FindFit[datatr, model /. {VH -> 10.*10^-3, cOH -> 0.0782, Kw -> 1.*10^-14}, {{Ka, 10^-5}, {cH, 0.1}}, pH] Show[Plot[Evaluate[model /. {VH -> 10.*10^-3, cOH -> 0.0782, Kw -> 1.*10^-14} /. fit], {pH, 0.00, 14.00}, PlotStyle -> {Red}], ListPlot[datatr, PlotRange -> All]] (* {Ka -> 0.0000145737, cH -> 0.107993} *)

enter image description here

Note that the given value for $K_a$ is approximately $10^{-4.84}$. Note also that I had to give Mathematica initial "guesses" for $K_a$ and $c_H$, but it then modified them to get what it considered a better fit.

This said, a word of warning: this method of using the inverse function is really not ideal because it reverses the role of the dependent and independent variable in your fit. Typically, when you do an experiment, you treat the measurements of the independent variable (the amount of OH ions added, in your case) as though they have negligible uncertainty compared to the errors on the dependent variable. The technique of least-squares fitting assumes, effectively, that the $x_i$ values are known pretty precisely while the $y_i$ variables have some random noise on them. By taking the inverse of the function and fitting that, you have effectively inverted the roles of the dependent and independent variables in the fit, and so the least-squares fit that FindFit is using may not be as reliable as it would be if you were doing things "the right way around".

Improve this answer
$\endgroup$
1
  • $\begingroup$ Thanks for your detailed explaining and warning! The reason I used the inverse function for fitting here is that the expression of its original function is much too complicated, which led to a bunch of issues. $\endgroup$ Commented Apr 18, 2021 at 8:13

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.