AdventOfCode2023/UBigInt.pas

{
  Solutions to the Advent Of Code.
  Copyright (C) 2022-2024  Stefan Müller

  This program is free software: you can redistribute it and/or modify it under
  the terms of the GNU General Public License as published by the Free Software
  Foundation, either version 3 of the License, or (at your option) any later
  version.

  This program is distributed in the hope that it will be useful, but WITHOUT
  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
  FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more details.

  You should have received a copy of the GNU General Public License along with
  this program.  If not, see <http://www.gnu.org/licenses/>.
}

unit UBigInt;

{$mode ObjFPC}{$H+}

interface

uses
  Classes, SysUtils, Math;

type
  TDigits = array of Cardinal;

  { TBigInt }

  // This is an abbreviated reimplementation in Freepascal of a C# class created in 2022.
  TBigInt = object
  private
    FDigits: TDigits;
    FIsNegative: Boolean;

    function GetSign: Integer;

    // Copies consecutive digits from this BigInt to create a new one. The result will be positive. Leading zeros are
    // removed from the result, but AIndex + ACount must not exceed the number of digits of this BigInt.
    // AIndex is the first (least significant) digit to be taken. The digit with this index will become the 0th digit of
    // the new BigInt.
    // ACount is the number of consecutive digits to be taken, and the number of digits of the result.
    function GetSegment(const AIndex, ACount: Integer): TBigInt;

    // Compares the absolute value of this TBigInt object to the absolute value of another one. Returns -1 if this
    // object is less than AOther, 1 if this object is greater than AOther, and 0 if they are equal.
    function CompareToAbsoluteValues(constref AOther: TBigInt): Integer;

    class function GetZero: TBigInt; static;

    // Adds A and B, ignoring their signs and using ReturnNegative instead. The result is
    //   Sign * (Abs(A) + Abs(B)),
    // where Sign is 1 for ReturnNegative = False and -1 otherwise.
    class function AddAbsoluteValues(constref AA, AB: TBigInt; const AReturnNegative: Boolean): TBigInt; static;

    // Subtracts B from A, ignoring their signs. However, the result might be negative, and the sign can be reversed by
    // setting ReturnNegative to True. The result is
    //   Sign * (Abs(A) - Abs(B)),
    // where Sign is 1 for ReturnNegative = False and -1 otherwise.
    class function SubtractAbsoluteValues(constref AA, AB: TBigInt; const AReturnNegative: Boolean): TBigInt; static;

    // Multiplies A and B, ignoring their signs and using ReturnNegative instead. This multiplication uses a recursive
    // implementation of the Karatsuba algorithm. See
    // https://www.geeksforgeeks.org/karatsuba-algorithm-for-fast-multiplication-using-divide-and-conquer-algorithm/
    // The result is
    //   Sign * (Abs(a) * Abs(b))
    // where Sign is 1 for ReturnNegative = False and -1 otherwise.
    class function MultiplyAbsoluteValues(constref AA, AB: TBigInt; const AReturnNegative: Boolean): TBigInt; static;

    class function FromHexOrBinString(const AValue: string; const AFromBase: Integer): TBigInt; static;
    class function ConvertDigitBlock(const AValue: string; var AStartIndex: Integer; const ACharBlockSize, AFromBase:
      Integer): Cardinal;
  public
    property IsNegative: Boolean read FIsNegative;
    property Sign: Integer read GetSign;
    class property Zero: TBigInt read GetZero;
    function CompareTo(constref AOther: TBigInt): Integer;
    function TryToInt64(out AOutput: Int64): Boolean;
    class function FromInt64(const AValue: Int64): TBigInt; static;
    class function FromHexadecimalString(const AValue: string): TBigInt; static;
    class function FromBinaryString(const AValue: string): TBigInt; static;
  end;

  { Operators }

  operator := (const A: Int64): TBigInt;
  operator - (const A: TBigInt): TBigInt;
  operator + (const A, B: TBigInt): TBigInt;
  operator - (const A, B: TBigInt): TBigInt;
  operator * (const A, B: TBigInt): TBigInt;
  operator shl (const A: TBigInt; const B: Integer): TBigInt;
  operator = (const A, B: TBigInt): Boolean;
  operator <> (const A, B: TBigInt): Boolean;
  operator < (const A, B: TBigInt): Boolean;
  operator <= (const A, B: TBigInt): Boolean;
  operator > (const A, B: TBigInt): Boolean;
  operator >= (const A, B: TBigInt): Boolean;

