[PoetryInMotion]

I have written a small MATLAB script that simulates a titration curve for the titration of a weak acid with NaOH (or similar) as titrator.

The result when I simulate the titration of 0.5 M acetic acid (pKa=4.76) with 0.5 M NaOH is this:



What bothers me is that I do not see any clear half-titration point. There ought to be one, according to with what seem to be experimental titration curves here and here. Hopefully someone here can help me understand why my script does not work. I have tried my best to eliminate approximation (for instance that the volume of the solution is approximately constant or that the self-dissociation of water does not affect the pH) but perhaps physical chemistry above my capacity is needed to get enough accuracy to see the half-titration point.

Below I include the code and the mathematical background, so that you get a better picture of what I have done this far.

Mathematical principle.

What my script does, is that it after each small addition of base, compute the pH in the solution. Both the protolysis equilibrium and the self-dissociation of water is taken into account. More precisely, I solve the following system numerically:

[tex]\begin{cases}K[HA]=[H^+][A^-]\\ Kw=[H^+][OH^-]\\ C_{a0} = [HA] + [A-]\\ C_{b0} = [OH^-]-[H^+] + [A-]\,,\end{cases}[/tex]

where K is the acid constant for my acid, Kw is the self-dissociation constant for water and C₀ represents the (imaginary) concentration of acid and base before the system is allowed to equilibriate.

My best guess for why

Substitutions give rise to
[tex]C_{a0}=\frac{[H^+][A^-]}{[K]}+[A^-]=(\frac{[H^]}{K}+1)(C_{b0}+[H^+]-\frac{K_w}{[H^+]})\,,[/tex]
[tex]K[H^+]C_{a0}=([H^+]+K)(C_{b0}[H^+]+[H^+]^2-K_w)\,,[/tex]
and finally
[tex][H^+]^3 + (K+C_{b0})[H^+]^2 + (C_{b0}K - C_{a0}K - K_w)[H^+] -KK_w=0\,,[/tex]
which is a polynomial which can easily be solve with a built in function in MATLAB.

My code.

Code: [Select]

Kw = 1e-14;
% Initial concentration of acid
C_acid = 0.5; % mol/dm^3
% Molarity of added base
C_base = 0.5; % mol/dm^3

% Initial volume of acid
vol = 25.0; % ml
% Stop graph when x ml base is added
max_added = 50.0; % ml

% Acid constant
K = exp(-4.76);

% Graph settings
dx = 0.0001;
x = 0.0:dx:max_added;
% pH = zeros(10 * max_added)

% Calculate the initial number of moles of acid
n_a = C_acid * vol * 0.001;

for i = 1:1/dx*max_added+1
    % Calculate the number of moles of base added
    n_b = x(i) * C_base * 0.001;
    % The total volume
    vol_tot = vol + x(i);
    % Initial concentrations of acid and base
    C_a0 = n_a / vol_tot * 1000;
    C_b0 = n_b / vol_tot * 1000;
    % Compute the concentration of hydrogen ions
    % Solve the polynomial equation
    R = roots([1, (K+C_b0), (C_b0*K - C_a0*K - 1e-14), -K*1e-14]);
    % Remove all chemically irrelevant solutions
    R = R(R==real(R) & R>0 & R<C_a0);
    % If there are more than one chemically relevant solution, display the number
    if not(length(R) == 1)
        disp(length(R))
    end
    H = min(R(R==real(R) & R>0 & R<C_a0));
    % Calculate the pH value
    pH(i) = -log10(H);
end
% Plot it
plot(x,pH,'red')
xlabel('TIllsatt volym bas (ml)')
%set(gca,'XTick',0:1:max_added);
%set(gca,'YTick',round(pH(1)):1:round(pH(length(pH)))+1);
%set(gca,{'xminorgrid' 'yminorgrid'}, {'on' 'on'});
ylabel('pH')
title('Titreing av ättiksprit')

Edit: For some reason, the LaTeX tool does not seem to work for me.

[mjc123]

Not sure how your code works, but when you write
K = exp(-4.76)
does that mean 10^-4.76 or e^-4.76?
If the latter, then your pKa, and pH at half titre, would be 2.1 - which looks close to what the graph says.

[Arkcon]

I don't understand what you claim to be missing, looking at your examples, there isn't much of a graphical structure associated with half-titration point.

[PoetryInMotion]

Not sure how your code works, but when you write
K = exp(-4.76)
does that mean 10^-4.76 or e^-4.76?
If the latter, then your pKa, and pH at half titre, would be 2.1 - which looks close to what the graph says.

Gah! Such an embarrassing mistake

Lesson of today: never underestimate your own stupidity when debugging your code.

New result:



I am incredibly thankful.

[samcarano]

I am trying to work through the derivations but I don't understand where the last two of the 4 equations come from.
Would somebody mind explaining to me where Ca0 and Cb0 are obtained from? 

The equations I see besides the ones mentioned here that seem relevant:
[HA]-->[H+]+[A-]
[HA]+[OH]-->[A-]+[H+]

I would think Ca0 the same thing as [HA] and Cb0 the same as [OH]

Clearly that is not true here but I can't seem to figure what I am not getting.

Any help would be much apprecitaed,

Thanks

[Borek]

Would somebody mind explaining to me where Ca0 and Cb0 are obtained from?

They aren't "obtained", they are defined.

They are analytical concentrations. When you put acid HA in the solution it dissociates, you have HA and A^- present - and even if you don't know the dissociation extent it is obvious that the sum of both forms must be identical to the amount of acid put in the solution. Say, you put 1 mole of acid in 1000 mL - no matter what the Ka is, [HA]+[A^-]=1.

[samcarano]

That makes some sense.  Thank you.  So you are adding the acid concentration and the acid conjugate concentration.

Cb0 still confuses me:

Why is Cb0= [OH-]-[H+]+[A-]?
I would think it would be just Cb0= [OH-]+[H+] to represent the total base in the solution.  Where does the [A-] come from and why not add the base and its conjugate again?

I might be getting throw off by what they poster means by "imaginary concentration"

Thanks for the help

[Borek]

Yeah, it can be confusing. Cb0 equation combines three things - mass balance for base, assumption that base is strong and fully dissociated, and charge balance.

BOH  B⁺ + OH^-

Start with

Cb0=[B⁺] + [BOH]

write charge balance and see what you get from these.