There is a 128 in my summation code, and I want to know who to thank for it.
The code adds numbers using a tree. Near the bottom, it switches to blocks of at most 128 values with eight running totals each. It's the arrangement used by NumPy's pairwise summation helper. I'm revising a paper on reproducible reductions and want one honest citation for the method I've borrowed.
It feels like a ten-minute job.
Inside that little function, I have three things to trace. First, the tree's improved error bound. Then the blocks that keep the function from spending too much time calling itself, and the eight running totals that give the processor independent additions to work on. Which citation covers them all? Finding a reference for one doesn't explain the other two.
The additions are the algorithm
A floating-point number has room for a limited number of significant digits. What happens when an addition needs more? Its result gets rounded. Moving the brackets changes which intermediate values get rounded, so it can change the final answer even when I feed in exactly the same numbers.
A running total keeps adding the next input to the previous result. An early input travels through nearly every addition, whereas in a balanced tree we add neighbouring inputs in pairs, then pair the results, and continue until one number remains.
Running total Balanced tree
+ d +
/ / \
+ c + +
/ / \ / \
a + b a b c d
Tree with larger leaves
+
/ \
leaf leaf
Inside each leaf, eight totals take turns receiving inputs:
lane 1: x0 + x8 + x16 + ...
lane 2: x1 + x9 + x17 + ...
...
lane 8: x7 + x15 + x23 + ...
Then combine the lane totals in a tree and add any tail values.
For a million inputs, a balanced binary tree has only about twenty levels. That shorter journey matters to rounding. Under the usual floating-point error model, the running total's worst-case bound grows roughly with the input count. The balanced tree's corresponding factor grows with the logarithm of that count.
Those bounds don't promise a tree wins on every input. Think of a ceiling in that error model. The actual error can be smaller.
Opposite-signed inputs can leave a relatively small result. That's cancellation. It can make particular examples behave differently. The tree limits how many rounded additions each input must pass through.
There is also a middle option. Split a million inputs into a thousand groups of a thousand, add each group with a running total, then add the thousand results. This is blocking. Why not stop there? It gives each input a much shorter journey without taking the tree all the way down to individual pairs.
I find the practical version appealing, because we can improve the error bound by moving brackets around while keeping the same inputs and ordinary additions. Every input still gets added. They don't all have to queue behind the same running total.
A textbook with no name attached
Nicholas Higham's chapter on summation says pairwise summation was “first discussed” by McCracken and Dorn, Babuška, and Linz. The references lead to 1964, 1969 and 1970.
Discussed. That's a careful word. I arrived wanting an inventor, and the standard reference is offering places to read.
McCracken and Dorn's Numerical Methods and FORTRAN Programming is the earliest statement I find in this search, with a tree for four inputs on printed page 62 and an error-bound comparison with a running total. They call the improvement “a result that is not intuitively obvious”.
What happens with more?
The next paragraph generalises by dividing n² positive numbers of similar magnitude into n groups of n, adding the values inside each group and then summing all n group totals. For four inputs, the picture is a tree. The general rule describes two levels of grouping. It doesn't say to recurse. Those are different instructions once the input gets larger.
That squared term deserves its own small celebration. My earlier quotations had lost it. Where did the other numbers go? With plain n the sentence doesn't account for all the numbers, while n² gives exactly the square-root blocking rule.
The useful part of the paragraph had disappeared into a superscript.
The grouping has no attribution. Does the surprised tone mean they discovered it themselves?
Actually, no, the tone can't tell me that by itself.
Two pages nearby make that argument fall apart. They also explain why adding small positive values first can help, in a similar voice without attribution, although John Carr had discussed ordering sums in his 1959 paper Error analysis in floating point arithmetic.
The textbook's tone tells me how its authors teach, not where each idea originated.
I like having that control case. A result already in print is firmer ground than interpreting how surprised two people sound on paper. Its tone is the same. I can't award an invention from that.
Another name, though. On page 56, McCracken and Dorn credit Kenneth M. King of Columbia University with suggesting the use of process graphs for error analysis. The diagrams trace intermediate results into later operations, letting us follow where errors accumulate. Six pages later, they compare the running total and tree that way.
Even his initials need checking. Stummel later calls him “H. M. King”, while the textbook spells out Kenneth M. King. I check the page image before following him.
Searching for the wrong initials would give me another reason to find nothing (and I already have plenty of those).
I can follow King too. The credit is for analysing the calculation, rather than grouping the additions. But it puts a person behind the pictures I've been studying and makes me wonder what else came out of those conversations, beyond the diagrams that reached the page.
Somebody must have used it
Since I can't find an inventor, I change the question. Somebody must have used this thing. If I can find an early user, their code might contain a comment or a reference that leads back to whoever showed it to them, so I start looking for a dated use I can follow. This involves a fair amount of dumpster diving through old source code. I feed huge quantities of FORTRAN, COBOL and other old code into AI, asking it to help find grouped sums or anything resembling the tree. The spelling needn't resemble the modern name. I'm looking for the calculation happening, not somebody announcing that they've invented pairwise summation.
NASA seems worth trying.
My thinking is that somebody doing all that numerical work must have had a reason to care about how a sum comes out. Surely there might be a routine tucked away somewhere, with an author and a date attached.
So I keep searching the old NASA code and related material.
Code from before the textbook would give me an earlier use. Shortly afterwards could help too.
Where did its author get the method? Did they cite the book, refer to an internal report, or describe it as something they already knew? Finding a user could give me another route towards an inventor.
I don't find a NASA example that gives me that trail.
It's frustrating, because I can explain why the search should help, but that doesn't produce the code I'm hoping for. Feeding more material to the AI mostly gets me more material.
The IBM library trail is no better. I can't find the grouping in the old FORTRAN maths-library material I search.
Was this a useful scribble in a textbook that nobody even picked up inside IBM?
The code scan finally turns up Graeme Fairweather's Algorithm 351, Modified Romberg Quadrature. It was submitted in September 1968 and published in June 1969. Blocked summation, at work. Its ROMINT routine estimates the area under a curve by grouping partial sums at three levels. The code has 32 and 512 in its loop bounds.
Here is the arrangement doing a job, with an author attached. Whom does he cite?
Fairweather calls it “a procedure proposed by Rutishauser to reduce the effect of rounding errors”. That sends me to Heinz Rutishauser's 1967 Description of ALGOL 60, pages 105–106, where he calls them presummations.
I have found a user and another possible inventor, but the source still comes after Dorn's textbook, so I have another trail that doesn't reach further back.
Linz's 1970 note is wonderfully brief. He explains the problem with one running sum, proposes multiple intermediate sums, and gives the pairwise recurrence. His bibliography contains Wilkinson's book on rounding errors.
I can follow the method perfectly well. Its origin remains much less tidy.
The footnote I nearly walked past
I try loads of Russian wordings, using AI to help crawl Russian-language sources. I'm hoping my English searches have missed its description.
The Soviet lead sends me round in a circle. First I dismiss it too quickly, then get excited about how much is there. I manage to be wrong about the same lead in opposite directions.
I. V. Viten'ko's 1968 paper examines optimal algorithms for addition and multiplication. It was submitted in April 1967. In its error model, minimising the depth of the additions gives an optimal error bound. How shallow can the calculation get? The minimum is the number of levels needed for a balanced binary tree.
Higham discusses that result too, separately from his list of early discussions. Optimal in what sense? In that error model. It establishes a precise relationship between the depth of the calculation and its error bound.
This is much richer than I expected from a citation hunt. Its 1967 submission still comes three years after the FORTRAN textbook. It cannot be the source of the 1964 passage, however much I want a neat trail leading backwards.
What did the search miss? Babuška's 1968 paper on numerical stability has a footnote on page 23 describing a summation error bound with logarithmic growth. The text search that dismissed the paper missed the very result I was looking for.
The page image makes it quite plain. Down at the bottom, beneath the main argument, there it is.
Searching scans is useful for finding a passage. It is a poor basis for saying a passage doesn't exist. Formulas lose their exponents, footnotes get mangled, and an idea can appear without the modern name I'm searching for. I now test an unsuccessful search against a passage I know exists. Can it find that? If not, I have a search problem.
At this point I'm less interested in choosing a winner. I want to understand why such a useful arrangement keeps appearing in remarks, examples and footnotes.
Further back, and further afield
I also wonder whether “first discussed by” hides a language barrier. An English-language author might know an earlier international source. Later readers might cite the English account and miss that connection. The inventor would disappear behind it. That would explain the discussions without a claim.
For a while, I'm convinced the answer is German, because another name for pairwise summation is cascade summation and “cascade” feels German to me. I have Zuse and his floating-point machines in mind too. I'm ready for this to turn into a German origin story.
Where does the word lead?
Following the word leads to Jack Wolfe's 1964 paper, with its “cascading accumulators”. Instead of pairing neighbouring inputs into the tree I'm chasing, his accumulators group values by magnitude, so the values themselves decide which group they belong in. The word names another method. A poor clue to my tree's birthplace.
Higham uses cascade as a name for pairwise summation too, so I still have a naming question. I can't join those uses into a history of the word. I was hoping a name would narrow the search. Instead, it has given me another thing to investigate.
A more substantial lead. Samelson and Bauer are studying floating-point accuracy in Munich in 1953. Compared with the older fixed-point papers I'm reading, that makes the German-speaking work look promising.
Which calculation are they analysing?
But their paper concerns individual operations and explicitly leaves error propagation through longer calculations aside. It doesn't give me the grouped sum.
Editions have histories too
Babuška leads me to his Czech colleagues, Milan Práger and Emil Vitásek. I follow the editions of their book, Numerical Processes in Differential Equations. The earlier the date, the more interested I become. I am getting quite easy to please at this point.
The earlier date is tempting, but it belongs to an earlier edition. What date can I attach to the passage? I find the summation material in the revised Russian edition of 1969. The English 1966 text I search doesn't turn it up. Earlier book, later passage. I can't put the earlier date beside material I've only located in the later edition, however much I'd like the dates on the books to do that work for me.
So finding an earlier edition hasn't given me an earlier statement. I still like the German theory, but I can't find the connection that would make it hold together.
I can't support the language-barrier explanation either. The Soviet result is later than the textbook, and the earlier German work doesn't supply the grouping. I never find an older passage supplying that connection.
One Czech lead stays open. Babuška's bibliography lists a 1963 paper in the Dresden journal, Numerische Stabilität von Rechenprozessen, by Babuška, Práger and Vitásek. I haven't obtained the paper. Its date is earlier than the textbook. Until I can read it, it stays a promising reference to chase.
A different arithmetic
In their 1947 paper on matrix inversion, von Neumann and Goldstine use fixed-point arithmetic. The fractional point stays put. Within the available range, adding two numbers on that grid introduces no rounding. Products need rounding while ordinary additions do not.
What addition rounding is there to repair?
Turing's 1948 matrix paper uses fixed-point calculations too. These authors haven't missed my sum's rounding problem. There is no addition-rounding error for regrouping to repair in their arithmetic. Range still needs careful attention, but that's a different job. I can stop expecting these authors to reach for my tree.
Or give the running total more digits than the inputs. Wilkinson discusses wider accumulators in Rounding Errors in Algebraic Processes. Those extra digits greatly reduce intermediate rounding before the final rounding to the required precision.
The cost depends on the machine and available operations. That helps explain why accuracy work might focus on the accumulator rather than grouping.
The textbook's additions use a fixed working precision. Regrouping can then improve the bound using ordinary additions. I can express its two-level scheme with nested loops. The change is accessible right where the book's reader is working.
It seems to me that this is why the search keeps resisting a grand invention story. Each author has a calculation to improve and a machine to improve it on. Regrouping needn't be their main result.
Dorn's other chain
There is one connection I particularly enjoy. In 1961, Dorn described a generalisation of Horner's rule for evaluating polynomials. He splits the polynomial's work into interleaved chains that can make progress together, rather than leaving the machine's parallel arithmetic units waiting behind one dependent chain.
What is the shared idea?
That resembles the eight lanes inside my summation leaf. The same calculation? No, Dorn is evaluating a polynomial. It has multiplications as well as additions. The connection is the dependency problem and the way he divides it, an older example of the same engineering thought without a documented chain of borrowing into NumPy.
So Dorn turns up in the speed discussion and, three years later, in the textbook's accuracy discussion. I came looking for a name to put beside one function. The same name appears beside two of its useful ideas, in different settings.
Did they know each other?
Once the papers stop naming predecessors, I start following the people, because if this came from a conversation then looking for an earlier publication might never find it. Right, who worked with whom? Where could somebody have passed the idea on?
General Electric looks promising. Both McCracken and Dorn worked there in the 1950s. McCracken's 1957 writing places him in the Computer Department in Phoenix while Dorn worked in the Aircraft Gas Turbine Division. Same company, different departments. The company name connects them but the records don't put them together. GE is a rather large place to go looking for a conversation.
At New York University, the connection gets closer. McCracken worked at its AEC Computing Center from 1958. The AEC is the US Atomic Energy Commission.
What about Dorn? His contemporary IBM biography places him at the same centre from 1957 to 1959 before he joined IBM. The dates and computing environment overlap. Here are the two authors before the textbook.
The funding records make me look twice. Dorn's 1960 paper on quadratic programming names the AEC contract AT(30-1)-1480. A decade later, Linz's summation paper carries the same contract number at Courant, NYU's institute for mathematical sciences.
I'm following a research community across time, with a government funding reference connecting papers on different subjects.
I have somehow made other people's employment histories part of adding up a list.
Invited to the same conference
The Czech trail reaches that community too. In June 1969, Vitásek and Widlund, a Courant researcher, are both invited speakers at the Dundee conference on the numerical solution of differential equations. Their papers appear in the same proceedings and Widlund's contribution acknowledges that AEC contract as well. Babuška's colleague and someone from Linz's institution are at the same event, a year before Linz's note appears. That's the sort of route I was hoping for. The idea could have travelled between people who weren't citing one another's papers. Did it? I don't have an account of anyone discussing this sum there. I don't have a record placing Linz at Dundee either. These connections make a conversation easier to imagine without telling me whether it happened.
The ten-minute job
By now, I've kept looking for much longer than is reasonable for the sentence I originally wanted to write. I have enough references to describe the method and explain why it helps. I could use them and get on with the paper, but instead I'm still trying to find out who first thought of moving the brackets.
I'm doing a remarkably thorough job of not writing that sentence.
The irritating part is that the search keeps giving me something worth following. A footnote changes what I think I know. An earlier paper answers a question I hadn't thought to ask. Even the failed searches start to make more sense once I understand the arithmetic those people were working with. There is always another reason to keep going.
But I really do want an answer. I want to stop opening another reference and finding the method sitting there, unexplained as a piece of history. Somebody wrote this down. Somebody drew that tree. Did they work it out themselves, or had somebody else shown it to them? By this point, leaving that unanswered bothers me far more than the original citation problem warrants.
The irony is that I can already use the result. Nothing in my code is waiting for me to settle its ancestry. I've turned a small piece of documentation into a separate research project, and I'm having trouble putting it down. The things I'm learning are interesting enough to keep me here, even while the answer I set out to find keeps getting further away.
Back to 128
I still want to know why the code stops splitting the input at 128 values.
The modern implementation history is easier to inspect. Steven G. Johnson's 2013 Julia proposal links to the Wikipedia explanation of pairwise summation. He says a base case of 128 makes the recursion overhead negligible in his comparison.
The NumPy change merged in December 2013 adds eight-way unrolling for accuracy and AVX vectorisation. Unrolling puts several iterations in one loop body, here using eight running totals. Vectorisation can handle several values per instruction. Its discussion moves comfortably between error, loop structure and speed. These are the concerns I recognise from my own code.
Johnson's 2017 explanation calls coarsening the base case “a standard technique to improve the performance of any recursive algorithm”. Here is someone who implemented it, describing the change as ordinary practice. I was hoping for a claim I could follow. He is explaining something programmers do to recursive functions.
Dorn's later paper cites IBM Research Report RC-647, Some mathematical aspects of parallel computation, by W. S. Dorn, N. C. Hsu and T. J. Rivlin, January 1962. Rivlin's bibliography lists it, but I haven't obtained its contents. If you have access to a copy, I'd love to read it. Please get in touch if you can share it or point me to an archive.
My current theory is that the first version comes out of that American circle, with Dorn at IBM, King at Columbia and the people around them. Drawing the operations made a useful rearrangement apparent, so they used it and carried on.
It also looks to me like five or so independent arrivals at versions of the same idea, apparently without knowing about one another's work. I keep looking for one inventor and finding more people who seem to have worked it out for themselves. It was so simple to its inventors that nobody bothered to claim it.
It is almost as if, in the 1960s, the world decided this must be invented. And so it was.