```Contributor: MIKE ANTTILA

(*
The two units used should come after this message. Uncomment several write-
commands to get a "fully" operational program rather than this benchmark
version. You then also can skip the Timer unit and the two commands from that
unit (TimerOn and TimerOff) to make the program much smaller (no float math
linked into the program).
*)

program PiCalc;  { The fastest PI calculator you'll ever find... :) }

{ From bits and pieces picked up mainly from the FidoNet PASCAL echo }
{ Collected, optimized, unitized, etc. by Bjorn Felten @ 2:203:208 }
{ Public Domain  --  Nov 1994 }

{ Units needed are at the end !! }

uses HugeUtil, Timer; { use Crt if you want fast printout on screen }
{ don't if you want to be able to redirekt o/p }

var
words, number   : longint;
nin, link, pii, a239    : HugePtr;

procedure ArcCoTan(n : integer; var angle : Huge);
var n2, del, remain : integer;
positive : boolean;

begin                               { corresp. integer operations }
ZeroHuge(angle,words);            { angle := 0 }
ZeroHuge(nin^,words);             { nin   := 0 }
ZeroHuge(link^,words);            { link  := 0 }
angle.dat[angle.len] := 1;        { angle := 1 }
DivHuge(angle,n,angle,remain);    { angle := angle div n }
n2 := n*n;                        { n2    := n * n }
del := 1;                         { del   := 1 }
positive := true;
CopyHuge(angle,nin^);             { nin   := angle }
repeat
DivHuge(nin^,n2,nin^,remain);   { nin   := nin div n2 }
inc(del, 2);                    { del   := del + 2 }
positive := not positive;
DivHuge(nin^,del,link^,remain); { link  := nin div del }
if positive then
AddHuge(angle,link^)          { angle := angle + link }
else
SubHuge(angle,link^);         { angle := angle - link }
{    write(#13,word(del)) } { uncomment to see that program is not dead }
until (link^.len <= 1) and (link^.dat[1] = 0);
{  writeln}                 { ... and this too }
end; { ArcCoTan }

begin
{  writeln('Program to get Pi (',pi:1:17,'...) with large precision.'); }
write('Digits(max 40.000): '); readln(number);
words := round(number / 4.7) + 3; { appr. 4.7 digits in one word }
write(number:6,#9);
TimerOn;
GetHuge(pii,  words+2);
GetHuge(a239, words+2);
GetHuge(link, words+2);
GetHuge(nin,  words+2);
ArcCoTan(5,   pii^);        { ATan(1/5)  }
AddHuge(pii^, pii^);
AddHuge(pii^, pii^);        { * 4        }
ArcCoTan(239, a239^);       { ATan(1/239)}
SubHuge(pii^, a239^);
AddHuge(pii^, pii^);
AddHuge(pii^, pii^);        { * 4        }
TimerOff;
{  WriteHuge(pii^, number)}     { uncomment if you want printout }
end.

unit HugeUtil;

interface

const HugeMax = \$8000-16;

type  Huge = record
len : word;
dat : array[1..HugeMax] of word;
end;
HugePtr = ^Huge;

procedure AddHuge  (var Answer, Add : Huge);
procedure MulHuge  (var A : Huge; Mul : integer; var Answer : Huge);
procedure DivHuge  (var A : Huge; Del : integer; var Answer : Huge;
var Remainder : integer);
procedure SubHuge  (var Answer, Sub : Huge);
procedure ZeroHuge (var L : Huge; Size : word);
procedure CopyHuge (var Fra,Til : Huge);
procedure GetHuge  (var P : HugePtr; Size : word);
procedure WriteHuge(var L : Huge; Size: word);

implementation

procedure AddHuge; assembler; asm
cld
push  ds
lds   di,Answer
les   si,Add
seges lodsw
mov   cx,ax
clc
@l1:
seges lodsw
adc   [si-2],ax
loop  @l1
jnb   @done
@l2:
add   word [si],1
inc   si
inc   si
jc    @l2
@done:
mov   si,di
lodsw
shl   ax,1
add   si,ax
lodsw
or    ax,ax
je    @d2
inc   word [di]
@d2:
pop   ds
end;

procedure MulHuge; assembler; asm
cld
push  ds
lds   si,A
mov   bx,Mul
les   di,Answer
mov   cx,[si]
mov   dx,si
inc   di
inc   di
clc
@l1:
mov   ax,[di]
pushf
mul   bx
popf
adc   ax,si
stosw
mov   si,dx
loop  @l1
adc   si,0
mov   es:[di],si
lds   di,A
mov   di,[di]
mov   ax,[di+2]
or    ax,ax
je    @l2
inc   di
inc   di
@l2:
lds   si,Answer
mov   [si],di
pop   ds
end;

procedure DivHuge; assembler; asm
std
push  ds
lds   si,A
mov   bx,Del
les   di,Answer
mov   cx,[si]
mov   di,cx
add   di,cx
xor   dx,dx
@l1:
mov   ax,[di]
div   bx
stosw
loop  @l1
lds   si,Remainder
mov   [si],dx
lds   si,A
mov   ax,[si]
lds   di,Answer
mov   [di],ax
mov   si,[di]
shl   si,1
@d3:
lodsw
or    ax,ax
jne   @d2
dec   word [di]
jne   @d3
inc   word [di]
@d2:
pop   ds
end;

procedure SubHuge; assembler; asm
cld
push  ds
lds   di,Answer
les   si,Sub
seges lodsw
mov   cx,ax
clc
@l1:
seges lodsw
sbb   [si-2],ax
loop  @l1
jnb   @done
@l2:
sub   word [si],1
inc   si
inc   si
jc    @l2
@done:
mov   si,[di]
shl   si,1
std
@d3:
lodsw
or    ax,ax
jne   @d2
dec   word [di]
jne   @d3
inc   word [di]
@d2:
pop   ds
end;

procedure WriteHuge;
var L1, L2, I, R, R1, X : integer;
begin
with L do begin
L1 := Len;
L2 := L1 - 1;
I := 1;
write(dat[L1],'.');
X := 0;
for I := 1 to Size div 4 do begin
Dat[L1] := 0;
Len := L2;
MulHuge(L,10000,L);
R := dat[L1];
R1 := R div 100;
R  := R mod 100;
write(chr(R1 div 10+48), chr(R1 mod 10+48),
chr(R  div 10+48), chr(R  mod 10+48));
inc(X);
write(' ');
if X > 14 then begin
writeln; write('  ');
X := 0
end
end
end;
writeln
end;                            { WriteHuge }

procedure ZeroHuge;
begin
fillchar(L.Dat, Size * 2, #0);
L.Len := Size
end;

procedure CopyHuge;
begin
move(Fra, Til, Fra.Len * 2 + 2)
end;

procedure GetHuge;
var D : ^byte;
Tries,
Bytes : word;
begin
Bytes := 2 * (Size + 1);
Tries:=0;
repeat
getmem(P,Bytes);

{ To make it possible to use maximally large arrays, and to increase
the speed of the computations, all records of type Huge MUST start
at a segment boundary! }

if ofs(P^) = 0 then begin
ZeroHuge(P^,Size);
exit
end;
inc(Tries);
freemem(P,Bytes);
new(D)
until Tries>10;   { if not done yet, it's not likely we ever will be }
writeln('Couldn''t get memory for array');
halt(1)
end;                                   { GetHuge }

end.

unit Timer;

interface

procedure TimerOn;
procedure TimerOff;

implementation

var
Time      : Longint absolute \$0040:\$006C;
WaitTime,
Temp      : Longint;

procedure TimerOn;
begin
WaitTime:=Time
end;

procedure TimerOff;
begin
Temp:=Time;
writeln('Done! It took ',(Temp-WaitTime)/18.2:6:2,'s.')
end;

end.
```