Sines and Cosines on the Datatrons - Part I

The reconstruction of SIN() and COS() for the Burroughs Algebraic Compiler on the Burroughs 205 was much easier once the SQRT() effort was complete and working.   So this might be a good place to talk about the structure of the subroutines and, more importantly, the invaluable assistance that Paul was able to locate for us.

I'll begin by noting that Paul Kimpel came into possession of copies of professor Don Knuth's initial "rough draft" source code for several of the Burroughs 205 subroutines back in 2015.

 

Paul now picks up the story of 

The Package from the Professor

 From Paul:

During 2014-15 I built a software emulator for the ElectroData/Burroughs 205. A major goal of that project was to recover the ALGOL-58 compiler that Professor Knuth wrote for that machine during 1960-61, a PDF scan of a listing for which was among the papers that he had donated to the Computer History Museum in 2005. I transcribed that listing, and after also recovering the EASY and MEASY assemblers that were required, eventually got the compiler to run. Alas, the listing had included only the compiler itself and not its run-time library, so while the compiler appeared to generate good code, there was no way for it to finish and produce an executable program.

Nonetheless, I was quite satisfied with the results of my effort, and wrote Professor Knuth in mid-November 2015, describing what I had done. He replied by email shortly thereafter, expressing enthusiasm for what I had accomplished, and kindly offered to send me some other items from his files. One of those was a copy of a listing of a non-trivial program he had written in 1960 for use by the Clevite Corporation of Cleveland, Ohio, where Professor Knuth was enrolled in undergraduate studies at what is now Case Western Reserve University. Clevite made bearings for automobile engines, and the program appears to be an analysis of forces on such bearings. Since the parts of a reciprocating engine articulate and rotate during operation, such a program would need to use a variety of math functions, particularly trigonometric functions.

I transcribed that listing and got the program to assemble using the EASY assembler, and started to puzzle out the format and units of the input data it would require, but that proved to be quite difficult. I set the program aside at some point and eventually forgot about it.

In July 2020, Tom Sawyer wrote me with a progress report on his effort to reconstruct the compiler's run-time library and raised an issue about the precision of the square root routine he had been working on. In the middle of replying to Tom, I suddenly remembered the Clevite program and its math routines. I suggested to Tom that he try extracting those Clevite routines for use in the ALGOL-58 compiler's library, which he subsequently did.

Once Paul jogged my memory, I immediately took a look at the Clevite program and discovered a set of math subroutines including SQRT, SIN, COS, ARCTAN, EXP and LOG.  This was wonderful; genuine "Knuth code" and they appeared to be about the right length for my existing 205 tape (see the final section of this post).

I implemented the SQRT, SIN and COS code from the Clevite program and was ready to move on to WRITE and READ which I knew were going to be difficult.  WRITE was around 300 words of rather complicated 205 code.  I thought that I would start with READ since there must be come similar code that I could borrow within the already working ALGOL compiler.  After all, it handled the processing of statements like "PI = 3.1415927$"

I started browsing the 80 page listing that professor Knuth had sent to Paul and was amazed to find two pages of EASY/MEASY assembler code devoted to the WRITE subroutine!  Just prior to this in the listing was source code for both SQRT and SIN.  What a wonderful gift to our project.

Here is an illustration of the code from both listings.  On the left is the original "rough draaft" source code for the sine calculation in the Burroughs 205 ALGOL-58 subroutines and on the right, the code in the Clevite program.  They are identical, except for two constants highlighted in blue.

In addition, the coefficients used in the calculation are omitted in the compiler listing but are present in the Clevite code.  We are fortunate to have both.  The Clevite listing even displays the coefficients twice, helping with a legibility problem.

The only difference between the two sets of code is that the Clevite subroutine expects the input to be expressed in degrees, and the ALGOL 58 program expects the input to be expressed in radians.  Interestingly both routines, as their first step, change the input value into fractions of a circle (or fractions of a sine cycle, if you prefer.)  

For a 45 degree input, the Clevite program multiplies by 1/360 = .0027777778 to get 0.125

The ALGOL58 program multiplies by 1/(2Pi) = 0.15915494 to also get 0.125

This step has some very convenient benefits.  For example, if the input were to be 7245 degrees (this is 20 x 360 +45)  the result becomes 20.125  By discarding the integer and only keeping the fraction, we have an input that has been scaled to a number between 0 and 1 representing a fraction of a circle.  The next 3 instructions do exactly that.  Let's call that number x.

The programs then proceed to evaluate the following polynomial in x:

 SIN(x) = Ax+Bx^3+Cx^5+Dx^7+Ex^9+Fx^11+Gx^13+Hx^15+Ix^17
where
A=   6.2831849
B= -41.341677
C=  81.604783
D= -76.701934
E=  42.040797
F= -15.046596
G=   3.7414651
H=  -0.64069755
I=   0.060675276

This polynomial results in quite accurate results at any point in the interval 0,1 (approximately 7 digits of accuracy.)  You can plug them into your favorite spreadsheet and check the results.

