BBC Micro Compendium by Jeremy Ruston

discussion of beeb/electron applications, languages, utils and educational s/w
User avatar
hoglet
Posts: 7495
Joined: Sat Oct 13, 2012 6:21 pm
Location: Bristol
Contact:

BBC Micro Compendium by Jeremy Ruston

Post by hoglet » Fri Apr 27, 2018 12:17 pm

Hi all,

I'm trying to track down a copy of the BBC Micro Compendium by Jeremy Ruston.
ISBN 10: 0907563333
ISBN 13: 9780907563334

I think CCH has one, but I'm rather far from there:
http://www.computinghistory.org.uk/det/ ... ompendium/

I'm particularly interested in the part that covers "Computer arithmetic, including full details of floating point algorithms"

Does anyone have a copy of this book?

Dave
Last edited by hoglet on Fri Apr 27, 2018 12:41 pm, edited 1 time in total.

User avatar
leenew
Posts: 3661
Joined: Wed Jul 04, 2012 3:27 pm
Location: Doncaster, Yorkshire
Contact:

Re: BBC Micro Compendium by Jeremy Rushton

Post by leenew » Fri Apr 27, 2018 12:26 pm

Hi Dave,
I don't have the book, but I think the author was on here a while back unless I am mistaken?
p.s. it is Jeremy Ruston and not Rushton.

Lee.

User avatar
hoglet
Posts: 7495
Joined: Sat Oct 13, 2012 6:21 pm
Location: Bristol
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by hoglet » Fri Apr 27, 2018 12:42 pm

Oops, fixed that now.

I think his was the first annotated disassembly of BBC Basic:
https://twitter.com/jermolene/status/50 ... lang=en-gb

Dave

RobC
Posts: 2300
Joined: Sat Sep 01, 2007 9:41 pm
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by RobC » Fri Apr 27, 2018 1:33 pm

Pretty certain I have a copy Dave. I'll have a look in my library (aka 'the loft') and get back to you.

You're welcome to borrow it or I can try to dig out the info you want.

EDIT: I've found a copy. From memory, I may well have another version of the book - the one without the BASIC ROM listing that I think caused some issues BITD.

User avatar
hoglet
Posts: 7495
Joined: Sat Oct 13, 2012 6:21 pm
Location: Bristol
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by hoglet » Fri Apr 27, 2018 1:54 pm

That's great Rob.

What BigEd and I are interested in is any algorithmic details of the floating point trig/maths functions. Specifically, we have been digging into the implementation of LN() in Basic 4.00 and Basic 4.32, and we now understand these quite well. We'd like to extend this to Basic II as well.

We're collecting our notes in a github repository:
https://github.com/hoglet67/BBCBasic4r3 ... isassembly

If you have a moment, it would be great if you could take a quick look at what Jeremy says about LN() which is at &A801 in Basic II. A quick photo or scan of just those pages would be brilliant.

Dave

RobC
Posts: 2300
Joined: Sat Sep 01, 2007 9:41 pm
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by RobC » Fri Apr 27, 2018 2:03 pm

Will do.

In the fp arithmetic chapter it seems he gives pseudo code for some operations (addition, multiplication and sqrt) but doesn't do this for the log functions. However, the BASIC ROM disassembly is well annotated.

Looking at it, it's probably simplest if I just grab the relevant part of the BASIC II disassembly and type up the notes.

User avatar
hoglet
Posts: 7495
Joined: Sat Oct 13, 2012 6:21 pm
Location: Bristol
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by hoglet » Fri Apr 27, 2018 2:19 pm

RobC wrote: Looking at it, it's probably simplest if I just grab the relevant part of the BASIC II disassembly and type up the notes.
Whatever you think is easiest.

Many thanks.

User avatar
BigEd
Posts: 2091
Joined: Sun Jan 24, 2010 10:24 am
Location: West
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by BigEd » Fri Apr 27, 2018 2:32 pm

Thanks Rob!

RobC
Posts: 2300
Joined: Sat Sep 01, 2007 9:41 pm
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by RobC » Fri Apr 27, 2018 2:51 pm

Here's the typed up disassembly with Jeremy's notes. I've added the comments in square brackets.
Attachments
BASIC 2 LN routine.zip
(1.29 KiB) Downloaded 28 times