implementation

const
  CBase = Cardinal.MaxValue + 1;
  CMaxDigit = Cardinal.MaxValue;
  CDigitSize = SizeOf(Cardinal);
  CBitsPerDigit = CDigitSize * 8;
  CHalfBits = CBitsPerDigit >> 1;
  CHalfDigitMax = (1 << CHalfBits) - 1;

  CZero: TBigInt = (FDigits: (0); FIsNegative: False);

{ TBigInt }

function TBigInt.GetSign: Integer;
begin
  if FIsNegative then
    Result := -1
  else if (Length(FDigits) > 1) or (FDigits[0] <> 0) then
    Result := 1
  else
    Result := 0;
end;

function TBigInt.GetSegment(const AIndex, ACount: Integer): TBigInt;
var
  trimmedCount: Integer;
begin
  trimmedCount := ACount;
  while (trimmedCount > 1) and (FDigits[AIndex + trimmedCount - 1] = 0) do
    Dec(trimmedCount);
  Result.FDigits := Copy(FDigits, AIndex, trimmedCount);
  Result.FIsNegative := False;
end;

function TBigInt.CompareToAbsoluteValues(constref AOther: TBigInt): Integer;
var
  i: Integer;
begin
  Result := Length(FDigits) - Length(AOther.FDigits);
  if Result = 0 then
  begin
    for i := High(FDigits) downto 0 do
      if FDigits[i] < AOther.FDigits[i] then
      begin
       Result := -1;
       Break;
      end
      else if FDigits[i] > AOther.FDigits[i] then
      begin
       Result := 1;
       Break;
      end;
  end;
end;

class function TBigInt.GetZero: TBigInt;
begin
  Result := CZero;
end;

class function TBigInt.AddAbsoluteValues(constref AA, AB: TBigInt; const AReturnNegative: Boolean): TBigInt;
var
  i, j, lenA, lenB, len, shorter: Integer;
  carry: Cardinal;
begin
  lenA := Length(AA.FDigits);
  lenB := Length(AB.FDigits);

  // Initializes the digits array of the result, with a simple test to try to predict a carry-over into a new digit. The
  // result could still carry into new digit depending on lower digits (carry over multiple digits), which would be
  // handled at the end.
  if lenA = lenB then
    if CMaxDigit - AA.FDigits[lenA - 1] < AB.FDigits[lenB - 1] then
      // Result will carry into new digit.
      SetLength(Result.FDigits, lenA + 1)
    else
      SetLength(Result.FDigits, lenA)
  else
    SetLength(Result.FDigits, Max(lenA, lenB));
  len := Length(Result.FDigits);

  // Calculates the new digits from less to more significant until the end of the shorter operand is reached.
  shorter := Min(Length(AA.FDigits), Length(AB.FDigits));
  i := 0;
  carry := 0;
  while i < shorter do
  begin
    if (AB.FDigits[i] = CMaxDigit) and (carry > 0) then
    begin
      Result.FDigits[i] := AA.FDigits[i];
      carry := 1;
    end
    else
      if CMaxDigit - AA.FDigits[i] < AB.FDigits[i] + carry then
      begin
        Result.FDigits[i] := AB.FDigits[i] + carry - 1 - (CMaxDigit - AA.FDigits[i]);
        carry := 1;
      end
      else begin
        Result.FDigits[i] := AA.FDigits[i] + AB.FDigits[i] + carry;
        carry := 0;
      end;
    Inc(i);
  end;

  // Copies the missing unchanged digits from the longer operand to the result, if any, before applying remaining
  // carry-overs. This avoids additional tests for finding the shorter digit array.
  if (i < lenA) or (i < lenB) then
    if lenA >= lenB then
      for j := i to len - 1 do
        Result.FDigits[j] := AA.FDigits[j]
    else
      for j := i to len - 1 do
        Result.FDigits[j] := AB.FDigits[j];

  // Applies the remaining carry-overs until the end of the prepared result digit array.
  while (carry > 0) and (i < len) do
  begin
    if Result.FDigits[i] = CMaxDigit then
      Result.FDigits[i] := 0
    else begin
      Inc(Result.FDigits[i]);
      carry := 0;
    end;
    Inc(i);
  end;

  // Applies the carry-over into a new digit that was not anticipated in the initialization at the top (carry over
  // multiple digits).
  if carry > 0 then
    Insert(1, Result.FDigits, len);
  Result.FIsNegative := AReturnNegative;
