Skip to content

Original ResearchVolume 1, Issue 3, Spring 2025

What a Buffer Is Actually Resisting

Buffer capacity computed exactly, not approximated. Where a buffer peaks, how fast it fails, which buffers have no peak at all, and why blood buffers nineteen times better than the same chemistry in a stoppered bottle.

Written jointly by the Science Journaling Club

Study
Computational study
Review
Peer-edited by the club review board
Length
12 sectionsabout 27 min

Abstract

Buffer capacity is introduced as a formula and usually abandoned there. We computed it instead, as the exact derivative β=dCb/dpH\beta = dC_b/d\mathrm{pH} of a full charge-balance solver, across dissociation constants from pKa 2 to pKa 12 and concentrations from 0.1 mM to 1 M. No solution was mixed at any point. The club has no laboratory, so the computation is the experiment. Three independent routes to β\beta agree: the analytic expression, a five-point finite difference, and a simulated titration that adds a measured dose of base and re-solves the nonlinear equilibrium, which match to 1.4×10−10 and 1.8×10−8 respectively. The maximum sits at pH equal to pKa, but not exactly: at 0.100 M and pKa 5 it sits 3.477×10−4 pH units below, a displacement the water term predicts to a ratio of 1.0008, and which is identically zero only at pKa 7.000. Peak height matches ln⁡(10) C/4\ln(10)\,C/4 to 8 parts per million at 0.100 M. Capacity falls to half over a window of 1.531103 pH units, so at one unit from pKa only 33.06% of peak survives and at two units 3.92%. Six of twenty-two tested buffers have no capacity maximum at all. A two-component buffer splits into twin peaks above a separation of 1.143895 pH units, which is 2 arccosh(2)/ln⁡102\,\mathrm{arccosh}(2)/\ln 10. For blood at 37 °C, pH 7.40 and pCO2 40 mmHg, the model produces 24.50 mmol/L bicarbonate without being told to, and gives an open-system capacity of 57.32 against a closed-system 2.940 mmol/L per pH unit, a ratio of 19.50. A Monte Carlo over five uncertain physiological inputs, 200,000 draws, gives 20.162 ± 0.003 with a standard deviation of 1.213. The open system at 1.3 pH units off its own optimum beats the same closed chemistry sitting exactly on that optimum by a factor of 3.85. Seed 20250321.

A Formula That Stops Too Soon

Buffer capacity has a definition, and it is a good one. Van Slyke wrote it down in 1922 as the slope of a titration curve [3], building on the mass-action treatment Henderson had published fourteen years earlier [1] and the logarithmic form Hasselbalch gave it in 1917 [2]: the number of moles of strong base you must add to a litre of solution to raise its pH by one unit. Call it β\beta. The larger it is, the harder the solution fights.

For a single weak acid HA at total concentration CC with dissociation constant KaK_a, the charge balance after adding CbC_b moles per litre of sodium hydroxide reads Cb+[H+]=[OH−]+[A−]C_b + [\mathrm{H}^+] = [\mathrm{OH}^-] + [\mathrm{A}^-], which rearranges to

Cb(h)  =  Kwh  −  h  +  CKaKa+hC_b(h) \;=\; \frac{K_w}{h} \;-\; h \;+\; \frac{C K_a}{K_a + h}

with h=[H+]h = [\mathrm{H}^+]. Differentiate with respect to pH, remembering that dh/dpH=−ln⁡(10) hdh/d\mathrm{pH} = -\ln(10)\,h, and the capacity falls out whole:

β  =  ln⁡(10)[Kwh  +  h  +  CKah(Ka+h)2]\beta \;=\; \ln(10)\left[\frac{K_w}{h} \;+\; h \;+\; \frac{C K_a h}{(K_a + h)^2}\right]

That is the formula. Textbooks print it, note that the middle term peaks when h=Kah = K_a, state that the peak value is ln⁡(10) C/4\ln(10)\,C/4, advise you to work within one pH unit of pKa, and move on to the next chapter. Urbansky and Schock set out the full derivation and complained, in 2000, that this is exactly where most treatments stop [11]. Every one of those statements is approximately true. Three of them are not exactly true, one of them is wrong in a way that matters for real laboratory buffers, and the last one is a rule of thumb that turns out to be a full pH unit too generous.