User avatar
BigEd
Posts: 2091
Joined: Sun Jan 24, 2010 10:24 am
Location: West
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by BigEd » Fri Apr 27, 2018 3:18 pm

Thanks Rob! What do we know about the "series evaluator" - is it a polynomial, maybe with a fixed number of terms?

Code: Select all

Load the address of the LN series table in A(lsb) and Y(msb):        
A838: LDA #&A873 AND 255         
A83A: LDY #&A873 DIV 256

Call the series evaluator:         
A83C: JSR &A897

User avatar
hoglet
Posts: 7495
Joined: Sat Oct 13, 2012 6:21 pm
Location: Bristol
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by hoglet » Fri Apr 27, 2018 3:24 pm

Ed,

It looks like the series evaluator for LN is driven off this table of values:

Code: Select all

.LA873
06
7a 1238a50b // 0.00892463796
88 790e9ff3 // 249.057128131
7c 2aac3fb5 // 0.04166817555
86 3401a27a // 45.0015963614
7f 638e37ec // 0.44444441562
82 3fffffc1 // 2.99999994133
Here's the code for it:

Code: Select all

.LA897
A897  85 4D              .M      STA &4D
A898  84 4E              .N      STY &4E
A89B  20 85 A3            ..     JSR &A385      // pack FWA into &46C to &470
A89E  A0 00              ..      LDY #&00
A8A0  B1 4D              .M      LDA (&4D),Y    // load count of fp table size
A8A2  85 48              .H      STA &48
A8A4  E6 4D              .M      INC &4D        // move to first fp var
A8A6  D0 02              ..      BNE &A8AA
A8A8  E6 4E              .N      INC &4E
.
A8AA  A5 4D              .M      LDA &4D
A8AC  85 4B              .K      STA &4B
A8AE  A5 4E              .N      LDA &4E
A8B0  85 4C              .L      STA &4C
A8B2  20 B5 A3            ..     JSR &A3B5      // unpack fp var into FWA
.
A8B5  20 F5 A7            ..     JSR &A7F5      // set 4B,4C to 046C (= value passed in)
A8B8  20 AD A6            ..     JSR &A6AD      // FWA = fp var / FWA normalised and rounded
A8BB  18                 .       CLC
A8BC  A5 4D              .M      LDA &4D
A8BE  69 05              i.      ADC #&05       // move to next fp variable
A8C0  85 4D              .M      STA &4D
A8C2  85 4B              .K      STA &4B
A8C4  A5 4E              .N      LDA &4E
A8C6  69 00              i.      ADC #&00
A8C8  85 4E              .N      STA &4E
A8CA  85 4C              .L      STA &4C
A8CC  20 00 A5            ..     JSR &A500      // FWA = FWA + fp var normalised and rounded
A8CF  C6 48              .H      DEC &48
A8D1  D0 E2              ..      BNE &A8B5
A8D3  60                 `       RTS
I'm just about to try to add this to our C Model:
https://github.com/hoglet67/BBCBasic4r3 ... log_test.c

Dave

User avatar
BigEd
Posts: 2091
Joined: Sun Jan 24, 2010 10:24 am
Location: West
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by BigEd » Fri Apr 27, 2018 3:27 pm

Thanks Dave!

RobC
Posts: 2300
Joined: Sat Sep 01, 2007 9:41 pm
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by RobC » Fri Apr 27, 2018 3:39 pm

BigEd wrote:Thanks Rob! What do we know about the "series evaluator" - is it a polynomial, maybe with a fixed number of terms?
Funnily enough, I was just reading about it...

The book says:
"Series evaluator
This routine is used by most of the transcendential [sic] functions to compute the series required. On entry, A(lsb) and Y(msb) must point to a series table. The table should contain a number (n) followed by n+1 constants.

This is a BASIC version of the series evaluator:

Code: Select all

10 REM Filename: SERDEMO
20
30 REM BASIC version of the series evaluator
40 REM
50 REM (c) 1983 Jeremy Ruston & Acorn
60
70 S=0
80 INPUT"Enter address:" A$
90 D%=EVAL(A$)
100 INPUT"Enter starting value:" FAC1
110
120 S46C=FAC1
130 I%=?D%:D%=D%+1
140 FAC1=FNnext
150 FAC1=S46C/FAC1
160 D%=D%+5
170 FAC1=FAC1+FNnext
180 I%=I%-1
190 IF I%<>0 THEN GOTO 150
200 PRINT FAC1
210 END
220
230 DEF FNnext
240 LOCAL T%
250 FOR T%=0 TO 3
260 ?(TOP +3 + T%)=D%?T%
270 NEXT T%
280 = S
In fact this can be simplified to: A/(A/(A/N1+N2)+N3)+N4

The above formula shows a case where there are 4 constants in the series. 'A' represents the contents of FAC#1 on entry."

So, I guess it's doing a continued fraction rather than a polynomial? According to the book, the seven constants used are:
8.92463796E-3, 249.057128, 4.16681756E-2, 45.0015964, 0.444444416, 2.99999994, -0.4999999999.

User avatar
BigEd
Posts: 2091
Joined: Sun Jan 24, 2010 10:24 am
Location: West
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by BigEd » Fri Apr 27, 2018 3:45 pm

Yes, looks like a continued fraction! Thanks. I think the continued fractions we see in Basic 4 versions are a slightly different form, presumably giving more accuracy with less calculation.

For log we saw this:
(z+z/[C1+(C2/(z^2)+C3/(C4+C5/(z^2)))])

Somehow I'd got the idea that early BBC Basic versions used the ratio of two polynomials, but that seems not to be the case. The Microsoft Basics, IIRC, use a polynomial, which Dave's experiments (and received wisdom) tell us is not the most efficient method.

RobC
Posts: 2300
Joined: Sat Sep 01, 2007 9:41 pm
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by RobC » Fri Apr 27, 2018 3:54 pm

I think that looks like Euler's formula (the 2z and z^2 ring a bell). I'll see if I can dig a reference out.

User avatar
BigEd
Posts: 2091
Joined: Sun Jan 24, 2010 10:24 am
Location: West
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by BigEd » Fri Apr 27, 2018 4:05 pm

We used wolframalpha.com to expand the expressions as a Taylor series, and one could then see how they come up with the right series. What that doesn't tell you is how to get from a Taylor series to a continued fraction! Nor how to tweak it, if necessary, to improve accuracy over a limited range.

Case 1: the continued fraction as used for ATN, as understood by wolframalpha:

"(z+z/[-1.8+(3/-(z^2)-.617142857*0/(-1.533333333+3/-(z^2)))]) series expansion"

returns the Taylor series
z - z^3/3 + 0.2 z^5 - 0.12 z^7 + 0.072 z^9 - 0.0432 z^11 + O(z^13)

(Notable that this is close to but not quite the Taylor series for Arctan. Perhaps the idea was to fit a truncated and adjusted polynomial??)

Case 2: as used for LN:

"-2(z+z/[-1.8+(3/z^2-.617142857*0/(-1.533333333+3/z^2))]) series expansion"

Answer:
-2 z - (2 z^3)/3 - 0.4 z^5 - 0.24 z^7 - 0.144 z^9 - 0.0864 z^11 + O(z^13)
(Taylor series)

Removing the -2 and asking again:

"(z+z/[-1.8+(3/z^2-.617142857*0/(-1.533333333+3/z^2))]) series expansion"

Answer:
z + z^3/3 + 0.2 z^5 + 0.12 z^7 + 0.072 z^9 + 0.0432 z^11 + O(z^13)
(Taylor series)

which is the same thing, but without the alternating signs - and that turns out to be what's needed, given the change of variables used for LN.

User avatar
hoglet
Posts: 7495
Joined: Sat Oct 13, 2012 6:21 pm
Location: Bristol
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by hoglet » Fri Apr 27, 2018 4:15 pm

Here's a quick C version of the Basic 2 LN() function in our test framework:
https://github.com/hoglet67/BBCBasic4r3 ... test.c#L55
(it's currently ignoring the exponent wrangling, so operates over a limited range of inputs)

And the results:

Code: Select all

BBC Basic 2                   : max_err = 7.599084e-11

BBC Basic 4r00                : max_err = 6.629041e-11

BBC Basic 4r32 with 1   coeffs: max_err = 2.105262e-07
BBC Basic 4r32 with 3   coeffs: max_err = 1.205916e-11
BBC Basic 4r32 with 1/3 coeffs: max_err = 1.205916e-11
Dave

User avatar
hoglet
Posts: 7495
Joined: Sat Oct 13, 2012 6:21 pm
Location: Bristol
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by hoglet » Fri Apr 27, 2018 4:38 pm

I've spotted a slight error in Jeremy's disassembly (or Rob's transposition of it):

Code: Select all

Set FAC#2 to one (notice that Y contains the exponent of FAC#2 after these instructions):        
A817: LDY #&80         
A819: STY &3B          
A81B: STY &3E          
A81D: INY              
A81E: STY &3D 
This is actually setting FAC#2 to minus one (-1) not one, which turns out to be quite important.

Dave

RobC
Posts: 2300
Joined: Sat Sep 01, 2007 9:41 pm
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by RobC » Fri Apr 27, 2018 4:59 pm

hoglet wrote:I've spotted a slight error in Jeremy's disassembly (or Rob's transposition of it)...
The error is in the original text.

User avatar
hoglet
Posts: 7495
Joined: Sat Oct 13, 2012 6:21 pm
Location: Bristol
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by hoglet » Fri Apr 27, 2018 5:09 pm

As pointed out earlier by Sydney, Jeremy is a member here:
http://www.stardot.org.uk/forums/member ... ile&u=9493

I've sent him a PM asking about the disassembly and the book, whether an electronic copy exists, and what his views are about making a copy available online somewhere. This looks like an excellent resource, so fingers crossed.

Dave

RobC
Posts: 2300
Joined: Sat Sep 01, 2007 9:41 pm
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by RobC » Fri Apr 27, 2018 5:41 pm

Pumping the BASIC II constants into Wolfram Alpha, I get a series expansion of: -0.5 + x/3 - (x^2)/4 + (x^3)/5 - (x^4)/6 + (x^5)/7 + O(x^6).

Multiplying by x^2 and adding x, as done in the code, gets us to the Taylor series for log(1+ x) (and shows why we add -1, not 1, to x at the start!).

So, the table is the continued fraction expansion of (log(1+x)-x)/(x^2).

User avatar
hoglet
Posts: 7495
Joined: Sat Oct 13, 2012 6:21 pm
Location: Bristol
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by hoglet » Fri Apr 27, 2018 6:02 pm

RobC wrote: So, the table is the continued fraction expansion of (log(1+x)-x)/(x^2).
So, knowing that, how would you actually generate the table values?

RobC
Posts: 2300
Joined: Sat Sep 01, 2007 9:41 pm
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by RobC » Fri Apr 27, 2018 6:16 pm

hoglet wrote:So, knowing that, how would you actually generate the table values?
I'm still trying to remember how you convert from one to another - I've got a text book somewhere that explains it but haven't found it yet!

However, Wolfram alpha will do it for you - do "((log(1+x)-x)/(x^2)) continued fraction". I've checked and the values it gives match but I could only get it to produce them using the Kn notation rather than as a fraction.

RobC
Posts: 2300
Joined: Sat Sep 01, 2007 9:41 pm
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by RobC » Fri Apr 27, 2018 8:10 pm

Think I may have figured out where the BASIC IV approximation came from. Looks like it's using the Pade approximation to log((1-x)/(1+x)).

Doing the third order Pade approximant of the function in Wolfram Alpha gives the rational function (30x - 8x^3)/(9x^2 - 15). Its series expansion is:
-2x - (2x^3)/3 - (2x^5)/5 - (6x^7)/25 - (18x^9)/125 - (54x^11)/625 + O(x^13).

This matches the series expansion for the BASIC IV continued fraction.
Last edited by RobC on Fri Apr 27, 2018 8:15 pm, edited 1 time in total.

User avatar
BigEd
Posts: 2091
Joined: Sun Jan 24, 2010 10:24 am
Location: West
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by BigEd » Fri Apr 27, 2018 8:13 pm

Well done!

User avatar
hoglet
Posts: 7495
Joined: Sat Oct 13, 2012 6:21 pm
Location: Bristol
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by hoglet » Sat Apr 28, 2018 8:38 am

Could one of the mods please enable PMs for JeremyRuston?

Cheers,

Dave

User avatar
Richard Russell
Posts: 503
Joined: Sun Feb 27, 2011 10:35 am
Location: Downham Market, Norfolk
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by Richard Russell » Wed Jun 13, 2018 2:35 pm

BigEd wrote:
Fri Apr 27, 2018 3:45 pm
Yes, looks like a continued fraction! Thanks. I think the continued fractions we see in Basic 4 versions are a slightly different form, presumably giving more accuracy with less calculation.
My attention was drawn to this thread by a post elsewhere, so apologies for the late reply. I can remember Sophie saying in the early days that she used a continued-fraction expansion for the transcendental functions. When I came to write my Z80 version of BBC BASIC I had to decide how I was going to do it; I initially considered using the same approach but there was one major snag: I didn't understand (then or now) how continued fractions work!

So I did some experiments to see how many terms I would need if I was to use a conventional (Chebychev-optimised) power-series expansion, and to my surprise it needed no more terms to achieve at least as good an accuracy as the 6502 version. So to this day I don't know what the attraction of the continued-fraction approach was, and all my versions of BBC BASIC (up to BB4W version 5) have used exactly the same power-series expansions. Of course in BB4W v6 and BBCSDL I've dropped support for 40-bit floats altogether so those routines are now gone.

An interesting historical anecdote is that the Chebychev-optimisation of the coeffecients was carried out by a very early version of BBC BASIC (Z80) which at that stage had no transcendental functions at all! Effectively it bootstrapped-itself, without needing any external assistance from another language or calculator. 8)

Richard.

User avatar
BigEd
Posts: 2091
Joined: Sun Jan 24, 2010 10:24 am
Location: West
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by BigEd » Sat Jun 16, 2018 5:05 pm

An interesting historical anecdote is that the Chebychev-optimisation of the coeffecients was carried out by a very early version of BBC BASIC (Z80) which at that stage had no transcendental functions at all! Effectively it bootstrapped-itself, without needing any external assistance from another language or calculator.
That is indeed very cool!

Also an interesting finding, that the Chebyshev approach is no more expensive and no less accurate. (My guess would be that his polynomials would come out more wiggly, even if they are in some sense optimal.) AFAICT the Padé approach is especially good when fitting something with a pole - I suppose LN does have a pole at zero, although I think we might be computing it near 1.

Dave has written a C program to mimic the LN computation, to help us understand what was going on. It's in the 4r32 repo.

I did run a short program to try to measure the accuracy of the various (6502) versions of BBC Basic:

Code: Select all

   10M=0:S=0
   20FORA=0.03125TO1.5STEP.0625/8
   30E=ABS(A-EXPLNA)/A
   40S=S+E
   50IF E>M M=E:L=A
   60NEXT
   70PRINT S
   80PRINT M,L
It prints the sum of absolute errors, and the worst error, in a roundtrip exp(ln(x)) for various x.

Here's what I got:
Basic 2: sum is 1.77639201E-8, worst is 1.24176343E-9
Basic 4: sum is 2.06907208E-8, worst is 8.27842289E-10
Basic 4.30 with fix: sum is 1.41910074E-8, worst is 4.38269447E-10
Basic 4.32 with fix: same as 4.30

So it seems to me they squeezed out an extra bit of accuracy each time they improved the algorithm. I believe Dave and I did count out the multiplications and divisions used, but I can't find the counts right now.

One of the tricks Acorn use is that they compute fewer terms of the Padé approximation for smaller arguments, which of course makes it a bit cheaper.

It would be mildly interesting to run a program like the above on one of your Z80 Basics.

User avatar
hoglet
Posts: 7495
Joined: Sat Oct 13, 2012 6:21 pm
Location: Bristol
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by hoglet » Sat Jun 16, 2018 5:14 pm

Ed,

It would be interesting to see the sum/max error without the fix for Basic 4.32.

Dave

User avatar
BigEd
Posts: 2091
Joined: Sun Jan 24, 2010 10:24 am
Location: West
Contact:

Re: BBC Micro Compendium by Jeremy Ruston

Post by BigEd » Sat Jun 16, 2018 5:25 pm

It's nothing to write home about:
Basic 4r32 as shipped: sum is 2.80552831E-4, worst is 2.61333737E-5

Post Reply