end;

class function TBigInt.SubtractAbsoluteValues(constref AA, AB: TBigInt; const AReturnNegative: Boolean): TBigInt;
var
  a, b: TBigInt;
  carry: Cardinal;
  compare, i, j, lastNonZeroDigitIndex, len: Integer;
begin
  // Establishes the operand order, such that Abs(a) is not less than Abs(b).
  compare := AA.CompareToAbsoluteValues(AB);
  if compare = 0 then
  begin
    Result := Zero;
    Exit;
  end;

  if compare > 0 then
  begin
    a := AA;
    b := AB;
    Result.FIsNegative := AReturnNegative;
  end
  else begin
    a := AB;
    b := AA;
    Result.FIsNegative := not AReturnNegative;
  end;

  // Calculates the new digits from less to more significant until the end of the shorter operand is reached and all
  // carry-overs have been applied.
  len := Length(a.FDigits);
  SetLength(Result.FDigits, len);
  carry := 0;
  // Tracks leading zeros for the trim below.
  lastNonZeroDigitIndex := 0;
  i := 0;
  while i < Length(b.FDigits) do
  begin
    if (a.FDigits[i] = b.FDigits[i]) and (carry > 0) then
    begin
      Result.FDigits[i] := CMaxDigit;
      carry := 1;
      lastNonZeroDigitIndex := i;
    end
    else begin
      if a.FDigits[i] < b.FDigits[i] then
      begin
        Result.FDigits[i] := CMaxDigit - (b.FDigits[i] - a.FDigits[i]) + 1 - carry;
        carry := 1;
      end
      else begin
        Result.FDigits[i] := a.FDigits[i] - b.FDigits[i] - carry;
        carry := 0;
      end;
      if Result.FDigits[i] > 0 then
        lastNonZeroDigitIndex := i;
    end;
    Inc(i);
  end;
  while carry > 0 do
  begin
    if a.FDigits[i] = 0 then
    begin
      Result.FDigits[i] := CMaxDigit;
      lastNonZeroDigitIndex := i;
    end
    else begin
      Result.FDigits[i] := a.FDigits[i] - carry;
      carry := 0;
      if Result.FDigits[i] > 0 then
        lastNonZeroDigitIndex := i;
    end;
    Inc(i);
  end;

  // Copies the missing unchanged digits from the longer operand to the result, if any. If there are none, then no trim
  // needs to occur because the most significant digit is not zero.
  if i < len then
    for j := i to len - 1 do
      Result.FDigits[j] := a.FDigits[j]
  else if (lastNonZeroDigitIndex + 1 < len) then
    // Trims leading zeros from the digits array.
    Delete(Result.FDigits, lastNonZeroDigitIndex + 1, len - lastNonZeroDigitIndex - 1);
end;

class function TBigInt.MultiplyAbsoluteValues(constref AA, AB: TBigInt; const AReturnNegative: Boolean): TBigInt;
var
  lenA, lenB, lenMax, floorHalfLength, ceilHalfLength: Integer;
  a1, a0, b1, b0, a1b1, a0b0: Cardinal;
  am, bm, middle, biga1, biga0, bigb1, bigb0, biga1b1, biga0b0: TBigInt;