The figure below is what the formula looks like when you draw it rather than reading it. Five buffers at the same concentration, five different dissociation constants. Each is a hill. The hills are identical in shape and height, and they slide along the pH axis with pKa. At the far edges the solvent takes over and the curves climb away.

Buffer capacity against pH for five 0.100 molar monoprotic buffers with pKa 3, 5, 7, 9 and 11, showing five identical peaks and the shaded half-capacity window around pKa 7 0 2 4 6 8 10 12 14 pH 0.000 0.015 0.030 0.045 0.060 0.075 beta (mol/L per pH) 3 5 7 9 11 0.100 M monoprotic buffers, 25 °C; peak ln(10)·C/4 = 0.05756
Figure 1. Capacity of five 0.100 M monoprotic buffers at 25 °C, computed from the exact charge balance. Every peak is the same height, 0.0575651 mol/L per pH unit, which is ln⁡(10) C/4\ln(10)\,C/4 plus a solvent correction of eight parts per million at pKa 7. The shaded band is the half-capacity window of the pKa 7 buffer, computed at 1.531152 pH units wide against an analytic 1.531103. The curves turn upward at both edges because water is itself a buffer once you are far enough from neutrality.

What We Actually Did

The club has computers and people. It does not have a burette, a pH meter, a fume cupboard or anywhere to put one. So the study is arithmetic, and the arithmetic is the experiment. Nothing below was observed in a beaker. When this article says "measured", it means measured from the output of a program written by club members, in the same sense that you measure the digits of π\pi by computing them.

The program does four things. It writes the full proton condition for a solution with any number of buffers in it, keeping the water terms. It inverts that condition by bisection, so that given an amount of added strong base it returns the pH to the last bit a double can hold. It differentiates the condition three separate ways so that the three answers can be compared. And it runs the carbon dioxide system twice, once stoppered and once with the lungs attached.

Nothing is approximated. There is no Henderson-Hasselbalch step anywhere in the code, no assumption that [HA]≈C[\mathrm{HA}] \approx C, no dropping of the water terms, and no linearisation of the titration curve. De Levie showed in 1993 that a titration curve has a single closed form needing no approximations and no segmentation [10], and the solver here follows that spirit. Where the literature's shortcuts appear, they appear as something to check the exact answer against, never as a step in getting it. Po and Senozan catalogued where the Henderson-Hasselbalch approximation goes wrong and by how much [12]. We did not set out to reproduce their tables, and we have not checked our numbers against them line by line, but the same kind of error appears in section 3 every time the exact solver is put beside a shortcut.

20250321master seed
4.7×10−13solver round trip
200,000Monte Carlo draws
0.3 stotal runtime

The assumptions are listed in full in the script's docstring and are worth stating here too. Activities are replaced by concentrations everywhere. Equilibrium is instantaneous, which ignores how slowly carbon dioxide hydrates, slowly enough that evolution invented carbonic anhydrase about it. Dissolved CO2 and H2CO3 are lumped together, as they always are. Temperature is fixed, at 25 °C for the generic work and 37 °C for the blood work. Adding base does not change the volume.

Three Routes to One Number

A derivative computed by a program is a claim, and claims need auditing. We computed β\beta three ways that share as little machinery as possible.

The first is the analytic expression above, evaluated directly. The second is a five-point central difference of Cb(pH)C_b(\mathrm{pH}), which uses the charge balance but never the derivative. The third adds a real dose of strong base to the modelled solution, re-solves the whole nonlinear equilibrium by bisection for the new pH, does the same with an equal dose of acid, and divides. That third route is a simulated titration. It never differentiates anything and never touches equation (2).

Buffer capacity against pH for 0.100 molar acetic acid with simulated titration points on the analytic curve, and a lower panel showing the relative difference of the two numerical routes 0.0001 0.001 0.01 0.1 1 beta (mol/L per pH) pKa water wall water wall 0.100 M acetic acid, pKa 4.756, 25 °C 0 2 4 6 8 10 12 14 pH 1e-13 1e-10 1e-7 |rel. diff| titration route finite difference
Figure 2. Validation. Upper panel: the analytic capacity of 0.100 M acetic acid as a line, with the simulated-titration result plotted as open circles every half pH unit. Lower panel: how far each numerical route sits from the analytic one, as a relative difference. The finite difference runs at 10−12 and is limited by cancellation in double precision. The titration route runs at 10−8 and is limited by the truncation error of a central difference over a finite 10−4 pH step.

