Umiestnenie súboru www.TrSek.com/pas/bignum.pas

{ BIGNUM.PAS                 Copyleft (c) Zdeno Sekerák, <etrsek@gmail.com> }
{                                                                           }
{ Tento program je slobodný software. Mo¾ete ho ïalej distribuova» a/alebo  }
{ upravova» pod podmienkou licence GNU General Public License vydanej       }
{ organizáciou Free Software Foundation, verzia licencie 3 alebo vy¹¹ej.    }
{                                                                           }
{ Tento program je distribuvaný v nádeji, ¾e bude u¾itoèným, ale            }
{ NEPOSKYTUJE ®IADNE ZÁRUKY. Bez akejkoµvek vyplývajúcej záruky na          }
{ OBCHODOVATE¥NOS« alebo VHODNOS« PRE KONKRÉTNE POU®ITIE. Pre viac          }
{ podrobností si preèítajte licenciu GNU General Public Licence.            }
{                                                                           }
{  http://www.gnu.org/copyleft/gpl.html                                     }
{                                                                           }
{ Author: Zdeno Sekerák (TrSek)                                             }
{ Datum : 08.01.2010                                                        }
{ Verzia: 0.91 RC1                                     http://www.trsek.com }
unit BigNum;

interface
uses Objects, Crt;

{ make testing computing in native format, use for debuging                 }
{ for release use UNDEF _BN_NATIVE                                          }
{ DEFINE _BN_NATIVE}

const BN_MAX_DIGIT     = 240;           { how many digital has number }
      BN_MAX_STRING    = 240;            { maximal 255-string limit    }
      BN_INFINITY      = BN_MAX_DIGIT;
      BN_NUMSYSTEM     = 10;            { numeral system }
      BN_EXPONENT_ZERO = 5;             { how many zero write when format 1E-BN_EXPONENT_ZERO }
      BN_DECIMAL_POINT = '.';

      BN_ARROW_LEFT    = #75;
      BN_BACKSPACE     = #8;
      BN_ENTER         = #13;
      BN_ESC           = #27;
      BN_CHAR_SIGN     = '-';

      { error }
      BNE_OK                = 0;
      BNE_USER_TERMINATE    = 1;
      BNE_MINIMAL_OVERFLOW  = 2;
      BNE_MAXIMAL_OVERFLOW  = 3;
      BNE_DIVIDE_ZERO       = 4;
      BNE_SQRT_MINUS        = 5;
      BNE_EXP_MINUS         = 6;
      BNE_ARCSINUS_BIG      = 7;
      BNE_LN_NONPOSITIVE    = 8;

      { for unit needs - do not change }
      BN_MANTISA_LESS  = -1;
      BN_MANTISA_EQUAL = 0;
      BN_MANTISA_HIGH  = 1;

      BN_MIN_EXPONENT  = -2147483647;   { min longint }
      BN_MAX_EXPONENT  =  2147483647;   { max longint }

type
  digit = set of 0..9;

  _bignumber_type = record
    sign     : boolean;          { sign of number plu/minus        }
    exponent : longint;          { exponent of number              }
    overflow : byte;             { something wrong during compute  }
    zero     : boolean;          { indicate zero }
{$IFDEF _BN_NATIVE}
    native   : extended;         { its for test in native format   }
{$ENDIF}
    mantisa  : array[1..BN_MAX_DIGIT] of byte;    { numbers }
  end;

  BigNumber = Object(TObject)
    _bn: _bignumber_type;

    public
    constructor init(s:string);
    destructor  done; virtual;

    { important and often use }
    private
    procedure ExponentLeft(step:integer);
    procedure ExponentRight(step:integer);
    procedure Normalize;
    procedure Clear;
    function  IsZero:boolean;

    { auxiliary }
    public
    {protected}
    function  compare_mantisa(bn: bignumber):integer;
    function  Min(a,b:integer):integer;
    function  ToInt(s:string):longint;
    function  ToStr(i:longint):string;
    procedure SplitOverFlow(add: byte);

    { input, output }
    public
    procedure Frac;
    procedure Trunc;
    function  GetErrorS:string;
    procedure Val(s:string);
    function  Str(n,d:integer):string;
    function  Read(ch:char):char;
    procedure Write(n,d:integer);

    { 3 (resp.4) very important mathematical operations }
    public
    procedure plus     (add: bignumber);    { addend   }
    procedure minus    (les: bignumber);    { less     }
    procedure multiply (sub: bignumber);    { miltiply }
    procedure divide   (del: bignumber);    { divide   }

    { don't work with Turbo Pascal 7.0
    function operator <+> (a,b: bignumber);
    function operator <-> (a,b: bignumber);
    function operator <*> (a,b: bignumber);
    function operator </> (a,b: bignumber);
    }
  end;

implementation
{---------------------------------------------------------------------------}

{ inicialize constructor }
constructor bignumber.init(s:string);
begin
  self.Clear;

  if(length(s)>0)then
     self.Val(s);
end;
{---------------------------------------------------------------------------}

{ destructor, nothing do }
destructor bignumber.done;
begin
end;
{---------------------------------------------------------------------------}

procedure bignumber.Clear;
var i:integer;
begin
  _bn.sign     := false;
  _bn.exponent := 0;
  _bn.overflow := BNE_OK;
  _bn.zero     := true;
{$IFDEF _BN_NATIVE}
  _bn.native   := 0;
{$ENDIF}

  for i:=1 to BN_MAX_DIGIT do
    _bn.mantisa[i]:=0;
end;
{---------------------------------------------------------------------------}

function bignumber.Read(ch:char):char;
var i: integer;
{  ch: char;}
  x,y: integer;
    s: string[BN_MAX_STRING];
begin
  x:=wherex;
  y:=wherey;
  s:='';

  if(ch in [' ',#0])then
     ch:=readkey;

  repeat
    { number, sign char on start }
    if((ch in ['0'..'9'])
    or ((ch=BN_CHAR_SIGN) and (length(s)=0)))then begin
       s:=s+ch;
       system.write(ch);
       { why, until }
       if(ch=BN_CHAR_SIGN)then
          ch:='0';
    end;

    { decimal point only 1 times }
    if((ch in [',','.',BN_DECIMAL_POINT])
    and (Pos(BN_DECIMAL_POINT, s)=0))then begin
       s:=s + BN_DECIMAL_POINT;
       system.write(ch);
    end;

    { change mind - arrow to left, backspace }
    if(((ch=BN_ARROW_LEFT) or (ch=BN_BACKSPACE))
    and (length(s)>0))then begin
       delete(s,length(s),1);
       gotoxy(x,y); system.write(s+' ');
       gotoxy(x,y); system.write(s);
    end;

    { next char }
    ch:=readkey;

{ until (ch in [BN_ENTER, BN_ESC]);}
  until not(ch in ['0'..'9', ',', '.',
                   BN_DECIMAL_POINT, BN_ARROW_LEFT, BN_BACKSPACE]);
{ writeln;}

  { delete }
  Read:=ch;
  if((ch=BN_ESC) or (length(s)=0))then begin
     self.Clear;
     exit;
  end;

  { convert to bignumber format }
  self.Val(s);
end;
{---------------------------------------------------------------------------}

{ CARE: moze sa stratit cislo }
procedure bignumber.ExponentLeft(step:integer);
var i:longint;
begin
  if(step=0)then
     exit;

  { podtecenie }
  if(_bn.exponent < (BN_MIN_EXPONENT + step))then
     _bn.overflow := BNE_MINIMAL_OVERFLOW;

  { exponent change }
  _bn.exponent := _bn.exponent-step;

  { mistake, lost value !? }
  for i:=1 to step do
    if(_bn.mantisa[i]<>0)then
       _bn.overflow := BNE_MINIMAL_OVERFLOW;

  for i:=1 to BN_MAX_DIGIT-step do
    _bn.mantisa[i] := _bn.mantisa[i+step];

  for i:=BN_MAX_DIGIT-step+1 to BN_MAX_DIGIT do
    _bn.mantisa[i]:=0;
end;
{---------------------------------------------------------------------------}

{ posunieme cislo a exponent o step doprava }
{ pozor moze zmiznut cislo }
procedure bignumber.ExponentRight(step:integer);
var i:longint;
begin
  if(step=0)then
     exit;

  { overflow }
  if(_bn.exponent > (BN_MAX_EXPONENT - step))then
     _bn.overflow := BNE_MAXIMAL_OVERFLOW;

  if(step > BN_MAX_DIGIT)then
     step := BN_MAX_DIGIT;

  { exponent change }
  _bn.exponent := _bn.exponent+step;

  { zero? - nothing do }
  if( _bn.zero )then
     exit;

  for i:=BN_MAX_DIGIT downto step+1 do
    _bn.mantisa[i] := _bn.mantisa[i-step];

  for i:=step downto 1 do
    _bn.mantisa[i]:=0;

  { number was lost }
  if( IsZero )then begin
     _bn.overflow := BNE_MAXIMAL_OVERFLOW;;
     _bn.zero := true;
  end;
end;
{---------------------------------------------------------------------------}

{ normalize, first digit will be other like zero }
procedure bignumber.Normalize;
var i:integer;
begin
  i:=0;
  while((_bn.mantisa[i+1]=0)
    and ((i+1) < BN_MAX_DIGIT)) do
    inc(i);

  { nema zmysel aj tak je nula }
  if((i+1) = BN_MAX_DIGIT)then begin
     _bn.exponent := 0;
     _bn.zero := true;
     exit;
  end else
     _bn.zero := false;

  self.ExponentLeft(i);
end;
{---------------------------------------------------------------------------}

{ compare 2 mantises - no compare sign or exponent }
function  bignumber.compare_mantisa(bn: bignumber):integer;
var i:integer;
begin
  for i:=1 to BN_MAX_DIGIT do
  begin
    { is less }
    if(_bn.mantisa[i] < bn._bn.mantisa[i])then
    begin
       compare_mantisa := BN_MANTISA_LESS;
       exit;
    end;

    { is bigger }
    if(_bn.mantisa[i] > bn._bn.mantisa[i])then
    begin
       compare_mantisa := BN_MANTISA_HIGH;
       exit;
    end;
  end;

  { is equal }
  compare_mantisa := BN_MANTISA_EQUAL;
end;
{---------------------------------------------------------------------------}

{ minimal of integers }
function bignumber.Min(a,b:integer):integer;
begin
  if(a<b)then Min:=a
         else Min:=b;
end;
{------------------------------------------------------------------------------}

function bignumber.ToInt(s:string):longint;
var i,odp:longint;
begin
  odp:=0;
  ToInt:=odp;

  for i:=1 to length(s) do
  begin
    if(s[i] in ['0'..'9'])then
       odp:=odp*10 + ord(s[i])-ord('0')
    else
       exit;

    ToInt:=odp;
  end;
end;
{---------------------------------------------------------------------------}

function bignumber.ToStr(i:longint):string;
var s:string;
begin
  system.Str(i,s);
  ToStr:=s;
end;
{---------------------------------------------------------------------------}

{ bignumber is zero }
function bignumber.IsZero:boolean;
var i:integer;
begin
  IsZero:=false;
  { when found digit (1-9) it isn't zero }
  for i:=1 to BN_MAX_DIGIT do
    if(_bn.mantisa[i]<>0)then
       exit;

  { all digits is zero }
  IsZero := true;
end;
{---------------------------------------------------------------------------}

function bignumber.GetErrorS:string;
begin
  case _bn.overflow of
    BNE_OK:                 GetErrorS:='';
    BNE_USER_TERMINATE:     GetErrorS:='User terminate operation';
    BNE_MINIMAL_OVERFLOW:   GetErrorS:='Overflow minimal';
    BNE_MAXIMAL_OVERFLOW:   GetErrorS:='Overflow maximal';
    BNE_DIVIDE_ZERO:        GetErrorS:='Divide by Zero';
    BNE_SQRT_MINUS:         GetErrorS:='Sqrt negative number';
    BNE_EXP_MINUS:          GetErrorS:='Exp negative number';
    BNE_ARCSINUS_BIG:       GetErrorS:='ArcSin out <-1,1>';
    BNE_LN_NONPOSITIVE:     GetErrorS:='Ln nonpositive';
  else
    GetErrorS:='Number overflow';
  end;
end;
{---------------------------------------------------------------------------}

procedure bignumber.SplitOverFlow(add:byte);
begin
  if(( add <> BNE_OK )
  and( _bn.overflow = BNE_OK ))then
      _bn.overflow := add;
end;
{---------------------------------------------------------------------------}

{ convert string to number }
procedure bignumber.Val(s:string);
var i: integer;
begin
  self.Clear;
  _bn.overflow:=BNE_OK;

  { to native format }
{$IFDEF _BN_NATIVE}
  system.val(s,_bn.native,i);
{$ENDIF}

  { sign minus }
  if(s[1] = BN_CHAR_SIGN)then begin
     _bn.sign:=true;
     delete(s,1,1);
  end
  else
     _bn.sign:=false;

  { compute exponent }
  _bn.exponent := pos(BN_DECIMAL_POINT,s);
  if(_bn.exponent > 0)then begin
    delete(s, _bn.exponent, 1);
    _bn.exponent := _bn.exponent-2;
  end else
  begin
    _bn.exponent := pos('E',s);
    if(_bn.exponent > 0)then
       _bn.exponent := _bn.exponent-2
    else
       _bn.exponent := length(s)-1;
  end;

  { exponent from E }
  if(pos('E',s)>0)then begin
    i:=pos('E',s);
    if(s[i+1]='-')then
       _bn.exponent:=_bn.exponent - self.ToInt( system.copy(s,i+2,length(s)-i-1))
    else
       _bn.exponent:=_bn.exponent + self.ToInt( system.copy(s,i+2,length(s)-i-1));
    delete(s,i,length(s)-i+1);
  end;

  { fill mantisa }
  for i:=1 to min(length(s),BN_MAX_DIGIT) do
    if(s[i] in ['0'..'9'])then
       _bn.mantisa[i]:=ord(s[i])-ord('0');

  { remainder clear }
  for i:=length(s)+1 to BN_MAX_DIGIT do
    _bn.mantisa[i]:=0;

  { move when need }
  self.Normalize;
end;
{------------------------------------------------------------------------------}

{ convert number to string }
function bignumber.Str(n,d:integer):string;
var  i: integer;
  elen: integer;
  _exp: longint;
  _max: integer;
     s: string;
begin
  s:='';
  { sing }
  if(_bn.sign)then
     s := BN_CHAR_SIGN;

  elen:=2 + length(ToStr(abs(_bn.exponent)));
  _exp:=_bn.exponent+1;

  { negative exponent }
  if(_bn.exponent < 0)then begin
     { je to v norme urobime tvar 0.000}
     if( abs(_bn.exponent) <= BN_EXPONENT_ZERO )then begin
        s:=s + '0' + BN_DECIMAL_POINT;
        for i:=1 to abs(_bn.exponent)-1 do s:=s+'0';
        _exp:=0;
     end else _exp:=1;
  end;

  { all digits }
  for i:=1 to BN_MAX_STRING do begin
    s:=s + ToStr(_bn.mantisa[i]);
    dec(_exp);
    dec(n);
    dec(d);
    if(_exp = 0)then begin
       s:=s + BN_DECIMAL_POINT;
       break;
    end;

    if(n<=elen) then
       break;
  end;

  { found finish }
  _max:=BN_MAX_STRING;
  while((_bn.mantisa[_max]=0)
    and (_max>1)) do
    dec(_max);

  { remainder }
  while((d>0) and (i<_max)) do begin
    dec(d);
    inc(i);
    s:=s + self.ToStr(_bn.mantisa[i]);
  end;

  { decimal point on finale clear }
  if(pos(BN_DECIMAL_POINT,s) = length(s))then
     delete(s,length(s),1);

  { exponent need }
  if(_exp > 0)then
     s:=s+'E+'+ToStr(_exp);

  { negative exponent }
  if(_exp = 0)then
   if( _bn.exponent < (-1)*BN_EXPONENT_ZERO )then
     s:=s+'E'+ToStr(_bn.exponent);

  { final }
  Str:=s;
end;
{------------------------------------------------------------------------------}

{ write number }
procedure bignumber.Write(n,d:integer);
begin
  { write number }
  system.write( self.Str(n,d));

  { when number wrong - announce }
  if(_bn.overflow <> BNE_OK)then
     system.write(' <', self.GetErrorS, '>' );

{$IFDEF _BN_NATIVE}
  { in bracket write native number }
  system.write(' (',_bn.native,')');
{$ENDIF}
end;
{------------------------------------------------------------------------------}

{ Alone decimal part. }
procedure bignumber.Frac;
var i:integer;
begin
{$IFDEF _BN_NATIVE}
  { for test }
  _bn.native := system.Frac(_bn.native);
{$ENDIF}

  for i:=1 to self.min(_bn.Exponent+1, BN_MAX_DIGIT) do
    _bn.mantisa[i]:=0;

  { if need }
  self.Normalize;
end;
{------------------------------------------------------------------------------}

{ Cut decimal part.}
procedure bignumber.Trunc;
var i:longint;
begin
{$IFDEF _BN_NATIVE}
  { for test }
  _bn.native := system.Trunc(_bn.native);
{$ENDIF}

  { its decimal number }
  if(_bn.Exponent < 0)then begin
     self.Clear;
     exit;
  end;

  for i:=_bn.Exponent+2 to BN_MAX_DIGIT do
    _bn.mantisa[i]:=0;

  { if need }
  self.Normalize;
end;
{------------------------------------------------------------------------------}

{ added two numbers, allow different sign }
procedure bignumber.plus(add: bignumber);
var carry:byte;
    i:integer;
begin
{$IFDEF _BN_NATIVE}
  { for test }
  _bn.native:=_bn.native + add._bn.native;
{$ENDIF}

  { normalize exponents }
  if(_bn.exponent > add._bn.exponent)then
     add.ExponentRight(_bn.exponent - add._bn.exponent);

  if(add._bn.exponent > _bn.exponent)then
     self.ExponentRight(add._bn.exponent - _bn.exponent);

  { when sign are equal then added it }
  if(_bn.sign = add._bn.sign)then
  begin
      carry:=0;
      for i:=BN_MAX_DIGIT downto 1 do begin
         carry          := _bn.mantisa[i] + add._bn.mantisa[i] + carry;
         _bn.mantisa[i] := carry mod BN_NUMSYSTEM;
         carry          := carry div BN_NUMSYSTEM;
      end;

      { remainder carry }
      if(carry > 0)then begin
         self.ExponentRight(1);
         _bn.mantisa[1] := carry;
      end;
  end
  else begin
  { it will be minus }
      { but only big of less }
      if( self.compare_mantisa(add) = BN_MANTISA_LESS)then begin
          { invert sign }
          _bn.sign:=not(_bn.sign);
          carry:=1;
          for i:=BN_MAX_DIGIT downto 1 do begin
             carry          := add._bn.mantisa[i]+BN_NUMSYSTEM - _bn.mantisa[i]+ carry-1;
             _bn.mantisa[i] := carry mod BN_NUMSYSTEM;
             carry          := carry div BN_NUMSYSTEM;
          end;
      end else
      begin
          carry:=1;
          for i:=BN_MAX_DIGIT downto 1 do begin
             carry          := _bn.mantisa[i]+BN_NUMSYSTEM - add._bn.mantisa[i]+ carry-1;
             _bn.mantisa[i] := carry mod BN_NUMSYSTEM;
             carry          := carry div BN_NUMSYSTEM;
          end;
      end;
  end;

  { normalize it? }
  self.Normalize;

  { split overflow }
  self.SplitOverFlow( add._bn.overflow );
end;
{------------------------------------------------------------------------------}

{ minus replace with plus }
procedure bignumber.minus(les: bignumber);
begin
{$IFDEF _BN_NATIVE}
  les._bn.native:=(-1)*les._bn.native;
{$ENDIF}

  les._bn.sign:=not(les._bn.sign);
  self.plus(les);
end;
{------------------------------------------------------------------------------}

procedure bignumber.multiply(sub: bignumber);
var _mantisa:array[1..2*BN_MAX_DIGIT] of byte;
    carry:byte;
    i:integer;
    i2:integer;
begin
{$IFDEF _BN_NATIVE}
  { for test }
  _bn.native:=_bn.native*sub._bn.native;
{$ENDIF}

  { it's was a simple }
  if( _bn.sign = sub._bn.sign )then
       _bn.sign:=false
  else _bn.sign:=true;

  { auxiliary array clear }
  for i:=1 to 2*BN_MAX_DIGIT do
    _mantisa[i]:=0;

  { multiplying }
  for i2:=BN_MAX_DIGIT downto 1 do
  begin
   carry:=0;

   { if non-zero }
   if(sub._bn.mantisa[i2]<>0)then
       for i:=BN_MAX_DIGIT downto 1 do
       begin
         carry          := (sub._bn.mantisa[i2] * _bn.mantisa[i]) + carry + _mantisa[i2+i];
         _mantisa[i2+i] := carry mod BN_NUMSYSTEM;
         carry          := carry div BN_NUMSYSTEM;
       end;

   { remainder of carry }
   _mantisa[i2]:=_mantisa[i2] + carry;
  end;

  { found start }
  if( _mantisa[1] = 0 )then
       carry:=1
  else carry:=0;

  { test for overflow exponents }
  if( sub._bn.exponent > 0)then
   if( _bn.exponent > (BN_MAX_EXPONENT - sub._bn.exponent-1))then
     _bn.overflow := BNE_MAXIMAL_OVERFLOW;

  { exponents adds }
  _bn.exponent:= _bn.exponent + sub._bn.exponent + 1 - carry;

  { move }
  for i:=1 to BN_MAX_DIGIT do
    _bn.mantisa[i] := _mantisa[i+carry];

  { normalize if need }
  self.Normalize;

  { split overflow }
  self.SplitOverFlow( sub._bn.overflow );
end;
{------------------------------------------------------------------------------}

procedure bignumber.divide(del: bignumber);
var _mantisa:array[1..BN_MAX_DIGIT] of byte;
    carry:byte;
    i:integer;
    i2:integer;
begin
  { ATTENTION!!! divide by zero }
  if(del.IsZero)then begin
    _bn.overflow := BNE_DIVIDE_ZERO;
    exit;
  end;

{$IFDEF _BN_NATIVE}
  { for test }
  _bn.native:=_bn.native/del._bn.native;
{$ENDIF}

  { sign }
  if( _bn.sign = del._bn.sign )then
       _bn.sign:=false
  else _bn.sign:=true;

  { auxiliary array clear }
  for i:=1 to BN_MAX_DIGIT do
    _mantisa[i]:=0;

  { move to right because we need move 10 }
  self.ExponentRight(1);
  del.ExponentRight(1);

  for i2:=1 to BN_MAX_DIGIT do
  begin
      while( compare_mantisa(del) <> BN_MANTISA_LESS )do
      begin
         { increment 1 }
        inc(_mantisa[i2]);

        { minus }
        carry:=1;
        for i:=BN_MAX_DIGIT downto 1 do begin
           carry          := _bn.mantisa[i]+BN_NUMSYSTEM - del._bn.mantisa[i]+ carry-1;
           _bn.mantisa[i] := carry mod BN_NUMSYSTEM;
           carry          := carry div BN_NUMSYSTEM;
        end;
      end;

      { multiply with 10 }
      for i:=1 to BN_MAX_DIGIT-1 do
        _bn.mantisa[i]:=_bn.mantisa[i+1];

      { this is zero }
      _bn.mantisa[BN_MAX_DIGIT]:=0;
  end;

  { move vysledok }
  for i:=1 to BN_MAX_DIGIT do
    _bn.mantisa[i] := _mantisa[i];

  { exponents minus }
  _bn.exponent:= _bn.exponent - del._bn.exponent;

  { normalize if need }
  self.Normalize;

  { split overflow }
  self.SplitOverFlow( del._bn.overflow );
end;
{------------------------------------------------------------------------------}

begin
  writeln('This program use unit bignumber.');
end.

Copyrigth by Zdeno Sekerak 2007, http://www.trsek.com