begin
  lenA := Length(AA.FDigits);
  lenB := Length(AB.FDigits);
  if (lenA <= 1) and (lenB <= 1) then
  begin
    if (AA.FDigits[0] <= CHalfDigitMax) and (AB.FDigits[0] <= CHalfDigitMax) then
      if (AA.FDigits[0] = 0) or (AB.FDigits[0] = 0) then
        Result := Zero
      else begin
        Result.FDigits := TDigits.Create(AA.FDigits[0] * AB.FDigits[0]);
        Result.FIsNegative := AReturnNegative;
      end
    else begin
      // a1, a0, b1, b0 use only the lower (less significant) half of the bits of a digit, so the product of any two of
      // these fits in one digit.
      a1 := AA.FDigits[0] >> CHalfBits;
      a0 := AA.FDigits[0] and CHalfDigitMax;
      b1 := AB.FDigits[0] >> CHalfBits;
      b0 := AB.FDigits[0] and CHalfDigitMax;
      a1b1 := a1 * b1;
      a0b0 := a0 * b0;

      if a1b1 > 0 then
        Result.FDigits := TDigits.Create(a0b0, a1b1)
      else
        Result.FDigits := TDigits.Create(a0b0);
      Result.FIsNegative := AReturnNegative;

      // The result of (a1 + a0) * (b1 + b0) might not fit in one digit, so one last recursion step is necessary.
      am := FromInt64(a1 + a0);
      bm := FromInt64(b1 + b0);
      middle := (MultiplyAbsoluteValues(am, bm, False) - a1b1 - a0b0) << CHalfBits;
      if AReturnNegative then
        Result := Result - middle
      else
        Result := Result + middle;
    end;
  end
  else begin
    // Calculates where to split the two numbers.
    lenMax := Max(lenA, lenB);
    floorHalfLength := lenMax >> 1;
    ceilHalfLength := lenMax - floorHalfLength;

    // Performs one recursion step.
    if ceilHalfLength < lenA then
    begin
      biga1 := AA.GetSegment(ceilHalfLength, lenA - ceilHalfLength);
      biga0 := AA.GetSegment(0, ceilHalfLength);

      if ceilHalfLength < lenB then
      begin
        bigb1 := AB.GetSegment(ceilHalfLength, lenB - ceilHalfLength);
        bigb0 := AB.GetSegment(0, ceilHalfLength);
        biga1b1 := MultiplyAbsoluteValues(biga1, bigb1, AReturnNegative);
        biga0b0 := MultiplyAbsoluteValues(biga0, bigb0, AReturnNegative);
        Result := (biga1b1 << (2 * ceilHalfLength * CBitsPerDigit))
          + ((MultiplyAbsoluteValues(biga1 + biga0, bigb1 + bigb0, AReturnNegative) - biga1b1 - biga0b0) << (ceilHalfLength * CBitsPerDigit))
          + biga0b0;
      end
      else begin
        biga0b0 := MultiplyAbsoluteValues(biga0, AB, AReturnNegative);
        Result := ((MultiplyAbsoluteValues(biga1 + biga0, AB, AReturnNegative) - biga0b0) << (ceilHalfLength * CBitsPerDigit)) + biga0b0;
      end;
    end
    else begin
      bigb1 := AB.GetSegment(ceilHalfLength, lenB - ceilHalfLength);
      bigb0 := AB.GetSegment(0, ceilHalfLength);
      biga0b0 := MultiplyAbsoluteValues(AA, bigb0, AReturnNegative);
      Result := ((MultiplyAbsoluteValues(AA, bigb1 + bigb0, AReturnNegative) - biga0b0) << (ceilHalfLength * CBitsPerDigit)) + biga0b0;
    end;
  end;
end;