Across pH 1 to 13 the three agree. The worst relative disagreement between the analytic expression and the finite difference is 1.395×10−10, which is what you get when you subtract two nearly equal double-precision numbers. The worst disagreement with the titration route is 1.769×10−8, and that one has a known cause: it is the truncation error of a central difference taken over a step of 10−4 pH units. Shrink the step and it shrinks as the square, exactly as a central difference should.

The solver itself was checked separately. Round-tripping pH to added base and back over the range 0.5 to 14.5 returns the input to 4.7×10−13. At half neutralisation the exact pH of 0.100 M acetic acid is 4.756304 against the pKa of 4.756000 that Harned and Ehlers measured in 1932 [4], a gap of 3×10−4 that is the water term and not an error. At the equivalence point the exact answer is 8.878022 against the textbook weak-base shortcut's 8.878000.

Plain Arithmetic

Three of the results in this article can be got with a pencil. Here they are, worked.

The peak height. Set h=Kah = K_a in the buffer term. The numerator becomes CKaKaC K_a K_a, the denominator (2Ka)2=4Ka2(2K_a)^2 = 4K_a^2, so the term is C/4C/4. Multiply by ln⁡10=2.302585\ln 10 = 2.302585.

2.302585 × 0.100 ÷ 4 = 0.05756463 mol/L per pH unit
computed by the program at pKa 7.000, 0.100 M = 0.05756509
difference = +4.605×10−7, which is exactly ln(10) × (10−7 + 10−7), the solvent

The half-capacity window. Write r=10pKa−pHr = 10^{\mathrm{pKa}-\mathrm{pH}}, so r=h/Kar = h/K_a. The buffer term relative to its peak is 4r/(1+r)24r/(1+r)^2. Set that to one half and you get 8r=(1+r)28r = (1+r)^2, so r2−6r+1=0r^2 - 6r + 1 = 0, so r=3+22=5.828427r = 3 + 2\sqrt{2} = 5.828427.

log10(5.828427) = 0.7655514 pH units either side
check: 4 × 5.828427 ÷ (6.828427)2 = 23.31371 ÷ 46.62742 = 0.500000
full window = 2 × 0.7655514 = 1.5311027 pH units
computed by the program at pKa 7.000 = 1.5311518, a gap of 4.9×10−5

What is left one unit out. At pH=pKa±1\mathrm{pH} = \mathrm{pKa} \pm 1, r=10r = 10 or r=0.1r = 0.1. Either way 4r/(1+r)2=40/1214r/(1+r)^2 = 40/121.

40 ÷ 121 = 0.330579, so a third of the peak, not half
two units out: 400 ÷ 10201 = 0.039212, so one part in 25.5
three units out: 0.003992, one part in 250.5

The last line is the one to keep. The standard advice is to stay within one pH unit of pKa, and that advice is quietly conceding two thirds of the buffer. The honest window, the one where you still have half the strength you paid for, is 0.766 units either side.

Notes from the Table

Tuesday, first session
Solver written. Bisection on pH rather than on h, because h spans fourteen decades and pH spans fourteen units. Round trip clean to 1e-13 straight away. Somebody asked whether we should use Newton. We tried. It is faster and it is not always right, and the whole point of this study is being right.

Wednesday
Peak-finder returns pH 8.0 for a pKa 3 buffer. Argued about this for twenty minutes. It is a minimum. Blind bisection on the derivative over the whole pH range cannot tell a maximum from a minimum, and beta has both once the water walls are in. Rewrote the finder to demand a bracket with the derivative positive on the left and negative on the right.

Wednesday, later
With the finder fixed, pKa 5 at 0.100 M gives a peak at pH 4.99965232, not 5. We assumed a bug. Hunted it for most of an afternoon. It is not a bug. Set d(beta)/d(pH) to zero and the condition is −Kw/h + h + CKa h(Ka−h)/(Ka+h)3 = 0, and at h = Ka the third term vanishes while the first two do not. The peak is exactly at pKa only when Ka2 = Kw, which is pKa 7.000. The displacement is −(8/(C ln10))(Ka − Kw/Ka), and the computed ratio of found to predicted is 1.0008 at pKa 5, 1.0081 at pKa 4, 1.0943 at pKa 3 where the expansion is giving up.