Note that the cosine subroutine simply consists of an alternate entry point one word previous to the sine entry point.  A single instruction there adds 90 degrees to the angle for the Clevite routine or Pi/4 radians to the angle in the ALGOL-58 language routine.  This takes advantage of the well-known relationship

COS(x)=SIN(x + 90 degrees)

 

This is as about as much that might be of interest to a programmer looking at the recovery of SIN(x) and COS(x) on the 205.  But what about the mathematics involved?  Where do those polynomial coefficients come from?  How do other people calculate sines and cosines.  And, how does this implementation compare to the one in BALGOL on the Datatron 220?  I will take that up in a separate post, Part II.

Meanwhile, lets look at the actual layout of the compiler and subroutines on the magnetic tape that was distributed as Burroughs Algebraic Compiler back in 1960. 

The Layout of the Compiler Tape

The full compiler occupies 276 blocks on a 10,000 block lane of a 205 tape.  That is 276 blocks times 20 words per block x 11 digits per word equals 60,720 digits.  If you count each word as roughly 5 bytes plus sign, it comes to 27,600 bytes.  Stuffing a pretty full implementation of ALGOL 58 into that compact an area was quite an accomplishment.

Here is breakdown of the compiler segments on the tape:

Loader:  This is a single block at the start which loads the next 200 blocks into memory and runs a checksum to see if everything is good to this point.  This takes two tape read commands since 100 blocks is the maximum to be read by a single Burroughs 205 command.  The tape reads require about ten seconds.  The loader then performs a checksum on the 4000 words that were loaded which takes an additional 35 seconds.  (A very useful hack when you are doing multiple compiles is to skip the checksum!  There might have been some hazard to this back in the 1960s but is pretty low risk today.  Even back then, you were much more likely to simply fail to read the tape than read it incorrectly.)

The Compiler:  The next two hundred blocks on the tape, shown in red, comprise the 4000 words that are the actual compiler.  You might notice a fair amount of empty space in the last quarter of that area.  That is where tables are going to be stored.  This code will read all of the source program and generate object code from the first card through your FINISH$ statement.  Since all of this work is done in a single pass, address resolution is left as a project for the object code loader later.  Some syntax checking takes place during compilation.  If you misspell a reserved word, it might be immediately caught - for example, GOTO instead of GO TO or a second OUTPUT statement with the same name as a previous one.  Other egregious errors like including a semicolon in front of an END will be caught on the spot and whenever errors are caught, compilation simply stops.  It can take multiple tries to get to the end of the program the first time and careful desk-checking is advised!

The Compiler Overlay:  Once the compiler has processed a program and detected a FINISH$ statement, the next 22 blocks, shown in orange, are read in overlaying 440 words beginning at location 2200.  A check is immediately run to see if the sum of four counters is correct.  If your program passes this test, you are on your way!  This overlay will then load format bands for the CARDATRON reader, printer and punch onto the object program tape and finally, load any necessary subroutines from the following segments.  If your program fails that first test, however, the compiler reports the dreaded ".1001" on the Flexowriter.  This means "Unmatched left parenthesis or missing right operand.  Probably goes back far into the program where the error actually occurred."

The Subroutine Segments:  There were initially six subroutine segments on the compiler tape:

1. Floating to Fixed and Fixed to Floating.  Calls to these subroutines are automatically generated when the compiler sees a statement like "X = 1$" and X has not been declared as an INTEGER.
2. LOG(X) and EXP(X)
3. SQRT(X), COS(X), SIN(X) and ARCTAN(X)
   
4. WRITE($$OUT,FORM)
5. READ($$IN)
6. PTREAD($$IN), PTWRITE($$OUT)

As can be seen in the picture above, my tape dump ran out of data a bit more than halfway through the LOG(X) and EXP(X) segment.

A seventh segment was added to the compiler tape later at the end of 1961 by Burroughs to incorporate magnetic tape access using the MTSEARCH(), MTREAD() and MTWRITE() calls.  We have the good fortune to have the actual subroutines and the instructions for adding theses three new subroutines to the compiler because of the initiative of  Paul Kimpel reaching out to Professor Knuth.  This proved invaluable for dealing with the missing READ and WRITE routines.  They also serve as an interesting template for anyone thinking about adding some new input or output subroutines.

Subroutine Lengths and Entry Points Overview

There are three tables within the compiler that relate information about the subroutine lengths.  One specifies the length of the tape segments in blocks, a second notes how many words will be loaded from each segment and a third defines the entry point of a specific subroutine relative to the segment.  For the SIN() function, it resides on a 7 block segment;  110 words will be loaded from that segment; and

SQRT()        begins at word 0
COS()          begins at word 26
SIN()           begins at word 27
ARCTAN() begins at word 57

all within the 110 word segment.  This also tells us that the COS/SIN code should be 31 words in length, extending from word 26 through 56.  It was very welcome news to discover that this was the length of the code found in the Knuth papers.