class function TBigInt.FromHexOrBinString(const AValue: string; const AFromBase: Integer): TBigInt;
var
  charBlockSize, offset, i, j, k, remainder: Integer;
  d: Cardinal;
begin
  // 2 ^ (32 / charBlockSize) = AFromBase
  case AFromBase of
    2: charBlockSize := 32;
    16: charBlockSize := 8;
  end;

  if AValue[1] = '-' then
  begin
    offset := 2;
    Result.FIsNegative := True;
  end
  else begin
    offset := 1;
    Result.FIsNegative := False;
  end;

  // Calculates the first (most significant) digit d of the result.
  DivMod(AValue.Length - offset + 1, charBlockSize, i, remainder);
  k := offset;
  d := 0;
  // Checks the first block of chars that is not a full block.
  if remainder > 0 then
    d := ConvertDigitBlock(AValue, k, remainder, AFromBase);
  // Checks full blocks of chars for first digit.
  while (d = 0) and (i > 0) do
  begin
    Dec(i);
    d := ConvertDigitBlock(AValue, k, charBlockSize, AFromBase);
  end;

  // Checks for zero.
  if (d = 0) and (i = 0) then
    Result := Zero
  else begin
    // Initializes the array of digits.
    SetLength(Result.FDigits, i + 1);
    Result.FDigits[i] := d;

    // Calculates the other digits.
    for j := i - 1 downto 0 do
      Result.FDigits[j] := ConvertDigitBlock(AValue, k, charBlockSize, AFromBase);
  end;
end;

class function TBigInt.ConvertDigitBlock(const AValue: string; var AStartIndex: Integer; const ACharBlockSize,
  AFromBase: Integer): Cardinal;
var
  part: string;
begin
  part := Copy(AValue, AStartIndex, ACharBlockSize);
  Inc(AStartIndex, ACharBlockSize);
  case AFromBase of
    2: part := '%' + part;
    8: part := '&' + part;
    16: part := '$' + part;
  end;
  Result := StrToDWord(part);
end;

function TBigInt.CompareTo(constref AOther: TBigInt): Integer;
begin
  if FIsNegative = AOther.FIsNegative then
    Result := CompareToAbsoluteValues(AOther)
  else
    Result := 1;
  if FIsNegative then
    Result := -Result;
end;

function TBigInt.TryToInt64(out AOutput: Int64): Boolean;
begin
  AOutput := 0;
  Result := False;
  case Length(FDigits) of
    0: Result := True;
    1: begin
      AOutput := FDigits[0];
      if FIsNegative then
        AOutput := -AOutput;
      Result := True;
    end;
    2: begin
      if FDigits[1] <= Integer.MaxValue then
      begin
        AOutput := FDigits[1] * CBase + FDigits[0];
        if FIsNegative then
          AOutput := -AOutput;
        Result := True;
      end
      else if (FDigits[1] = Integer.MaxValue + 1) and (FDigits[0] = 0) and FIsNegative then
      begin
        AOutput := Int64.MinValue;
        Result := True;
      end;
    end;
  end;
end;

class function TBigInt.FromInt64(const AValue: Int64): TBigInt;
var
  absVal: Int64;
begin
  if AValue <> Int64.MinValue then
  begin
    absVal := Abs(AValue);
    if absVal >= CBase then
      Result.FDigits := TDigits.Create(absVal mod CBase, absVal div CBase)
    else
      Result.FDigits := TDigits.Create(absVal);
    Result.FIsNegative := AValue < 0;
  end
  else begin
    Result.FDigits := TDigits.Create(0, 1 << 31);
    Result.FIsNegative := True;
  end;
end;

class function TBigInt.FromHexadecimalString(const AValue: string): TBigInt;
begin
  Result := FromHexOrBinString(AValue, 16);
end;

class function TBigInt.FromBinaryString(const AValue: string): TBigInt;
begin
  Result := FromHexOrBinString(AValue, 2);
end;

{ Operators }

operator := (const A: Int64): TBigInt;
begin
  Result := TBigInt.FromInt64(A);