Thursday
pKa 2 returns nothing at all. No peak. Checked the bracket logic twice, then plotted it and saw why: at pKa 2 the capacity of a 0.100 M buffer is smaller than the strong-acid wall it is sitting on, so the curve falls monotonically through the buffer's own region. Six of the twenty-two grid cases behave this way. Nobody in the room had ever seen this stated anywhere.

Friday
Blood. Gave the model pH 7.40 and pCO2 40 mmHg with the clinical apparent pK1' of 6.10 and asked what bicarbonate it implied. It said 24.50 mmol/L. That is what a blood gas analyser reports. Nothing was tuned. Long pause in the room.

The Buffers With No Peak

This one surprised us. Sweep pKa from 2 to 12 at two concentrations and ask the program where each buffer's capacity peaks. Six of the twenty-two answers are that it does not.

The reason is the solvent. Near pH 2 the term ln⁡(10) h\ln(10)\,h is 0.0230 mol/L per pH unit and falling steeply. A 0.100 M buffer can contribute at most 0.0576 at the very top of its hill, and the hill is sitting on a slope steeper than itself. The sum has no turning point. What the buffer produces is a shoulder, a place where the descent briefly slows, and not a maximum.

Push the concentration down to 0.010 M and the same thing starts happening at pKa 3. Go the other way, to pKa 11 and 12, and it happens in mirror image against the hydroxide wall.

So a bottle labelled as a pH 2 buffer is not doing what the label implies. Most of its resistance to added base comes from the free acid already in solution, which is to say from the solvent chemistry rather than from the conjugate pair you paid for. That is not a criticism of the bottle. It works. It is a criticism of the mental picture. Good and colleagues made proximity of pKa to the working pH their first selection criterion when they designed the biological buffers still in use today [6], and this result says the criterion has teeth at the edges of the pH range as well as in the middle.

Two in One Beaker

Put two buffers in the same solution and their capacities add, because each conjugate pair takes up protons without asking the other one. That is equation (3) in the script, and it is the whole basis of the mixed buffers people actually use, and the reason a shelf of Good's buffers covers a continuous range rather than a set of isolated points [6].

Fix the total concentration at 0.100 M, split it evenly, and slide the two pKa values apart by a separation DD either side of pH 7. For small DD you get one peak with a flattened top. For large DD you get two peaks with a valley between. Where does it change over?

Set the second derivative of the sum at the midpoint to zero. Writing s(u)s(u) for the shape of one buffer term with u=pH−pKau = \mathrm{pH} - \mathrm{pKa}, the shape is 1/(2+2cosh⁡(uln⁡10))1/(2 + 2\cosh(u\ln 10)), and the condition reduces to cosh⁡(ln⁡(10) D/2)=2\cosh(\ln(10)\,D/2) = 2. So

D∗  =  2 arccosh(2)ln⁡10  =  1.143895 pH unitsD^{*} \;=\; \frac{2\,\mathrm{arccosh}(2)}{\ln 10} \;=\; 1.143895 \text{ pH units}

The numerical scan flips from a single central peak to a split pair between D=1.1439D = 1.1439 and D=1.2000D = 1.2000, bracketing that number. Below the threshold, the pair is one broad hill. Above it, two.

Buffer capacity against pH for a two-component buffer at four separations, showing a single flat peak below the threshold and twin peaks above it 4 5 6 7 8 9 10 pH 0.000 0.013 0.026 0.039 0.052 0.065 beta (mol/L per pH) D = 0.000 D = 1.144 D = 2.000 D = 3.000 0.100 M total, split between two buffers about pH 7
Figure 3. A 0.100 M total buffer split evenly between two components whose pKa values straddle pH 7 by ±D/2\pm D/2. At D=0D = 0 the pair is one buffer of double strength. At D=1.143895D = 1.143895, the analytic split threshold, the top is as flat as it can be without dividing. At D=2D = 2 and D=3D = 3 the valley between the two peaks deepens toward the single-component floor.

The engineering consequence is a trade. Measure the width over which the mixture still delivers half of what a single 0.100 M buffer manages at its own peak. A single buffer holds that for 1.531 pH units. The pair at the threshold holds it for 1.792. Push the separation to 3.0 and the width reaches 3.107, but the middle has sagged to 0.0068 against a single buffer's 0.0576, so what you have is two buffers with a hole between them rather than one wide buffer. Width is bought with height, and the exchange rate is not generous.

