Pi: Difference between revisions

Content added Content deleted
(Added solutions for Delphi and EDSAC.)
Line 1,459: Line 1,459:
4287554687311595628638823537875937519577818577805321712268066130019278
4287554687311595628638823537875937519577818577805321712268066130019278
76611195909216420198</pre>
76611195909216420198</pre>

=={{header|Delphi}}==
{{trans|BBC BASIC}}
250,000 digits of pi on a BBC micro is pretty impressive, so a translation of the BBC code has been used for this Delphi 7 solution. A Delphi TMemo replaces the BBC screen, and the output is copied from the TMemo to a disk file at the end of the program.

The number of digits written depends on the variable M, and is found by experiment to be very close to 2*M/7. The value of M is limited by the danger of overflow in D, which can approach (but doesn't exceed) 400*M.

With M = floor( (2^31 - 1)/400 ), the Delphi version writes 1,533,913 decimal places, which have been verified against [sorry, defeated by the captcha as usual]

introcs dot cs dot princeton dot edu slash java slash data slash pi-10million.txt
<lang delphi>
unit Pi_BBC_Main;

interface

uses
Classes, Controls, Forms, Dialogs, StdCtrls;

type
TForm1 = class(TForm)
btnRunSpigotAlgo: TButton;
memScreen: TMemo;
procedure btnRunSpigotAlgoClick(Sender: TObject);
procedure FormCreate(Sender: TObject);
private
fScreenWidth : integer;
fLineBuffer : string;
procedure ClearText();
procedure AddText( const s : string);
procedure FlushText();
end;

var
Form1: TForm1;

implementation

{$R *.dfm}

uses SysUtils;

// Button clicked to run algorithm
procedure TForm1.btnRunSpigotAlgoClick(Sender: TObject);
var
// BBC Basic variables. Delphi longint is 32 bits.
B : array of longint;
A, C, D, E, I, L, M, P : longint;
// Added for Delphi version
temp : string;
h, j, t : integer;
begin
fScreenWidth := 80;
ClearText();
M := 5368709; // floor( (2^31 - 1)/400 )

// DIM B%(M%) in BBC Basic declares an array [0..M%], i.e. M% + 1 elements
SetLength( B, M + 1);
for I := 0 to M do B[I] := 20;
E := 0;
L := 2;

// FOR C% = M% TO 14 STEP -7
// In Delphi (or at least Delphi 7) the step size in a for loop has to be 1.
// So the BBC Basic FOR loop has been replaced by a repeat loop.
C := M;
repeat
D := 0;
A := C*2 - 1;
for P := C downto 1 do begin
D := D*P + B[P]*$64; // hex notation copied from BBC version
B[P] := D mod A;
D := D div A;
dec( A, 2);
end;

// The BBC CASE statement here amounts to a series of if ... else
if (D = 99) then begin
E := E*100 + D;
inc( L, 2);
end
else if (C = M) then begin
AddText( SysUtils.Format( '%2.1f', [1.0*(D div 100) / 10.0] ));
E := D mod 100;
end
else begin
// PRINT RIGHT$(STRING$(L%,"0") + STR$(E% + D% DIV 100),L%);
// This can't be done so concisely in Delphi 7
SetLength( temp, L);
for j := 1 to L do temp[j] := '0';
temp := temp + SysUtils.IntToStr( E + D div 100);
t := Length( temp);
AddText( Copy( temp, t - L + 1, L));
E := D mod 100;
L := 2;
end;
dec( C, 7);
until (C < 14);
FlushText();

// Delphi addition: Write screen output to a file for checking
h := SysUtils.FileCreate( 'C:\Delphi\PiDigits.txt'); // h = file handle
for j := 0 to memScreen.Lines.Count - 1 do
SysUtils.FileWrite( h, memScreen.Lines[j][1], Length( memScreen.Lines[j]));
SysUtils.FileClose( h);
end;

{=========================== Auxiliary routines ===========================}

// Form created
procedure TForm1.FormCreate(Sender: TObject);
begin
fScreenWidth := 80; // in case not set by the algotithm
ClearText();
end;

// This Delphi version builds each screen line in a buffer and puts
// the line into the TMemo when the buffer is full.
// This is faster than writing to the TMemo a few characters at a time,
// but note that the buffer must be flushed at the end of the program.
procedure TForm1.ClearText();
begin
memScreen.Lines.Clear();
fLineBuffer := '';
end;

procedure TForm1.AddText( const s : string);
var
nrChars, nrLeft : integer;
begin
nrChars := Length( s);
nrLeft := fScreenWidth - Length( fLineBuffer); // nr chars left in line
if (nrChars <= nrLeft) then fLineBuffer := fLineBuffer + s
else begin
fLineBuffer := fLineBuffer + Copy( s, 1, nrLeft);
memScreen.Lines.Add( fLineBuffer);
fLineBuffer := Copy( s, nrLeft + 1, nrChars - nrLeft);
end;
end;

procedure TForm1.FlushText();
begin
if (Length(fLineBuffer) > 0) then begin
memScreen.Lines.Add( fLineBuffer);
fLineBuffer := '';
end;
end;

end.
</lang>
{{out}}
<pre>
3.141592653589793238462643383279502884197169399375105820974944592307816406286208
99862803482534211706798214808651328230664709384460955058223172535940812848111745
02841027019385211055596446229489549303819644288109756659334461284756482337867831
...
</pre>

=={{header|EDSAC order code}}==
The task is impossible if taken literally, because a program with finite memory cannot generate a non-repeating infinite sequence.

Like many of the solutions posted here, this program implements the spigot algorithm of Rabinowitz and Wagon. The code takes up 192 EDSAC locations, leaving 832 for storage, which is enough for 252 correct digits of pi.

With only a few changes (noted in the comments), the same code can be used to
output the first 2070 digits of e.
<lang edsac>
[EDSAC program, Initial Orders 2.
Calculates digits of pi by spigot algorithm.
Easily edited to calculate digits of e (see "if pi" and "if e" in comments).
Based on http://pi314.net/eng/goutte.php
See also https://www.cut-the-knot.org/Curriculum/Algorithms/SpigotForPi.shtml

Uses 17-bit values throughout.
Array index and counters are stored in the address field,
i.e. with the least significant bit in bit 1.
For integer arithmetic, the least significant bit is bit 0.

Variables don't need to be initialized at load time, so they overwrite
locations 6..12 of initial orders to save space.]
[6] P F [index into remainder array]
[7] P F [carry in the spigot algorithm]
[8] P F [negative count of digits]
[9] P F [pending digit, always < 9]
[10] P F [negative count of pending 9's]
[11] P F [9 or 0 in top 5 bits, for printing]
[12] P F [negative count of characters in current line]

[Array corresponding to Remainder row on http://pi314.net/eng/goutte.php]
T 53 K [refer to array via B parameter]
P 13 F [start array immediately after variables]

[Subroutine for short (17-bit) integer division.
Input: dividend at 4F, divisor at 5F.
Output: remainder at 4F, quotient at 5F.
Working locations 0F, 1F. 37 locations.]
T 987 K
GKA3FT34@A5FUFT35@A4FRDS35@G13@T1FA35@LDE4@T1FT5FA4FS35@G22@
T4FA5FA36@T5FT1FAFS35@E34@T1FA35@RDT35@A5FLDT5FE15@EFPFPD

T 845 K
G K
[Constants]
[Enter the editable numbers as addresses, e.g. P 100 F for 100.
Reducing the maximum array index will make the program take less time.
A maximum index of 831 (i.e. using all available memory) will give
252 correct digits of pi, or 2070 correct digits of e.]
[0] P 831 F [maximum array index <-- EDIT HERE, don't exceed 831]
[1] P 252 F [number of digits <-- EDIT HERE]
[2] P 72 F [digits per line <-- EDIT HERE]
[3] P D [short-value 1]
[4] P 5 F [short-value 10]
[5] J F [10*(2^12)]
[6] X F [add to T order to make V order]
[7] M F [used to convert T order to A order]
[8] # F [figures shift (for printer)]
[9] @ F [carriage return]
[10] & F [line feed]

[Main routine. Enter with acc = 0.]
[11] O 8 @ [set printer to figures; also used to print '9']
S 2 @ [load negative characters per line]
T 12 F [initialize character count]
S 1 @ [load negated number of digits]
T 8 F [initialize digit count]
T 10 F [clear negative count of 9's]
S 2 F [load -2 (any value < -1 would do)]
T 9 F [initialize digit buffer]
[Start algorithm: fill the remainder array with 2's (or 1's for e)
The code is a bit neater if we work backwards.]
A @ [maximum index]
[20] A 81 @ [make T order for array entry]
T 23 @ [plant in code]
[if pi] A 2 F [acc := 2]
[if e] [A 3 @] [acc := 1]
[23] T F [store in array entry]
A 23 @ [dec address in array]
S 2 F
S 81 @ [finished array?]
E 20 @ [loop back if not

Outer loop. Here for next digit.]
[28] T F [clear acc]
[Multiply remainder array by 10.
NB To preserve integer scaling, we need product times 2^16.]
H 5 @ [mult reg := 10*(2^12)]
A @ [acc := maximum index]
[31] A 81 @ [make T order for array entry]
U 37 @ [plant in code]
A 6 @ [convert to V order, same address]
T 35 @ [plant in code]
[35] V F [acc := array entry * 10*(2^12)]
L 4 F [shift to complete mult by 10*(2^16)]
[37] T F [store result in array]
A 37 @ [load T order]
S 2 F [dec address]
S 81 @ [test for done]
E 31 @ [loop back if not]
T F [clear acc]
T 7 F [clear carry]
A @ [acc := maximum index]
T 6 F [initialize array index]
[Inner loop to get next digit.
Work backwards through remainder array.]
[46] T F [clear acc]
A 6 F [load index]
A 81 @ [make T order for array entry]
U 61 @ [plant in code]
A 7 @ [convert to A order]
T 52 @ [plant in code]
[52] A F [load array element]
A 7 F [add carry from last time round loop]
T 4 F [sum to 4F for division routine]
A 6 F [acc := index as address = 2*(index as integer)]
[if pi] A 3 @ [plus 1]
[if e] [R D] [shift right, address --> integer]
T 5 F [to 5F for division routine (for e, 5F = index/2)]
[58] A 58 @ [call routine to divide 4F by 5F]
G 987 F
A 4 F [load remainder]
[61] T F [update element of remainder array]
[if pi: 4 orders]
H 5 F [mult reg := quotient]
V 6 F [multiply by index
NB need to shift 15 left to preserve integer scaling]
L F [shift 13 left]
L 1 F [shift 2 more left (for e, just use quotient)
[if e: 1 order, plus 3 no-ops to keep addresses the same]
[A 5 F] [load quotient]
[XF XF XF]
T 7 F [update carry for next time round loop]
A 6 F [load index]
S 2 F [decrement]
U 6 F [store back]
[We want to terminate after doing index = 1]
S 2 F [dec again]
E 46 @ [jump back if index >= 1]
[Treatment of index = 0 is different]
T F [clear acc]
A B [load rem{0)]
A 7 F [add carry]
T 4 F [sum to 4F for division routine]
A 4 @ [load 10]
T 5 F [to 5F for division routine]
[78] A 78 @ [call division routine]
G 987 F
A 4 F [load remainder]
[81] T B [store in rem{0}; also used to manufacture orders]
[82] A 82 @ [call subroutine to deal with quotient (clears acc)]
G 93 @
A 8 F [load negative digit count]
A 2 F [increment]
U 8 F [store back]
G 28 @ [if not yet 0, loop for next digit]
[Fake a zero digit to flush the last genuine digit(s)]
T 5 F [store fake digit 0 in 5F]
[89] A 89 @ [call subroutine to deal with digit]
G 93 @
O 8 @ [set figures: dummy character to flush print buffer]
Z F [stop]

[Subroutine to handle buffering and printing of digits.
Here with quotient from spigot algorithm still in 5F.
The quotient at 5F is usually the new decimal digit.
But the quotient can be 10, in which case we must treat it as 0
and ripple a carry through the previously-computed digits.
Hence the need for buffering.]
[93] A 3 F [make and plant return link]
T 130 @
A 5 F [load quotient]
S 4 @ [subtract 10]
E 105 @ [jump if quotient >= 10]
A 3 @ [add 1]
G 109 @ [jump if quotient < 9]
[Here if quotient = 9. Update count of 9's,
don't do anything with the buffer.]
T F [clear acc]
A 10 F [load negative count of 9's]
S 2 F [subtract 1]
T 10 F [update count]
E 130 @ [exit with acc = 0]
[Here if quotient >= 10. Take digit = quotient - 10,
and ripple a carry through the buffered digits.]
[105] T 5 F [store (quotient - 10) formed above]
T 11 F [store 0 to print '0' not '9']
A 3 @ [add 1 to buffered digit]
E 112 @ [join common code]
[Here if quotient < 9. Flush the stored digits.]
[109] T F [clear acc]
A 11 @ [load any O order (code for O is 9)]
T 11 F [store to print '9']
[112] A 9 F [load buffered digit, plus 1 if quotient >= 10]
G 118 @ [skip printing if buffer is empty]
L 1024 F [shift digits to top 5 bits]
T 1 F [store in 1F for printing]
[116] A 116 @ [call print routine]
G 131 @
[118] T F [clear acc]
A 5 F [load quotient (possibly modified as above)]
T 9 F [store in buffer]
[121] A 10 F [load negative count of 9's]
E 130 @ [if none, exit with acc = 0]
A 2 F [inc count]
T 10 F
A 11 F [load 9 (or 0 if there's a carry)]
T 1 F [to 1F for printing]
[127] A 127 @ [call print routine (clears acc)]
G 131 @
E 121 @ [jump back (always)]
[130] E F [return to caller]

[Subroutine to print character at 1F.
Also prints CR LF if necessary.]
[131] A 3 F [make and plant link for return]
T 141 @
A 12 F [load negative character count]
G 138 @ [jump if not end of line]
S 2 @ [reset character count]
O 9 @ [print CR LF]
O 10 @
[138] O 1 F [print character]
A 2 F [add 1]
T 12 F
[141] E F
E 11 Z [define entry point]
P F [acc = 0 on entry]
</lang>
{{out}}
<pre>
314159265358979323846264338327950288419716939937510582097494459230781640
628620899862803482534211706798214808651328230664709384460955058223172535
940812848111745028410270193852110555964462294895493038196442881097566593
344612847564823378678316527120190914
</pre>


=={{header|Elixir}}==
=={{header|Elixir}}==