end;

operator - (const A: TBigInt): TBigInt;
var
  len: Integer;
begin
  len := Length(A.FDigits);
  if (len > 1) or (A.FDigits[0] > 0) then
  begin
    Result.FDigits := Copy(A.FDigits, 0, len);
    Result.FIsNegative := not A.FIsNegative;
  end
  else
    Result := TBigInt.Zero;
end;

operator + (const A, B: TBigInt): TBigInt;
begin
  if A.IsNegative = B.IsNegative then
    Result := TBigInt.AddAbsoluteValues(A, B, A.IsNegative)
  else
    Result := TBigInt.SubtractAbsoluteValues(A, B, A.IsNegative);
end;

operator - (const A, B: TBigInt): TBigInt;
begin
  if A.IsNegative = B.IsNegative then
    Result := TBigInt.SubtractAbsoluteValues(A, B, A.IsNegative)
  else
    Result := TBigInt.AddAbsoluteValues(A, B, A.IsNegative);
end;

operator * (const A, B: TBigInt): TBigInt;
begin
  if  (A = TBigInt.Zero) or (B = TBigInt.Zero) then
    Result := TBigInt.Zero
  else
    Result := TBigInt.MultiplyAbsoluteValues(A, B, A.IsNegative <> B.IsNegative);
end;

// TODO: Shift operator could be implemented with a single Move call, but I do not want to change it without test cases.
operator shl(const A: TBigInt; const B: Integer): TBigInt;
var
  i, j, digitShifts, bitShifts, reverseShift, len, newLength: Integer;
  lastDigit: Cardinal;
begin
  // Handles shift of zero.
  if A = 0 then
    Result := TBigInt.Zero
  else begin
    // Determines full digit shifts and bit shifts.
    DivMod(B, CBitsPerDigit, digitShifts, bitShifts);

    if bitShifts > 0 then
    begin
      reverseShift := CBitsPerDigit - bitShifts;
      len := Length(A.FDigits);
      lastDigit := A.FDigits[len - 1] >> reverseShift;
      newLength := len + digitShifts;

      if lastDigit = 0 then
        SetLength(Result.FDigits, newLength)
      else
        SetLength(Result.FDigits, newLength + 1);

      // Performs full digit shifts by shifting the access index j for A.FDigits.
      Result.FDigits[digitShifts] := A.FDigits[0] << bitShifts;
      j := 0;
      for i := digitShifts + 1 to newLength - 1 do
      begin
        // Performs bit shifts.
        Result.FDigits[i] := A.FDigits[j] >> reverseShift;
        Inc(j);
        Result.FDigits[i] := Result.FDigits[i] or (A.FDigits[j] << bitShifts);
      end;

      if Length(Result.FDigits) > newLength then
        Result.FDigits[newLength] := lastDigit;
    end
    else begin
      // Performs full digit shifts by copy if there are no bit shifts.
      len := Length(A.FDigits);
      SetLength(Result.FDigits, len + digitShifts);
      for i := 0 to digitShifts - 1 do
        Result.FDigits[i] := 0;
      for i := 0 to len - 1 do
        Result.FDigits[i + digitShifts] := A.FDigits[i];
    end;

    Result.FIsNegative := A.IsNegative;
  end;
end;

operator = (const A, B: TBigInt): Boolean;
begin
  Result := A.CompareTo(B) = 0;
end;

operator <> (const A, B: TBigInt): Boolean;
begin
  Result := A.CompareTo(B) <> 0;
end;

operator < (const A, B: TBigInt): Boolean;
begin
  Result := A.CompareTo(B) < 0;
end;

operator <= (const A, B: TBigInt): Boolean;
begin
  Result := A.CompareTo(B) <= 0;
end;

operator > (const A, B: TBigInt): Boolean;
begin
  Result := A.CompareTo(B) > 0;
end;

operator >= (const A, B: TBigInt): Boolean;
begin
  Result := A.CompareTo(B) >= 0;
end;

end.