Real pairs behave as predicted. Acetate and the second phosphate step sit 2.442 pH units apart, well past the threshold, and the mixture sags to 42.1% of its peak in the middle. Citrate's own second and third dissociations sit 1.635 apart and sag only to 83.5%. Phosphate with tris, 0.874 apart and inside the threshold, shows no dip at all.

Blood Should Not Work

Here is the problem, stated as baldly as possible. The dominant buffer in blood plasma is the carbon dioxide system. Its apparent first dissociation constant at 37 °C and plasma ionic strength corresponds to pK1' = 6.10 [15]. Blood is held at pH 7.40. That is 1.30 pH units away from the optimum.

Section 4 gives the penalty for being 1.30 units off. Retained fraction 4r/(1+r)24r/(1+r)^2 with r=101.3=19.95r = 10^{1.3} = 19.95 comes to 0.1818. A buffer sitting that far from its own pKa keeps 18% of its strength. Blood's total carbonate is about 25.8 mmol/L, so the closed system's capacity at pH 7.4 should be, and in our model is, 2.940 mmol/L per pH unit. Whole blood measures 38.5 and plasma 16.1 mEq/L per pH [5], and Siggaard-Andersen built the clinical acid-base chart on a non-bicarbonate buffer value of about that size [7]. The carbonate system, on this accounting, is contributing almost nothing.

Something is missing from that accounting, and the missing thing is that you are breathing.

In a stoppered bottle, adding acid to a bicarbonate solution converts bicarbonate into dissolved carbon dioxide, which stays put. The total carbonate is conserved and the buffer eats itself. In a body, the dissolved carbon dioxide is removed by the lungs as fast as it appears, and pCO2 barely moves. The right boundary condition is not "total carbonate fixed". It is "dissolved carbon dioxide fixed". Rewriting the proton condition with S=α pCO2S = \alpha\,p_{\mathrm{CO}_2} held constant gives

βopen  =  ln⁡(10)[Kwh+h+K1Sh+4K1K2Sh2]\beta_{\text{open}} \;=\; \ln(10)\left[\frac{K_w}{h} + h + \frac{K_1 S}{h} + \frac{4K_1K_2S}{h^2}\right]

and the third term is simply ln⁡(10) [HCO3−]\ln(10)\,[\mathrm{HCO}_3^-]. Our model, given only pH 7.40, pCO2 40 mmHg and the apparent constants, produces 24.50 mmol/L of bicarbonate and 57.32 mmol/L per pH unit of capacity. Against the closed system's 2.940, that is a factor of 19.50.

Buffer capacity against pH for the carbon dioxide system, open and closed, on a logarithmic scale, with the two values at blood pH 7.4 marked and the ratio between them 6.0 6.4 6.8 7.2 7.6 8.0 8.4 pH 0.001 0.01 0.1 1 beta (mol/L per pH) open, pCO2 held closed, CT held blood pH 7.40 ×19.50 carbon dioxide system, 37 °C, pCO2 40 mmHg
Figure 4. The same chemistry, two boundary conditions. The closed curve has a maximum near its own pK1' and falls away. The open curve does not turn over anywhere above pH 4.506, because at fixed pCO2 the bicarbonate concentration rises by a factor of ten for every pH unit and carries the capacity up with it. At blood pH the gap is a factor of 19.50.

Now the part that changes how the whole thing reads. The open system has no optimum. Not a broad one, not a displaced one. None. Above pH 4.506 the capacity rises without turning back, because bicarbonate at fixed pCO2 is proportional to 1/h1/h. Being 1.30 units above pK1' is not a penalty in an open system. It is the reason there is so much bicarbonate in the first place.

Put a number on that. A closed carbonate system with the same 25.8 mmol/L of total carbonate, sitting exactly on its own optimum at pH 6.10, is worth 14.881 mmol/L per pH unit. The open system at pH 7.40, a full 1.30 units off that optimum, is worth 57.317. The off-optimum open system beats the on-optimum closed one by 3.85 times. To match it with an ordinary closed buffer perfectly tuned to pH 7.40 you would need 99.6 mmol/L of it, nearly four times the carbonate blood actually carries.

The vividness comes out in the dose table. Add 10 mmol/L of strong acid to a litre and re-solve the equilibrium exactly. The open system goes from 7.400 to 7.176. The closed one goes to 6.220. One of those is a bad afternoon and the other is a death certificate, and the chemistry in the bottle is identical. The only difference is whether anything is carrying the carbon dioxide away.

The Objection We Take Seriously

The strongest objection to this study is that the open system is not a buffer and calling it one is a category error.

The argument goes like this. A buffer resists pH change using only what is dissolved in it. The open carbon dioxide system resists pH change by continuously exporting matter to a reservoir of effectively infinite size. Of course it looks better. You have compared a closed box to a box with a hole in it and a pump attached, declared the second one nineteen times the buffer, and learned nothing about buffering. Worse, our own parameterisation makes the point for us: hold the pump still and the capacity collapses. What is being measured is the pump.

We think this objection is largely correct and that it improves the result rather than destroying it.

It is correct that βopen\beta_{\text{open}} is not a property of the solution alone. It is a property of the solution plus its boundary condition, and the boundary condition is doing most of the work. Our own section 10 numbers say so plainly: let pCO2 drift by one tenth of a decade per pH unit and you lose 10% of the capacity immediately. Let it drift by 0.956 decades, which is what a cork does, and you are back to the closed value.

Where we part company with the objection is on what follows. Van Slyke's definition is a derivative of a titration curve, and a titration curve is defined by whatever is held fixed while you titrate [3]. Nothing in the definition requires the system to be closed. If it did, the definition would be unable to describe the single most important pH regulator in vertebrate physiology, and a definition with that hole in it is the wrong definition. Roos and Boron treat open-system buffering power as a legitimate quantity for exactly this reason [8], Boron makes the same case in his teaching account of intracellular pH regulation [13], and Stewart's reformulation makes pCO2 an independent variable rather than a dependent one [9].

The honest version of the finding is therefore narrower than the headline, and better. It is not that bicarbonate is a wonderful buffer. It is that bicarbonate is a mediocre buffer at pH 7.4, and that attaching it to a controlled reservoir converts a mediocre buffer into an excellent one, by a factor this model puts at 19.50 and a Monte Carlo puts at 20.16. The chemistry contributes the 2.94. The lungs contribute the other 54.38.

One further concession. Our open model holds pCO2 perfectly constant, which no lung does. That choice sets an upper bound and we have not modelled the ventilatory control loop at all. A model with realistic respiratory dynamics would give a number between our two, closer to the open value for a slow acid load and closer to the closed value for a fast one.

Where a Different Choice Would Have Changed the Answer

Five inputs go into the blood calculation and none of them is known exactly. We varied them one at a time, and then all at once.

QuantityAnalytic or acceptedThis studyDifference
Single monoprotic buffer, 0.100 M, 25 °C
Peak capacity0.057564630.05756509+8.0×10−6 rel.
Peak position, pKa 7.000pH 7.000000pH 7.000000000 exactly
Peak position, pKa 5.000pH 5.000000pH 4.99965232−3.477×10−4
Predicted displacement, pKa 5−3.474×10−4−3.477×10−4ratio 1.00080
Half-capacity window1.53110271.5311518+4.9×10−5
Retained at pKa ± 140/121 = 0.3305790.330579exact
Retained at pKa ± 2400/10201 = 0.0392120.039212exact
Two components, 0.100 M total
Split threshold D∗D^{*}1.143895between 1.1439 and 1.2000brackets it
Half-width, single buffer1.5311031.531162+5.9×10−5
Half-width at D∗D^{*}no closed form1.792408–
Carbon dioxide system, 37 °C, pH 7.40, pCO2 40 mmHg
Bicarbonate24–26 mmol/L [16]24.50 mmol/Linside range
Dissolved CO21.2 mmol/L [15]1.228 mmol/L+2.3%
Open-system capacityln10 × 24 = 55.2657.317+3.7%
Closed-system capacity~4 mEq/L per pH [16]2.940−26.5%
Open / closed ratio1 + 101.3 = 20.95319.498−6.9%
Same, Monte Carlo mean–20.162 ± 0.003sd 1.213
Open off-optimum vs closed on-optimum–3.85×–

The sensitivity run says which choices matter. Moving pK1' from 6.10 to 6.03 takes the ratio from 19.498 to 22.461. Moving pH from 7.40 to 7.25 takes it to 14.585. Those two dominate everything else, which is what the approximation 1+10pH−pK1′1 + 10^{\mathrm{pH}-\mathrm{pK}_1'} predicts they should. Changing pCO2 from 40 to 55 mmHg moves the ratio by 0.0025, because it scales the open and closed capacities together. Switching the water constant from 37 °C back to 25 °C changes the ratio in the fifth decimal place, which is as clear a statement as you could want that water contributes nothing at blood pH.

Two modelling choices would have changed the published number. First, apparent pK2'. We used 9.80. The thermodynamic value near 10.33 gives a ratio of 20.487 instead of 19.498, a 5% shift, because the carbonate ion enters the open capacity with a weight of four. There is no consensus figure for plasma, so we put a uniform prior across the published range and carried it through. The same sensitivity turns up in ocean carbonate work, where the equivalent buffer factors depend sharply on which constants you adopt [14]. Second, the decision to derive bicarbonate from pH and pCO2 rather than to fix it at the clinical 24 mmol/L. Fixing it would have given a slightly smaller open capacity and a story about how well the model agrees with a number we had put in by hand.

Running mean of the open over closed capacity ratio against number of Monte Carlo trials, with a one standard error band narrowing as trials accumulate 10 100 1000 10000 100000 trials 19.25 19.50 19.75 20.00 20.25 20.50 20.75 21.00 open / closed ratio final mean 20.1620 ±1 standard error Monte Carlo running mean, seed 20250321
Figure 5. Monte Carlo convergence, 200,000 draws over pH, pCO2, pK1', pK2' and the CO2 solubility. The running mean settles by about 10,000 trials and the band is ± one standard error computed from the draws themselves. Final value 20.1620 ± 0.0027. The band is narrow because we took many draws. It says nothing about how well plasma chemistry is known, which is the standard deviation of 1.213, four hundred times wider.

That distinction is worth being blunt about. The Monte Carlo mean of 20.1620 sits 245 standard errors from the central-value calculation of 19.4975, and for about ten minutes we treated that as a disagreement. It is not. The central value uses pK2' = 9.800, the bottom edge of the uniform prior, whose own mean is 10.100. Recomputing the central value at 10.100 gives 20.1848, which sits two percent of one standard deviation from the Monte Carlo mean. The remainder is genuine curvature, since the ratio is nonlinear in pH and the average of a function is not the function of the average.

Reproducing This, and the Thing We Keep Coming Back To

Everything above comes from one file. From the repository root:

cd analysis
python buffer-capacity.py > buffer-capacity-output.txt

Python 3.12.3 and NumPy 2.4.2 were used here. NumPy supplies the PCG64 generator for the Monte Carlo and nothing else; the equilibrium work is standard-library arithmetic. The master seed is 20250321, stated in the first line of the output. Expected runtime is under a second on a school laptop, and the ten-minute budget the club sets for its studies was never in danger. Running the file twice produces byte-identical output. Passing --svg <directory> writes the five figures above as SVG fragments from the same arrays that produced the tables, so no figure in this article was drawn by hand or adjusted afterwards.

And then there is the thing none of us can put down.

We started this because buffer capacity looked like a formula worth computing properly, and it was. The peak really is at pKa, to four decimal places, and the small residual turns out to be water. The window really is 1.531 pH units wide, and 3+223 + 2\sqrt{2} really does fall out of a quadratic. Those are satisfying in the ordinary way that correct arithmetic is satisfying.

The carbon dioxide result is not like that. Blood is buffered by a weak acid whose dissociation constant is in the wrong place by a factor of twenty, at a concentration that ought to be far too low, and it works anyway, because the buffer is plumbed into a gas exchanger that runs about twelve times a minute and holds one of its two components at a constant value. The chemistry on its own is worth 2.94. The chemistry with the plumbing is worth 57.32. Every textbook diagram of the bicarbonate buffer that stops at the equilibrium arrows has drawn the less interesting half of the mechanism.

I keep thinking about the 3.85. Being 1.3 pH units away from your own optimum is supposed to be a design flaw. In an open system it is the whole point, because the distance is what fills the solution with bicarbonate. Sitting in the wrong place turns out to be the trick, provided someone is carrying the carbon dioxide away, which at this moment happens to be me.

References

Reproduce these numbers

From the analysis folder, run:

python buffer-capacity.py > buffer-capacity-output.txt
Analysis code
buffer-capacity.py71 kB