Code-Beispiel
Einfache und schnelle Bigint-Erweiterung
| Lizenz: | Erster Autor: | Letzte Bearbeitung: |
| k. A. |  stephanbrunker |  29.07.2018 |
Wer mit Ganzzahlen > 64 bit arbeiten will, kommt um eine Bigint / Largeint -Erweiterung nicht herum. Die im Package enthaltene big_int ist leider verwaist, nur als complierte dll vorhanden und deshalb schwer zu verstehen. Alternativ gibt es aus dem Math-Package die largeint.bas, die aber in der Anwendung auch alles andere als einfach ist. Im freebasic.net - Forum bin ich auf die Big_Integer von Richard gestoßen, die dank eines eigenen Datentyps und überladener Operatoren ganz einfach einzubinden ist. Leider war dieses Programm sehr viel langsamer als die largeint.bas. Das konnte ich lösen, indem ich das ganze Programm optimiert habe, dazu musste ich auch die Division und die Komparatoren komplett neu schreiben. Dann habe ich die Struktur / den Typ komplett überarbeitet. Die Funktionen sind jetzt alle Member des Typs, alle Datentypen sind vollständig überladen.
Das bedeutet, man kann jetzt den Datentyp BIGINT genauso verwenden wie einen INTEGER. Lediglich die Funktion CAST (bigint, (irgendwas) ) geht nicht, weil man dazu den CAST-Teil von irgendwas ändern müsste. Darum funktionieren auch die Makros Bit, BitSet und BitReset nicht, weil diese intern einen CAST verwenden.
Eine kurze Anleitung steht am Anfang des Codes. Die gebräuchlichsten Formate zur Konvertierung in einen Bigint sind naturgemäß Strings, und zwar entweder dezimal, hexadezimal oder als binär mit den gleichen Konvertierungsfunktionen wie sonst auch (Bin, Hex, CV, MK ...).
Mit der Ursprungsversion dauerte die Suche nach einer 1024-bit Sophie-Germain-Primzahl drei Tage auf fünf Rechnern, jetzt bekommt das ein Rechner in einer Stunde hin. Die Geschwindigkeit ist aber nicht vordergründig, das können in anderen Sprachen programmierte Programme / Bibliotheken besser hin, FB ist sicherlich nicht für Geschwindigkeit optimiert. Dafür aber für einfach geschriebenen, leicht verständlichen Code und Bigint ist jetzt zu benutzen wie ein Integer, einfacher kann es nicht mehr werden.
Hier ist der Code am Stück, 2000 Zeilen. Weil das etwas unübersichtlich ist und weil ich die weitere Mitarbeit vereinfachen will, habe ich den Code aufgeteilt und ein Github-Repository erstellt: 
https://github.com/stephanbrunker/big_integer Außerdem habe ich eine ausführlichere Readme erstellt und ein Beispiel programmiert (Miller-Rabin-Primzahltest).
' version 2.8.1 29 July 2018
'================================================================
' by Stephan Brunker (stephanbrunker at web punkt de)
' based on version 1.2 by Richard @ freebasic.net/forum
' with credits to: Hartmut Ruske and frisian
'================================================================
' The easy way to use - analog to the "hello world" the 1+1:
'------------------------------
' #include "big_integer.bas"
' Dim as Bigint a,b,c,d
' a = 1234
' b = "1357239875892745982784975230987590287"
' c = VALBigint("&h3a84458e9bc47")
' d = a + b * c
' print d
'------------------------------
' hex and string input can be so long as you want - up to 4.29E09 bytes.
' input can be any data type (ulong, ulongint)
' also strings coded in decimal, hex or binary(endian reversed)
' all the operators are overloaded
' also all the comparators: = < > <> >= <=
' for additional functions look in the code, it should be
' well commented to see for what the functions are for
' you can also input and output in various formats, see section
' "Conversion Functions"
'===============================================================
' This package only handles integers encoded in a two's complement format.
' The first bit in the two's complement format is always the sign which
' takes the value of 1 = negative, 0 = positive or zero. Byte examples are;
' +127 = 0111 1111 most positive
' +8 = 0000 1000
' +7 = 0000 0111
' +4 = 0000 0100
' +2 = 0000 0010
' +1 = 0000 0001
' 0 = 0000 0000 zero
' -1 = 1111 1111
' -2 = 1111 1110
' -4 = 1111 1100
' -7 = 1111 1001
' -8 = 1111 1000
' -127 = 1000 0001
' -128 = 1000 0000 most negative
'----------------------------------------------------------------
' Each successive 32 bits are packed into a 4 byte block.
' Each block is stored in 4 successive bytes of a string.
'----------------------------------------------------------------
' Strings in FreeBASIC are indexed and printed from left to right. Big
' integers appear to be stored in strings backwards which can be confusing.
' The Left side of the string is the Right side of the number and vice versa.
' Beware: Shift_Left moves a string to the right, Shift_Right moves it left.
'----------------------------------------------------------------
' The first byte in the string is the least significant byte of the number.
' The last block in the string is the most significant block of the number.
' String s indexing has s[0] as the LS byte and s[len(s)-1] as the MS byte.
' These big integer strings are always multiples of 4 bytes long.
' The msb of a number is sign extended to the MS bit of the MS block.
'----------------------------------------------------------------
' A number is always stored in a way that correctly represents that number.
' Where an operation would overflow into the MSB, a sign block is
' appended to the number so as to prevent overflow or sign change.
' Unnecessary leading zero or all ones blocks are removed by prune.
'----------------------------------------------------------------
' String pointers may change if a string is created or any length is changed.
'================================================================
Type Bigint ' a little endian, two's complement binary number
s As String ' packed into a string, in blocks four bytes long
' Constructors
'-------------
Declare Constructor () ' default constructor
Declare Constructor (ByRef a As Bigint) ' copy constructor
Declare Constructor (ByRef a As Byte)
Declare Constructor (ByRef a As UByte)
Declare Constructor (ByRef a As Short)
Declare Constructor (ByRef a As UShort)
Declare Constructor (ByRef a As Long)
Declare Constructor (ByRef a As ULong)
Declare Constructor (ByRef a As Integer)
Declare Constructor (ByRef a As UInteger)
Declare Constructor (ByRef a As LongInt)
Declare Constructor (ByRef a As ULongInt)
Declare Constructor (ByRef a As Single)
Declare Constructor (ByRef a As Double)
Declare Constructor (ByRef a As String)
' let
'----
Declare Operator Let (ByRef a As Bigint)
' cast to other datatypes
'------------------------
Declare Operator Cast() As Byte ' CByte
Declare Operator Cast() As UByte ' CUByte
Declare Operator Cast() As Short ' CShort
Declare Operator Cast() As UShort ' CUShort
Declare Operator Cast() As Long ' CLng
Declare Operator Cast() As ULong ' CULng
Declare Operator Cast() As LongInt ' CLngint
Declare Operator Cast() As ULongInt ' CULngint
Declare Operator Cast() As Integer ' Cint
Declare Operator Cast() As UInteger ' CUInt
Declare Operator Cast() As Single ' CSng
Declare Operator Cast() As Double ' CDbl
Declare Operator Cast() As String ' Str
' functions
'----------
' private: (cannot make them really protected because the operators access them)
'using them directly on your own risk, because everything is passed ByRef and changed internally
Declare Static Function compare (ByRef a As Bigint, ByRef b As Bigint) As Long
Declare Static Function square(ByRef aa As Bigint) As Bigint
Declare Static Function mul2(ByRef a As Bigint) As Bigint
Declare Static Function div2(ByRef a As Bigint) As Bigint
Declare Static Sub div(ByRef aa As Bigint, ByRef bb As Bigint,ByRef q As Bigint, ByRef r As Bigint)
Declare Static Sub prune(ByRef a As Bigint)
' public:
Declare Static Function modpower(ByRef bb As Bigint, ByRef ee As Bigint, ByRef m As Bigint) As Bigint
Declare Static Function factorial (ByRef a As Bigint) As Bigint
Declare Static Function msbit(ByRef a As Bigint) As Long
' Overloading Bit / Bitset / BitReset is not possible, because these are macros
' containing a "cast(bigint,1)" operator, and you can only overload a cast from bigint to other
' and not to the bigint itself (this had to be done in the type integer for example)
' Couriously, #undef bit and declare a new Function bit(a as bigint, b as Ulongint)
' would work in this file, but not if the file is included in another project
Declare Static Function Bit_Value(ByRef v As Bigint, ByVal b As ULongInt) As Long
Declare Static Function Bit_Set(ByRef vv As Bigint, ByVal b As ULongInt) As Bigint
Declare Static Function Bit_Reset(ByRef vv As Bigint, ByVal b As ULongInt) As Bigint
' implicit step versions
Declare Operator For ()
Declare Operator Step()
Declare Operator Next(ByRef end_cond As Bigint) As Integer
' explicit step versions
Declare Operator For (ByRef step_var As Bigint)
Declare Operator Step(ByRef step_var As Bigint)
Declare Operator Next(ByRef end_cond As Bigint, ByRef step_var As Bigint ) As Integer
' operate and assign
Declare Operator += (ByRef rhs As Bigint)
Declare Operator -= (ByRef rhs As Bigint)
Declare Operator *= (ByRef rhs As Bigint)
Declare Operator \= (ByRef rhs As Bigint)
Declare Operator Mod= (ByRef rhs As Bigint)
Declare Operator ^= (ByRef rhs As Bigint)
Declare Operator Shl= (ByRef rhs As LongInt)
Declare Operator Shr= (ByRef rhs As LongInt)
Declare Operator And= (ByRef rhs As Bigint)
Declare Operator Or= (ByRef rhs As Bigint)
Declare Operator Xor= (ByRef rhs As Bigint)
Declare Operator Imp= (ByRef rhs As Bigint)
Declare Operator Eqv= (ByRef rhs As Bigint)
End Type
'================================================================
' OVERLOADED OPERATORS:
'================================================================
' Comparators
'-----------
Declare Operator < (ByRef As Bigint, ByRef As Bigint) As Integer
Declare Operator > (ByRef As Bigint, ByRef As Bigint) As Integer
Declare Operator = (ByRef As Bigint, ByRef As Bigint) As Integer
Declare Operator <> (ByRef a As Bigint, ByRef b As Bigint) As Integer
Declare Operator <= (ByRef a As Bigint, ByRef b As Bigint) As Integer
Declare Operator >= (ByRef a As Bigint, ByRef b As Bigint) As Integer
' Mathematical Operators
'-----------------------
Declare Operator + (ByRef x As Bigint) As Bigint
Declare Operator - (ByRef x As Bigint) As Bigint
Declare Operator + (ByRef x As Bigint, ByRef y As Bigint) As Bigint
Declare Operator - (ByRef x As Bigint, ByRef y As Bigint) As Bigint
Declare Operator * (ByRef x As Bigint, ByRef y As Bigint) As Bigint
Declare Operator \ (ByRef x As Bigint, ByRef y As Bigint) As Bigint
Declare Operator Mod (ByRef x As Bigint, ByRef y As Bigint) As Bigint
Declare Operator ^ (ByRef x As Bigint, ByRef y As LongInt) As Bigint
Declare Operator Abs (Byref x As Bigint) As Bigint
Declare Operator Sgn (Byref x As Bigint) As Integer
' Bitwise Operations
'-------------------
Declare Operator Not (ByRef x As Bigint) As Bigint
Declare Operator And (ByRef x As Bigint, ByRef y As Bigint) As Bigint
Declare Operator Or (ByRef x As Bigint, ByRef y As Bigint) As Bigint
Declare Operator Xor (ByRef x As Bigint, ByRef y As Bigint) As Bigint
Declare Operator Imp (ByRef x As Bigint, ByRef y As Bigint) As Bigint
Declare Operator Eqv (ByRef x As Bigint, ByRef y As Bigint) As Bigint
Declare Operator Shl (Byref a As Bigint, ByVal n As LongInt) As Bigint
Declare Operator Shr (Byref a As Bigint, ByVal n As LongInt) As Bigint
' Conversion Functions
'---------------------
Declare Function CBig Overload(a as Byte) as Bigint
Declare Function CBig Overload(a as UByte) as Bigint
Declare Function CBig Overload(a as Short) as Bigint
Declare Function CBig Overload(a as UShort) as Bigint
Declare Function CBig Overload(a as Integer) as Bigint
Declare Function CBig Overload(a as UInteger) as Bigint
Declare Function CBig Overload(a as Long) as Bigint
Declare Function CBig Overload(a as ULong) as Bigint
Declare Function CBig Overload(a as LongInt) as Bigint
Declare Function CBig Overload(a as ULongInt) as Bigint
Declare Function CBig Overload(a as Single) as Bigint
Declare Function CBig Overload(a as Double) as Bigint
Declare Function CBig Overload(a as String) as Bigint
Declare Function Bin(ByRef s As Bigint) As String
Declare Function Hex(ByRef s As Bigint) As String
Declare Function Hex(ByRef s As Bigint, ByRef n As Ulong) As String
Declare Function Uhex(ByRef s As Bigint) As String
Declare Function UHexT(ByRef s As Bigint,ByRef n As ULong) As String
Declare Function Oct(ByRef s As Bigint) As String
Declare Function MkBigint(ByRef a As Bigint) As String
Declare Function MkUBigint(ByRef a As Bigint) As String
Declare Function ValBigint(ByRef a As String) As Bigint
Declare Function ValUBigint(ByRef a As String) As Bigint
Declare Function CVBigint(ByRef aa As String) As Bigint
Declare Function CVUBigint(ByRef aa As String) As Bigint
'================================================================
' CONSTANTS
'================================================================
' pre-defining the strings for certain values is much faster
' then generating them on the fly, also operations with them
Dim Shared As String Bigint_s0
Bigint_s0 = Chr(0,0,0,0)
Dim Shared As String Bigint_s00
Bigint_s00 = Chr(0,0,0,0,0,0,0,0)
Dim Shared As String Bigint_s1
Bigint_s1 = Chr(1,0,0,0)
Dim Shared As String Bigint_s2
Bigint_s2 = Chr(2,0,0,0)
Dim Shared As String Bigint_s_1
Bigint_s_1 = Chr(-1,-1,-1,-1)
'================================================================
' CONSTRUCTORS
'================================================================
Constructor Bigint () ' default constructor
this.s = Bigint_s0 ' zero
End Constructor
Constructor Bigint (ByRef a As Bigint) ' copy constructor
this.s = a.s
End Constructor
Constructor Bigint (ByRef a As Byte)
If (128 And a) Then
this.s = Chr(a,-1,-1,-1)
Else
this.s = Chr(a,0,0,0)
End If
End Constructor
Constructor Bigint (ByRef a As UByte)
this.s = Chr(a,0,0,0)
End Constructor
Constructor Bigint (ByRef a As Short)
If (32768 And a) Then
this.s = Chr(LoByte(a), HiByte(a), -1, -1 )
Else
this.s = Chr(LoByte(a), HiByte(a), 0, 0 )
End If
End Constructor
Constructor Bigint (ByRef a As UShort)
this.s = Chr(LoByte(a), HiByte(a), 0, 0 )
End Constructor
Constructor Bigint (ByRef a As Long)
this.s = Bigint_s0
Dim As Long Ptr bip = CPtr(Long Ptr, StrPtr(this.s))
Dim As Long Ptr ip = CPtr(Long Ptr, @a)
*bip = *ip
End Constructor
Constructor Bigint (ByRef a As ULong)
this.s = Bigint_s0
Dim As ULong Ptr bip = CPtr(ULong Ptr, StrPtr(this.s))
Dim As ULong Ptr uip = CPtr(ULong Ptr, @a)
*bip = *uip
If (128 And this.s[3]) Then this.s &= Bigint_s0 ' make it positive
End Constructor
Constructor Bigint (ByRef a As Integer)
If a < -2147483648 Or a > 2147483647 Then ' integer<64>
this.s = Bigint_s00
Dim As LongInt Ptr bip = CPtr(LongInt Ptr, StrPtr(this.s))
Dim As LongInt Ptr lip = CPtr(LongInt Ptr, @a)
*bip = *lip
Else ' integer<32>
this.s = Bigint_s0
Dim As Long Ptr bip = CPtr(Long Ptr, StrPtr(this.s))
Dim As Long Ptr ip = CPtr(Long Ptr, @a)
*bip = *ip
End If
End Constructor
Constructor Bigint (ByRef a As UInteger)
If a > 4294967295 Then ' uinteger<64>
this.s = Bigint_s00
Dim As ULongInt Ptr bip = CPtr(ULongInt Ptr, StrPtr(this.s))
Dim As ULongInt Ptr ulip = CPtr(ULongInt Ptr, @a)
*bip = *ulip
If (128 And this.s[7]) Then this.s &= Bigint_s0 ' make it positive
Else ' integer<32>
this.s = Bigint_s0
Dim As ULong Ptr bip = CPtr(ULong Ptr, StrPtr(this.s))
Dim As ULong Ptr uip = CPtr(ULong Ptr, @a)
*bip = *uip
If (128 And this.s[3]) Then this.s &= Bigint_s0 ' make it positive
End If
End Constructor
Constructor Bigint (ByRef a As LongInt)
If a < -2147483648 Or a > 2147483647 Then '8Byte-Bigint
this.s = Bigint_s00
Dim As LongInt Ptr bip = CPtr(LongInt Ptr, StrPtr(this.s))
Dim As LongInt Ptr lip = CPtr(LongInt Ptr, @a)
*bip = *lip
Else '4Byte-Bigint
this.s = Bigint_s0
Dim As Long Ptr bip = CPtr(Long Ptr, StrPtr(this.s))
Dim As Long Ptr ip = CPtr(Long Ptr, @a)
*bip = *ip
End If
End Constructor
Constructor Bigint (ByRef a As ULongInt)
If a > 4294967295 Then ' 8Byte-Bigint
this.s = Bigint_s00
Dim As ULongInt Ptr bip = CPtr(ULongInt Ptr, StrPtr(this.s))
Dim As ULongInt Ptr ulip = CPtr(ULongInt Ptr, @a)
*bip = *ulip
If (128 And this.s[7]) Then this.s &= Bigint_s0 ' make it positive
Else ' 4Byte-Bigint
this.s = Bigint_s0
Dim As ULong Ptr bip = CPtr(ULong Ptr, StrPtr(this.s))
Dim As ULong Ptr uip = CPtr(ULong Ptr, @a)
*bip = *uip
If (128 And this.s[3]) Then this.s &= Bigint_s0 ' make it positive
End If
End Constructor
Constructor Bigint (ByRef a As Single)
Const As ULong implicit_bit = 2^23
Const As ULong mant_mask = implicit_bit - 1
Dim As ULong u, mant
Dim As Long negative, expo
Dim As Bigint x
'----------------------------------------------------
If a < 0 Then negative = -1 ' remember sign
a = Int(Abs(a) + 0.5) ' rectify and round to closest integer
'----------------------------------------------------
u = *(CPtr(ULong Ptr, @a)) ' copy Single into a Ulong frame
expo = (u Shr 23) And 255 ' the 8 bit exponent
mant = (u And mant_mask ) ' 23 mantissa bits
'----------------------------------------------------
If expo = 0 Then ' the double has zero value or is de-normalized
this.s = Bigint_s0 ' de-normalized is very very small
Else
mant = mant + implicit_bit ' insert the missing implicit bit
expo = expo - 127 ' remove the exponent bias
If expo < 24 Then
mant = mant Shr (23 - expo) ' reduce it to integer
x.s = Bigint_s0
*(CPtr(ULong Ptr, StrPtr(x.s))) = mant ' make Bigint
If negative Then x = - x
Else
Print "WARNING: Single argument was unreliable."
Sleep : End
End If
This = x
End If
End Constructor
Constructor Bigint (ByRef a As Double)
Const As ULongInt implicit_bit = 2^52 ' 4503599627370496
Const As ULongInt mant_mask = implicit_bit - 1 ' 4503599627370495
Dim As ULongInt u, mant
Dim As Long negative, expo
Dim As Bigint x
'----------------------------------------------------
If a < 0 Then negative = -1 ' remember sign
a = Int(Abs(a) + 0.5) ' rectify and round to closest integer
'----------------------------------------------------
u = *(CPtr(ULongInt Ptr, @a)) ' copy Double into a Ulongint frame
expo = (u Shr 52) And 2047 ' the 11 bit exponent
mant = (u And mant_mask ) ' 52 mantissa bits
'----------------------------------------------------
If expo = 0 Then ' the double has zero value or is de-normalized
this.s = Bigint_s0 ' de-normalized is very very small
Else
mant = mant + implicit_bit ' insert the missing implicit bit
expo = expo - 1023 ' remove the exponent bias
If expo < 53 Then
mant = mant Shr (52 - expo) ' reduce it to integer
x.s = Bigint_s00
*(CPtr(ULongInt Ptr, StrPtr(x.s))) = mant ' make Bigint
If negative Then x = - x
Else
Print "WARNING: Double argument was unreliable."
Sleep : End
End If
This = x
End If
End Constructor
' pack ascii numeric, octal, dual or hexadecimal string to a Bigint
' the output Bigint will have a length that is a multiple of 4 bytes
' works equivalent to VAL: ignore leading space, accept + and - at
' the beginning and stop if a non-number appears in the string
' 'u' enforces Unsigned conversion
Constructor Bigint (ByRef aa As String)
Dim as String a = LTrim(aa) ' remove space
If Len(a) = 0 Then This.s = Bigint_s0 : Exit Constructor
#Define decimal 1
#Define hexadecimal 2
#Define dual 3
#Define octal 4
Dim as Byte mode, sign, usgn
'------------------------------------------------------------
' check first character and trim to numbers only
Select case a[0]
Case Asc("&")
If Len(a) > 1 Then
Select case a[1]
Case Asc("h"), Asc("H")
mode = hexadecimal
a = Right(a,Len(a)-2)
Case Asc("b"), Asc("B")
mode = dual
a = Right(a,Len(a)-2)
Case Asc("o"), Asc("O")
mode = octal
a = Right(a,Len(a)-2)
Case Else
If a[1] > 47 And a[1] < 58 Then ' octal is &3333 and &o3333
mode = octal
a = LTrim(a,"&")
Else
This.s = Bigint_s0 : Exit Constructor
End If
End Select
Else
This.s = Bigint_s0 : Exit Constructor
End If
Case Asc("-")
a = LTrim(a,"-")
sign = -1
mode = decimal
Case Asc("+")
a = LTrim(a,"+")
mode = decimal
Case Else
If a[0] > 47 And a[0] < 58 Then
mode = decimal
Else
This.s = Bigint_s0 : Exit Constructor
End If
End Select
'------------------------------------------------------------
' check for suffix
If a[Len(a)-1] = Asc("u") Then
usgn = -1
a = RTrim(a,"u")
End If
'------------------------------------------------------------
' stop conversion if a nonconform type occurs
Dim i as Long
For i = 0 to Len(a)-1
Select Case mode
Case hexadecimal
If (a[i] < 48 And a[i] > 57 ) And (a[i] < 65 And a[i] > 70) And (a[i] < 97 And a[i] > 102) Then
a = Left(a,i)
Exit For
End If
Case octal
If a[i] < 48 And a[i] > 56 Then
a = Left(a,i)
Exit For
End If
Case dual
If a[i] < 48 And a[i] > 49 Then
a = Left(a,i)
Exit For
End If
Case decimal
If a[i] < 48 And a[i] > 57 Then
a = Left(a,i)
Exit For
End If
End Select
Next i
Select Case Mode
Case hexadecimal
' pad leading zeroes for positive hex
If (Len(a) Mod 8) <> 0 Then a = String(8-(Len(a) Mod 8),"0") & a
Dim as Ulong b
Dim as Long d
Dim as Bigint c
' get positive or negative value for the first block
If usgn = -1 Then ' unsigned suffix
b = ValUInt("&h" & Mid(a,1,8))
c = b
Else
d = ValInt("&h" & Mid(a,1,8)) ' can be positive or negative
c = d
End If
' add all following blocks (freakish, but correct also for negative values)
For i = 9 To Len(a) Step 8
c = c Shl 32
b = ValUInt("&h" & Mid(a,i,8))
c += b
Next i
This = c
Case octal
' pad leading zeroes
If (Len(a) Mod 10) <> 0 Then a = String(10-(Len(a) Mod 10),"0") & a
Dim as Ulong b
Dim as Bigint c
For i = 1 To Len(a) Step 10
b = ValUInt("&o" & Mid(a,i,10))
c = c Shl 30
c += b
Next i
This = c
Case dual
' pad leading zeroes if blocksize <> 32bit - positive
If (Len(a) Mod 32) <> 0 Then a = String(32-(Len(a) Mod 32),"0") & a
Dim as Ulong b
Dim as Long d
Dim as Bigint c
If usgn = -1 Then
b = ValUInt("&b" & Mid(a,1,32))
c = b
Else
d = ValInt("&b" & Mid(a,1,32))
c = d
End If
For i = 33 To Len(a) Step 32
c = c Shl 32
b = ValUInt("&b" & Mid(a,i,32))
c += b
Next i
This = c
Case decimal
Dim As Long p, j, blocks
' extend to next multiple of 9 digits
i = Len(a)
blocks = i \ 9 ' number of 9 digit blocks needed
If i Mod 9 <> 0 Then blocks += 1
p = 9 * blocks
a = String(p - i, "0") & a ' pad to next multiple of 9 digits
'------------------------------------------------------------
' decimal to binary conversion
i = ( 8 + Len(a) * 3.32192809488) \ 8 ' bytes needed for binary
blocks = 1 + (i \ 4) ' adjust to multiple of 4
this.s = String(blocks * 4, 0 ) ' binary destination string
'------------------------------------------------------------
Dim As ULong Ptr bp, bpz, bpcarry, bpdata
bpz = Cast(ULong Ptr, StrPtr(this.s)) ' binary output string[0]
Dim As ULongInt product, carry, multiplier = 1e9
bpdata = Cast(ULong Ptr, @product) ' bottom half of product
bpcarry = bpdata + 1 ' top half of product
'------------------------------------------------------------
blocks = 1 ' blocks will be advanced as required by carry
For i = 1 To Len(a)-8 Step 9 ' msd to lsd in blocks of 9
bp = bpz ' point back to the start
carry = ValULng(Mid(a, i, 9)) ' take the next 9 digit block
For j = 1 To blocks
product = multiplier * *bp + carry
*bp = CULng(*bpdata)
carry = CULngInt(*bpcarry)
bp += 1
Next j
' advancing blocks only as needed doubles the speed of conversion
If Carry Then
*bp = carry
blocks += 1 ' an exact count of the blocks used
End If
Next i
this.s = Left(this.s, blocks * 4) ' keep only used blocks
'-------------------------------------------------------------
If this.s[Len(this.s)-1] And 128 Then this.s &= Bigint_s0 ' MSB must be 0
If sign Then
This = - This
' That one doesn't make sense: negative Input to unsigned:
If usgn Then this.s &= Bigint_s0
End If
End Select
End Constructor
Operator Bigint.let (ByRef a As Bigint)
this.s=a.s
End Operator
'=================================================================
' OVERLOADED OPERATORS AND FUNCTIONS
'=================================================================
' unary plus +
'================================================================
Operator + (ByRef x As Bigint) As Bigint
Return x
End Operator
' unary minus -
'================================================================
Operator - (ByRef a As Bigint) As Bigint
' Negate the twos complement binary number in a Bigint
Dim As Bigint s = a
Dim As Long blocks = Len(s.s) \ 4
Dim As ULongInt sum
Dim As ULong carry
Dim As ULong Ptr ps
ps = Cast( ULong Ptr, StrPtr(s.s))' the Uinteger data in Bigint
carry = 1 ' set carry
Do ' slow ahead until clear of the carry
*ps = Not *ps
sum = CULngInt(*ps) + carry
*ps = sum
carry = sum Shr 32
ps += 1
blocks -= 1
Loop Until (carry = 0) Or (blocks = 0)
If blocks > 0 Then
Do ' no carry, so full speed ahead
*ps = Not *ps
ps +=1
blocks -= 1
Loop Until blocks = 0
End If
' Negating the most negative integer is a problem because carry propagation
' flips the sign which should have become positive. But negation of zero
' does not change the sign either, so we need to differentiate between zero
' and one by quickly examining the final carry bit from the two's complement.
If carry = 0 Then ' carry was not generated by the most negative number
If (128 And a.s[Len(a.s)-1]) = (128 And s.s[Len(s.s)-1]) Then s.s &= Bigint_s0
End If ' this prevents a negated zero being extended by an extra zero block
Return s
End Operator
' addition
'================================================================
Operator + (ByRef aa As Bigint, ByRef bb As Bigint) As Bigint
If aa = bb Then
Return Bigint.mul2(aa)
ElseIf aa.s = Bigint_s0 Then
Return bb
ElseIf bb.s = Bigint_s0 Then
Return aa
End If
Dim As Bigint a = aa, b = bb
Dim As Long blocks, i, j, lena, lenb, sa, sb, delta
'------------------------------------------------------------
' test to see if the two most significant digits differ which
lena = Len(a.s) - 1 ' might change the sign without carry
If a.s[lena] And 128 Then sa = 255 ' sign as a byte
i = a.s[lena] And 192 ' if MSBs differ then extend the sign
If (i = 64) Or (i = 128) Then a.s = a.s + String(4, Chr(sa) )
'------------------------------------------------------------
lenb = Len(b.s) - 1
If b.s[lenb] And 128 Then sb = 255
i = b.s[lenb] And 192
If (i = 64) Or (i = 128) Then b.s = b.s + String(4, Chr(sb) )
'------------------------------------------------------------
' align the two Bigints to be added
delta = Len(a.s) - Len(b.s) ' new values
If delta <> 0 Then ' sign extend the shorter
If delta > 0 Then
' a = a
If b.s[Len(b.s)-1] And 128 Then i = 255 Else i = 0
b.s = b.s + String(delta, Chr(i) ) ' extend b
Else
If aa.s[Len(aa.s)-1] And 128 Then i = 255 Else i = 0
a.s = a.s + String(-delta, Chr(i) ) ' extend a
' b = b
End If
End If ' a and b are now the same length
'------------------------------------------------------------
' accumulate b into a
blocks = Len(a.s) \ 4
Dim As ULongInt sum = 0 ' clear carry
Dim As ULong carry
Dim As ULong Ptr pa, pb
pa = Cast(ULong Ptr, StrPtr(a.s) )
pb = Cast(ULong Ptr, StrPtr(b.s) )
For i = 0 To blocks-1
sum = CULngInt(pa[i]) + pb[i] + carry
pa[i] = sum
carry = sum Shr 32
Next i
Bigint.prune(a)
Return a
End Operator
' subtraction
'================================================================
Operator - (ByRef aa As Bigint, ByRef bb As Bigint) As Bigint
Dim As Bigint cc = aa + (-bb)
Return cc
End Operator
' multiply
'================================================================
Operator * (ByRef aa As Bigint, ByRef bb As Bigint) As Bigint
If aa.s = Bigint_s0 Or bb.s = Bigint_s0 Then
Return 0
ElseIf aa.s = Bigint_s1 Then
Return bb
ElseIf bb.s = Bigint_s1 Then
Return aa
ElseIf aa = bb Then
Return bigint.square(aa) ' squaring is faster
ElseIf aa.s = Bigint_s2 Then
Return Bigint.mul2(bb)
ElseIf bb.s = Bigint_s2 Then
Return Bigint.mul2(aa)
Else
' sort out the signs and rectify the inputs
Dim As Bigint a = aa, b = bb, c
Dim As Long sign_a, sign_b, sign_c
sign_a = a.s[Len(a.s)-1] And 128
sign_b = b.s[Len(b.s)-1] And 128
sign_c = sign_a Xor sign_b
If sign_a Then a = -a
If sign_b Then b = -b
'------------------------------------------------------------
' find the dimensions of the problem
Dim As Long i, j, asize, bsize
asize = Len(a.s) ' number of bytes in a
bsize = Len(b.s) ' number of bytes in b
c.s = String(asize + bsize, Chr(0)) ' initialise accumulator
asize = asize \ 4 - 1 ' number of blocks in a
bsize = bsize \ 4 - 1
'------------------------------------------------------------
' pointers into all the Bigints
Dim As ULong Ptr ia, ib, ic
ia = Cast(ULong Ptr, StrPtr(a.s) )
ib = Cast(ULong Ptr, StrPtr(b.s) )
ic = Cast(ULong Ptr, StrPtr(c.s) )
Dim As ULongInt product
Dim As ULong carry
'------------------------------------------------------------
For i = 0 To asize
carry = 0 ' clear carry
For j = 0 To bsize
product = CULngInt(ia[i]) * ib[j] + ic[i+j] + carry
ic[i+j] = product
carry = product Shr 32
Next j
ic[i+j] = carry
Next i
'------------------------------------------------------------
If sign_c = 128 Then c = - c
Bigint.prune(c)
Return c
End If
End Operator
Operator \ (ByRef x As Bigint, ByRef y As Bigint) As Bigint
Dim As Bigint a,b
bigint.div(x,y,a,b)
Return a
End Operator
Operator / (ByRef x As Bigint, ByRef y As Bigint) As Bigint
Print "no floating point division in Big_Integer"
Sleep: End
Return 0
End Operator
Operator Mod (ByRef x As Bigint, ByRef y As Bigint) As Bigint
Dim As Bigint c,d
bigint.div(x,y,c,d)
Return d
End Operator
' exponentiation
'================================================================
Operator ^ (ByRef x As Bigint, ByRef n As LongInt) As Bigint
If n = 2 Then
Return bigint.square(x)
Else
If n < 0 Then
Print "Cannot raise a big integer to a negative power."
Sleep : End
End If
Dim As Long i = 64
Do ' find first set bit
i = i - 1
Loop Until Bit(n, i) Or (i = 0)
i = i + 1
Dim As Bigint pwr
pwr = 1
Do
i = i - 1
pwr = bigint.square(pwr) ' this was a multiply but square is faster
If Bit(n, i) Then pwr = pwr * x
Loop Until i = 0
Return pwr ' pwr was pruned by square and by multiply
End If
End Operator
' NOT. Invert all the bits in a Bigint
'================================================================
Operator Not (ByRef aa As Bigint) As Bigint
Dim As Bigint a = aa
For i As Long = 0 To Len(a.s)-1
a.s[i] = 255 - a.s[i]
Next i
Return a
End Operator
' shift Bigint n bits left
'================================================================
Operator Shl (ByRef a As Bigint, ByVal n As LongInt) As Bigint
If n = 0 Then Return a
If n < 0 Then Return 0
If a.s = Bigint_s0 Then Return 0
Dim As LongInt nblocks = n \ 32
Dim As Byte nbits = n Mod 32
Dim As Bigint s
Dim As ULong m = BitSet(CLng(0), nbits)
s.s = String(nblocks * 4, Chr(0)) + a.s ' put zeros on the rhs
s = m * s
Return s
End Operator
' shift Bigint n bits right, by shifting left nbits and right nblocks
'================================================================
Operator Shr (ByRef a As Bigint, ByVal n As LongInt) As Bigint
If n = 0 Then Return a
If n < 0 Then Return 0
If n > (8 * Len(a.s)) Then Return 0
Dim As LongInt nblocks = n \ 32
Dim As Byte nbits = n Mod 32
Dim As Bigint s = a
Dim As ULongint m = BitSet(CLngint(0), 32 - nbits )
s = m * s ' move bits left
s.s = Right(s.s, Len(s.s) - (nblocks+1)*4 )
If Len(s.s) = 0 Then Return 0
Return s
End Operator
' bitwise AND
'================================================================
Operator And (ByRef aa As Bigint, ByRef bb As Bigint) As Bigint
Dim As Bigint a = aa, b = bb, c
Dim As Long lena, lenb, i
lena = Len(a.s)
lenb = Len(b.s)
If lena > lenb Then
b.s = b.s + String(lena - lenb, Bit(b.s[lenb - 1], 7))
Else
a.s = a.s + String(lenb - lena, Bit(a.s[lena - 1], 7))
End If
c = b
For i = 0 To Len(c.s) - 1
c.s[i] = c.s[i] And a.s[i]
Next i
Return c
End Operator
' bitwise Or
'================================================================
Operator Or (ByRef aa As Bigint, ByRef bb As Bigint) As Bigint
Dim As Bigint a = aa, b = bb, c
Dim As Long lena, lenb, i
lena = Len(a.s)
lenb = Len(b.s)
If lena > lenb Then
b.s = b.s + String(lena - lenb, Bit(b.s[lenb - 1], 7))
Else
a.s = a.s + String(lenb - lena, Bit(a.s[lena - 1], 7))
End If
c = b
For i = 0 To Len(c.s) - 1
c.s[i] = c.s[i] Or a.s[i]
Next i
Return c
End Operator
' bitwise Xor
'================================================================
Operator Xor (ByRef aa As Bigint, ByRef bb As Bigint) As Bigint
Dim As Bigint a = aa, b = bb, c
Dim As Long lena, lenb, i
lena = Len(a.s)
lenb = Len(b.s)
If lena > lenb Then
b.s = b.s + String(lena - lenb, Bit(b.s[lenb - 1], 7))
Else
a.s = a.s + String(lenb - lena, Bit(a.s[lena - 1], 7))
End If
c = b
For i = 0 To Len(c.s) - 1
c.s[i] = c.s[i] Xor a.s[i]
Next i
Return c
End Operator
' bitwise Imp, implication
'================================================================
Operator Imp (ByRef aa As Bigint, ByRef bb As Bigint) As Bigint
Dim As Bigint a = aa, b = bb, c
Dim As Long lena, lenb, i
lena = Len(a.s)
lenb = Len(b.s)
If lena > lenb Then
b.s = b.s + String(lena - lenb, Bit(b.s[lenb - 1], 7))
Else
a.s = a.s + String(lenb - lena, Bit(a.s[lena - 1], 7))
End If
c = b
For i = 0 To Len(c.s) - 1
c.s[i] = (c.s[i] Imp a.s[i])
Next i
Return c
End Operator
' bitwise Eqv, equivalence is the complement of Xor
'================================================================
Operator Eqv (ByRef aa As Bigint, ByRef bb As Bigint) As Bigint
Dim As Bigint a = aa, b = bb, c
Dim As Long lena, lenb, i
lena = Len(a.s)
lenb = Len(b.s)
If lena > lenb Then
b.s = b.s + String(lena - lenb, Bit(b.s[lenb - 1], 7))
Else
a.s = a.s + String(lenb - lena, Bit(a.s[lena - 1], 7))
End If
c = b
For i = 0 To Len(c.s) - 1
c.s[i] = (c.s[i] Eqv a.s[i])
Next i
Return c
End Operator
Operator Abs (ByRef a As Bigint) As Bigint
Dim As Bigint b = a
If 128 And b.s[Len(b.s)-1] Then b = - b
Return b
End Operator
Operator Sgn (ByRef a As Bigint) As Integer
Dim As Long i = 128 And a.s[Len(a.s)-1]
If i Then Return -1 ' is negative
If a.s = Bigint_s0 Then Return 0 ' is zero
Return 1 ' is positive
End Operator
'----------------------------------------------------------------
' operate and assign
Operator Bigint.+= (ByRef rhs As Bigint)
This = This + rhs
End Operator
Operator Bigint.-= (ByRef rhs As Bigint)
This = This - rhs
End Operator
Operator Bigint.*= (ByRef rhs As Bigint)
This = This * rhs
End Operator
Operator Bigint.\= (ByRef rhs As Bigint)
Dim As Bigint c = This, d
bigint.div(c,rhs,This,d)
End Operator
Operator Bigint.mod= (ByRef rhs As Bigint)
Dim As Bigint c, d = This
bigint.div(d,rhs,c,This)
End Operator
Operator Bigint.^= (ByRef rhs As Bigint)
This = This ^ rhs
End Operator
Operator Bigint.shl= (ByRef rhs As LongInt)
This = This Shl rhs
End Operator
Operator Bigint.shr= (ByRef rhs As LongInt)
This = This Shr rhs
End Operator
Operator Bigint.and= (ByRef rhs As Bigint)
This = This And rhs
End Operator
Operator Bigint.or= (ByRef rhs As Bigint)
This = This Or rhs
End Operator
Operator Bigint.xor= (ByRef rhs As Bigint)
This = This Xor rhs
End Operator
Operator Bigint.imp= (ByRef rhs As Bigint)
This = This Imp rhs
End Operator
Operator Bigint.eqv= (ByRef rhs As Bigint)
This = This Eqv rhs
End Operator
'----------------------------------------------------------------
' comparate
Operator = (ByRef a As Bigint, ByRef b As Bigint) As Integer
Dim As Integer c = bigint.compare(a, b)
If c=0 Then Return -1 Else Return 0
End Operator
Operator <> (ByRef a As Bigint, ByRef b As Bigint) As Integer
Dim As Integer c = bigint.compare(a, b)
If c=0 Then Return 0 Else Return -1
End Operator
Operator < (ByRef a As Bigint, ByRef b As Bigint) As Integer
Dim As Integer c = bigint.compare(a,b)
If c = 1 Then Return -1 Else Return 0
End Operator
Operator > (ByRef a As Bigint, ByRef b As Bigint) As Integer
Dim As Integer c = bigint.compare(a,b)
If c = -1 Then Return -1 Else Return 0
End Operator
Operator <= (ByRef a As Bigint, ByRef b As Bigint) As Integer
Dim As Integer c = bigint.compare(a,b)
If c = 1 Or c = 0 Then Return -1 Else Return 0
End Operator
Operator >= (ByRef a As Bigint, ByRef b As Bigint) As Integer
Dim As Integer c = bigint.compare(a, b)
If c = -1 Or c = 0 Then Return -1 Else Return 0
End Operator
'----------------------------------------------------------------
' FOR Bigint, NEXT Bigint and STEP Bigint
' implicit step versions. Implicit step is 1
Operator Bigint.for( )
End Operator
Operator Bigint.step( )
This += 1
End Operator
Operator Bigint.next( ByRef end_cond As Bigint ) As Integer
Return This <= end_cond
End Operator
'----------------------------------------------------------------
' explicit step versions
Operator Bigint.for( ByRef step_var As Bigint )
End Operator
Operator Bigint.step( ByRef step_var As Bigint )
This += step_var
End Operator
Operator Bigint.next( ByRef end_cond As Bigint, ByRef step_var As Bigint ) As Integer
If step_var < 0 Then
Return This >= end_cond
Else
Return This <= end_cond
End If
End Operator
'================================================================
' FUNCTIONS
'================================================================
' remove unnecessary leading blocks, prune time easily earns it's keep
Sub Bigint.prune(ByRef a As Bigint)
a.s = Left(a.s, (((bigint.msbit(a) + 1) \ 32 ) + 1 ) * 4)
' Times and space are improved through the sensible use of prune.
' Addition of opposite signs or subtraction can cancel many blocks.
' A redundant block can appear during multiplication. Square or div2.
' Mul2, Complement, Negate and Absolute do not generate unnecessary blocks.
' Power is pruned internally by the prune in multiply and square.
End Sub
'=======================================================================
' square a number, approaches twice the speed of multiply for big numbers
Function Bigint.square(ByRef aa As Bigint) As Bigint
If aa.s = Bigint_s0 Then
Return aa
ElseIf aa.s = Bigint_s1 Then
Return aa
End If
Dim As Bigint a = aa, c
If (128 And a.s[Len(a.s)-1]) Then a = -a
'------------------------------------------------------------
' find the dimension of the problem
Dim As Long i, j, asize = Len(a.s) ' number of bytes in a
c.s = String(2 * asize, Chr(0) ) ' initialise accumulator
asize = (asize \ 4) - 1 ' highest block number in a
'------------------------------------------------------------
' pointers into all the Bigints
Dim As ULong Ptr pa, pc
pa = Cast(ULong Ptr, StrPtr(a.s) ) ' the base addresses of Bigints
pc = Cast(ULong Ptr, StrPtr(c.s) )
Dim As ULongInt product ' bottom half is data, top half will be carry
Dim As ULong carry, sum
'------------------------------------------------------------
' multiply one triangular corner only
For i = 1 To asize
' pa starts at 1 not zero because 0,0 is on the diagonal
' the second element in a starts at zero
' pc is the accumulator ic = icz + i + j
carry = 0 ' clear carry
For j = 0 To i - 1
product = CULngInt(pa[i]) * pa[j] + pc[j+i] + carry
pc[j+i] = product
carry = product Shr 32
Next j
pc[j+i] = carry ' line of blocks gets one longer each loop
Next i
'------------------------------------------------------------
' double the triangle, cannot do it at same time as squaring diagonal
' because it can cause the product to overflow
carry = 0 ' clear carry
For i = 0 To (2 * asize) + 1
sum = pc[i]
product = CULngInt(sum) + sum + carry
pc[i] = product
carry = product Shr 32
Next i
'------------------------------------------------------------
' square and accumulate the diagonal elements
carry = 0 ' clear carry
For i = 0 To asize
' square the diagonal entry, while propagating carry
sum = pa[i]
product = CULngInt(sum) * sum + pc[i+i] + carry
pc[i+i] = product
carry = product Shr 32
' propagate carry through the following odd block of C
product = CULngInt(pc[i+i+1]) + carry
pc[i+i+1] = product
carry = product Shr 32
Next i
'------------------------------------------------------------
If 128 And c.s[Len(c.s)-1] Then c.s = c.s & Bigint_s0 ' sign is positive
Bigint.prune(c)
Return c
End Function
'=======================================================================
' shift up one bit, low towards high
Function Bigint.mul2(ByRef a As Bigint) As Bigint
If a.s = Bigint_s0 Then Return a
Dim As Long i, sign, sum, carry
Dim As Bigint b = a
sign = b.s[Len(b.s) - 1] And 128 ' snag the msb of highest byte
For i = 0 To Len(b.s) - 1
sum = b.s[i] + b.s[i] + carry
carry = HiByte(sum)
b.s[i] = LoByte(sum)
Next i
If ( b.s[Len(b.s) - 1] And 128 ) <> sign Then
carry = carry * 255
b.s = b.s + Chr(carry,carry,carry,carry) ' sign extend four bytes
End If
Return b
End Function
'=======================================================================
' shift down one bit, high towards low
Function Bigint.div2(ByRef a As Bigint) As Bigint
If a.s = Bigint_s0 Then Return a
Dim As Long i, carry
Dim As Bigint b = a
For i = 0 To Len(a.s)-2 ' all except the top byte of four
b.s[i] = ( b.s[i] \ 2 ) + (128 * (b.s[i+1] And 1))
Next i
i = Len(b.s) - 1
b.s[i] = b.s[i] \ 2
b.s[i] = b.s[i] Or (2 * ( b.s[i] And 64)) ' sign extend the msb
Bigint.prune(b)
Return b
End Function
'=======================================================================
' integer divide, a \ b, return div and mod
Sub Bigint.div(_
ByRef a As Bigint, ByRef bb As Bigint,_
ByRef q As Bigint, ByRef r As Bigint)
If bb.s = Bigint_s0 Then
Print " Division by zero. "
Sleep : End
End If
Dim As Long lena = Len(a.s), lenbb = Len(bb.s)
'------------------------------------------------------------
' Test if Longint Division possible
If (lena <= 8) And (lenbb <= 8) Then ' arguments are one or two blocks
Dim As Longint va = a, vb = bb
q = va \ vb
r = va Mod vb
Exit Sub
End If
'------------------------------------------------------------
' Test if divisor is bigger than dividend
If Abs(bb) > Abs(a) Then
r.s = a.s
q.s = Bigint_s0
Exit Sub
End If
'------------------------------------------------------------
' Read Signs
Dim As Long sa, sb, sq
sa = 128 And a.s[lena-1]
sb = 128 And bb.s[lenbb-1]
sq = sa Xor sb ' sign of the result
'---------------------------------------------------------------------
' Setup variables and pointers
' r=dividend and remainder
' b=divisor
' q=quotient
' sum=interim result for subtraction
Dim As Bigint b = bb
If sb Then b = -b
r.s = a.s
If sa Then r = -r
Dim As UShort Ptr pr, pb, pq
Dim As ULong sum, offset=&hFFFF0000
Dim As ULong rblocks,bblocks,blockshift,rounds,blocks
Dim As ULong lenr = Len(r.s), lenb = Len(b.s)
Dim As UShort qi, carry
Dim As Long substract, i
Dim As ULongInt high48b
Dim As ULongInt high48r = 0
pr = Cast(UShort Ptr, StrPtr(r.s))
pb = Cast(UShort Ptr, StrPtr(b.s))
rblocks = (lenr \ 2) - 1
bblocks = (lenb \ 2) - 1
'------------------------------------------------------------
' convert to 16bit blocklength for effective testdivison
' append exactly two zero blocks and ignore the second leading zeroblock
If pr[rblocks] <> 0 Then
r.s &= Bigint_s0
rblocks += 1
ElseIf pr[rblocks] = 0 And pr[rblocks-1] = 0 Then
rblocks -= 1
End If
If pb[bblocks] <> 0 Then
b.s &= Bigint_s0
bblocks += 1
ElseIf pb[bblocks] = 0 And pb[bblocks-1] = 0 Then
bblocks -= 1
End If
' new pointer after changing string
pr = Cast(UShort Ptr, StrPtr(r.s))
pb = Cast(UShort Ptr, StrPtr(b.s))
blockshift = rblocks - bblocks
rounds = blockshift + 1
'------------------------------------------------------------
' setup quotient to usual 32bit blocklength
If Bit(rounds,0) = 0 Then
q.s = String(rounds * 2, 0) ' quotient
Else
q.s = String((rounds + 1) * 2, 0)
End If
pq = Cast(UShort Ptr, StrPtr(q.s))
pq += blockshift ' start at msb
'------------------------------------------------------------
' start calculation
' most significant bits of divisor are constant
If bblocks = 1 Then
' shift left if divisor is less than 16 bits
high48b = CULng(pb[bblocks-1] Shl 16)
Else
' most significant 32bits of divisor, rounded up
high48b = CULngInt(pb[bblocks-1]) Shl 16 + pb[bblocks-2] + 1
End If
For blocks = 1 To rounds
' msbits: remainder of previous step and following 32bits (= ms48bits)
high48r = CULngInt(pr[rblocks]) Shl 32 + CULng(pr[rblocks-1] Shl 16) + pr[rblocks-2]
' testdivision: because of rounding up result may be 1 too low
qi = (high48r \ high48b )
If qi > 0 Then
' r= r - q*b for every block, begin at lsblock, with carry
carry = 0
For i = blockshift To rblocks
If carry = 0 Then ' don't apply the offset!
sum = pr[i] - (pb[i-blockshift] * qi)
Else ' offset because the overflow 16->32 bit
sum = pr[i] - (pb[i-blockshift] * qi) + carry + offset
End If
pr[i] = sum
carry = sum Shr 16
Next i
End If
'------------------------------------------------------------
' if testdivision was too low, additional substraction needed
substract = 0
' test if remainder > shifted quotient
For i = rblocks To blockshift Step -1
substract = pr[i] - pb[i-blockshift]
If substract <> 0 Then Exit For ' not equal
Next i ' equal, check next block
' if higher, then substract quotient once and thus fix rounding error
If substract >= 0 Then
carry = 1
qi += 1
For i = blockshift To rblocks
sum = pr[i] + CUShort(Not(pb[i-blockshift])) + carry
pr[i] = sum
carry = sum Shr 16
Next i
End If
'------------------------------------------------------------
rblocks -= 1
blockshift -= 1
*pq = qi ' write quotient with pointer
pq -= 1 ' next block of quotient
Next blocks
'------------------------------------------------------------
' finalisation
If Bit(q.s[Len(q.s) - 1], 7) Then q.s &= Bigint_s0
Bigint.prune(r) ' trim result
Bigint.prune(q)
If sq Then q = -q ' Xor of input signs
If sa Then r = -r ' sign of original input A
'------------------------------------------------------------
End Sub
'=======================================================================
Function Bigint.factorial(ByRef a As Bigint) As Bigint
Dim As Bigint f = 1, n = a
Do Until n < 2
f = f * n
n = n - 1
Loop
Return f
End Function
'=======================================================================
'exponentiation modulus
Function Bigint.modpower(ByRef bb As Bigint, ByRef ee As Bigint, ByRef m As Bigint) As Bigint
If m.s = Bigint_s0 Then
Print " Division by zero. "
Sleep : End
End If
If (ee.s[Len(ee.s)-1] And 128) Then
Print "Cannot raise a Bigint to a negative power"
Sleep : End
ElseIf ee.s = Bigint_s0 Then
If Abs(m) = 1 Then Return 0 Else Return 1
End If
Dim As Bigint r = 1, b, y = bb, x ' x is dump for the quotient
Dim As Long bitlen,i,z
Dim As ULong Ptr pee
pee = Cast(ULong Ptr, StrPtr(ee.s))
Dim As ULong spee
spee = *pee ' load first block of exponent in variable
bitlen = bigint.msbit(ee) ' the highest set bit
bigint.div(y,m,x,b) ' initial reduction
' i counts from the lsb to msb, z counts the 32 bits in a block
For i=0 To bitlen-1
If Bit(spee,z) Then ' if bit is set then multiply
y = r * b
bigint.div(y,m,x,r)
End If
y = bigint.square(b)
bigint.div(y,m,x,b)
If z = 31 Then ' reset z, next block
z = 0
pee += 1
spee = *pee
Else
z += 1
End If
Next i
y = r * b ' bitvalue for highest bit=1
bigint.div(y,m,x,r)
Return r
End Function
'=======================================================================
' BIT FUNCTIONS
'=======================================================================
' find the bit position of the first bit that differs from the sign bit
Function Bigint.msbit(ByRef a As Bigint) As Long
Dim As Long i, j, k = 0
i = Len(a.s)
If 128 And a.s[Len(a.s)-1] Then ' negative
Do ' find the highest non-255 byte in the string
i = i - 1
j = a.s[i]
Loop Until (j < 255) Or (i = 0)
j = 255 - j
Else ' positive
Do ' find the highest non-zero byte in the string
i = i - 1
j = a.s[i]
Loop Until (j > 0) Or (i = 0)
End If
' find the highest non-sign bit in the byte
If j And 1 Then k = 1 ' straight code is faster than a loop
If j And 2 Then k = 2
If j And 4 Then k = 3
If j And 8 Then k = 4
If j And 16 Then k = 5
If j And 32 Then k = 6
If j And 64 Then k = 7
If j And 128 Then k = 8
k = k + (i * 8) - 1 ' if no bits differ (-1 or 0) then return -1
Return k
End Function
'=======================================================================
' get the value of a specified bit in a big integer
Function Bigint.Bit_Value(ByRef v As Bigint, ByVal b As ULongInt) As Long
Dim As Long bitval, by = b \ 8
If by < Len(v.s) Then
bitval = Bit(v.s[by], b Mod 8)
Else
If v.s[Len(v.s)-1] And 128 Then bitval = -1 ' the extended sign bit
End If
Return bitval
End Function
'================================================================
' set a specified bit in a big integer
Function Bigint.Bit_Set(ByRef vv As Bigint, ByVal b As ULongInt) As Bigint
Dim As Bigint v = vv
Dim As Long by, bi, delta, sign
by = b \ 8 ' byte number
bi = b Mod 8 ' bit number
delta = by - Len(v.s) + 1
If bi = 7 Then delta = delta + 1 ' protect the sign bit
If delta > 0 Then ' lengthen the number
delta = ((delta + 3)\ 4) * 4
If v.s[Len(v.s)-1] And 128 Then sign = -1 ' the extended sign bit
v.s = v.s + String(delta, sign)
End If
v.s[by] = BitSet(v.s[by], bi)
Return v
End Function
'================================================================
' clear a specified bit in a big integer
Function Bigint.Bit_Reset(ByRef vv As Bigint, ByVal b As ULongInt) As Bigint
Dim As Bigint v = vv
Dim As Long by, bi, delta, sign
by = b \ 8 ' byte number
bi = b Mod 8 ' bit number
delta = by - Len(v.s) + 1
If bi = 7 Then delta = delta + 1 ' protect the sign bit
If delta > 0 Then ' lengthen the number
delta = ((delta + 3)\ 4) * 4
If v.s[Len(v.s)-1] And 128 Then sign = -1 ' the extended sign bit
v.s = v.s + String(delta, sign)
End If
v.s[by] = BitReset(v.s[by], bi)
Return v
End Function
'=====================================================================
'Main Comparation function - faster than substract ...
'=====================================================================
Function Bigint.compare(ByRef a As Bigint, ByRef b As Bigint) As Long
' return -1 for a>b; 1 for a<b and 0 for equal
Dim As Byte signa, signb
signa = 128 And a.s[Len(a.s)-1] ' -1 (true) negative, 0 (false) positive
signb = 128 And b.s[Len(b.s)-1]
'-------------------------------------------
' sign is different - easy:
If signa = 0 And signb = -128 Then
Return -1
ElseIf signa = -128 And signb = 0 Then
Return 1
End If
'-------------------------------------------
' len is different - easy:
If Len(a.s) > Len(b.s) Then
If signa = 0 Then
Return -1
Else
Return 1
End If
ElseIf Len(a.s) < Len(b.s) Then
If signa = 0 Then
Return 1
Else
Return -1
End If
End If
'-------------------------------------------
' compare block for block:
Dim As Long i
Dim as Ulong Ptr pa = CPtr(Ulong Ptr, Strptr(a.s))
Dim as Ulong Ptr pb = CPtr(Ulong Ptr, Strptr(b.s))
For i = (Len(a.s) \ 4 ) - 1 To 0 Step -1
If pa[i] > pb[i] Then
Return -1
ElseIf pa[i] < pb[i] Then
Return 1
End If
Next i
Return 0
End Function
'================================================================
' CAST AND CONVERSION FUNCTIONS
'================================================================
Operator Bigint.cast() As Byte ' CByte
If This > 127 Or This < -128 Then
Print " Overflow in BigInt to Byte conversion. "
Print This
Sleep : End
End If
Return *CPtr(Byte Ptr, StrPtr(this.s))
End Operator
'----------------------------------------------------------------
Operator Bigint.cast() As UByte ' CUByte
If This > 255 Or This < 0 Then
Print " Overflow in BigInt to Ubyte conversion. "
Print This
Sleep : End
End If
Return *CPtr(UByte Ptr, StrPtr(this.s))
End Operator
'----------------------------------------------------------------
Operator Bigint.cast() As Short ' CShort
If This > 32767 Or This < -32768 Then
Print " Overflow in BigInt to Short conversion. "
Print This
Sleep : End
End If
Return *CPtr(Short Ptr, StrPtr(this.s))
End Operator
'----------------------------------------------------------------
Operator Bigint.cast() As UShort ' CUShort
If This > 65535 Or This < 0 Then
Print " Overflow in BigInt to UShort conversion. "
Print This
Sleep : End
End If
Return *CPtr(UShort Ptr, StrPtr(this.s))
End Operator
'----------------------------------------------------------------
Operator Bigint.Cast() As Long ' CLng
If Len(this.s) <> 4 Then
Print " Overflow in BigInt to Long conversion. "
Print This
Sleep : End
End If
Return *CPtr(Long Ptr, StrPtr(this.s))
End Operator
'----------------------------------------------------------------
Operator Bigint.Cast() As ULong ' CULng
If This > 4294967295 Or This < 0 Then
Print " Overflow in BigInt to ULong conversion. "
Print This
Sleep : End
End If
Return *CPtr(ULong Ptr, StrPtr(this.s))
End Operator
'----------------------------------------------------------------
Operator Bigint.Cast() As LongInt ' CLongInt
Dim As String s = this.s
If Len(s) > 8 Then
Print " Overflow in BigInt to LongInteger conversion. "
Print This
Sleep : End
End If
If Len(s) = 4 Then ' sign extend 4 byte integer to 8 byte LongInt
If Bit(s[3], 7) Then
s &= Bigint_s_1
Else
s &= Bigint_s0
End If
End If
Return *CPtr(LongInt Ptr, StrPtr(s))
End Operator
'----------------------------------------------------------------
Operator Bigint.Cast() As ULongInt ' CULongInt
Dim As String s = this.s
If This > 18446744073709551615 Or This < 0 Then
Print " Overflow in BigInt to LongInteger conversion. "
Print This
Sleep : End
End If
If Len(s) = 4 Then s &= Bigint_s0 'extend to len 8
Return *CPtr(ULongInt Ptr, StrPtr(s))
End Operator
'----------------------------------------------------------------
Operator Bigint.Cast() As Integer ' CInt
Dim a As Integer = 2147483647
If (a + 1) > a Then '64bit
Return CLngInt(This)
Else
Return CLng(This)
End If
End Operator
'----------------------------------------------------------------
Operator Bigint.Cast() As UInteger ' CUint
Dim a As UInteger = 4294967295
If (a + 1) > a Then '64bit
Return CULngInt(This)
Else
Return CULng(This)
End If
End Operator
'----------------------------------------------------------------
Operator Bigint.Cast() As Single ' CSng
Dim As Bigint b = This
Dim As ULong ul, Sign_Bit
If Bit( b.s[Len(b.s) - 1], 7) Then ' extract the sign bit
Sign_Bit = CULng(1) Shl 31 ' remember the sign
b = -b ' rectify
End If ' b is now a positive BigInt
' overflows single?
If Len(b.s) > 32 Then ' 32 bytes * 8 bits per byte = 256 bits
' special case of sign bit = lead block &FFFF
If Len(b.s) = 36 And ( b.s[35] And b.s[34] And b.s[33] And b.s[32] ) = &hFF Then
ul = Sign_Bit Or &hFF000000 ' all ones exponent
Return *Cast(Single Ptr, @ul) ' bit pattern is a double
End If
Print " Overflow in BigInt to Single conversion. "
Print This
Sleep : End
End If
If b = 0 Then Return 0 ' test for simple zero
Dim As LongInt expo = 8 * Len(b.s) + 126 ' = bits + expo_bias - 1
' if needed for the conversion, extend tail with two LS blocks of zero
If Len(b.s) = 4 Then b.s = Bigint_s0 + b.s
' the ms block still contains the data, so no change to expo
Dim As UByte Ptr ubp = StrPtr(b.s) + Len(b.s) - 1 ' point to the MSbyte
Dim As Long i
For i = 0 To 4 ' find the leading non-zero byte, MS block may be zero
If *ubp > 0 Then Exit For
ubp = ubp - 1
expo = expo - 8 ' expo reduction of 8 bits per zero byte skipped
Next i ' ubp now points to the MS non-zero byte
ul = *Cast(ULong Ptr, ubp - 3) ' normalize bytes, left justify
For i = 31 To 24 Step -1 ' find the MS set bit
If Bit(ul, i) Then Exit For
expo = expo - 1
Next i ' i now points to MSbit
ul = ul Shr (i - 23) ' shift right to put MSbit in bit 23
ul = BitReset(ul, 23) ' kill only the implicit bit now in bit 52
ul = Sign_Bit Or (expo Shl 23) Or ul ' build the single
Return *Cast(Single Ptr, @ul) ' return the bit pattern as a double
End Operator
'----------------------------------------------------------------
Operator Bigint.Cast() As Double ' CDbl
Dim As Bigint b = This
Dim As ULongInt uli, Sign_Bit
If Bit( b.s[Len(b.s) - 1], 7) Then ' extract the sign bit
Sign_Bit = CULngInt(1) Shl 63 ' remember the sign
b = -b ' rectify
End If ' b is now a positive BigInt
' overflows double? if mag > 1.797693134862310e308 = signed infinity
If Len(b.s) > 128 Then ' 128 bytes * 8 bits per byte = 1024 bits
' special case of sign bit = entire block
If Len(b.s) = 132 And ( b.s[131] And b.s[130] And b.s[129] And b.s[128] ) = &hFF Then
uli = Sign_Bit Or &h7FF0000000000000 ' all ones exponent
Return *Cast(Double Ptr, @uli) ' bit pattern is a double
End If
Print " Overflow in BigInt to Double conversion. "
Print This
Sleep : End
End If
If Len(b.s) = 4 Then ' test for simple zero
If ( b.s[3] Or b.s[2] Or b.s[1] Or b.s[0] ) = 0 Then Return 0
End If
Dim As LongInt expo = 8 * Len(b.s) + 1022 ' = bits + expo_bias - 1
' if needed for the conversion, extend tail with two LS blocks of zero
If Len(b.s) < 12 Then b.s = Chr(0,0,0,0, 0,0,0,0) + b.s
' the ms block still contains the data, so no change to expo
Dim As UByte Ptr ubp = StrPtr(b.s) + Len(b.s) - 1 ' point to the MSbyte
Dim As Long i
For i = 0 To 4 ' find the leading non-zero byte, MS block may be zero
If *ubp > 0 Then Exit For
ubp = ubp - 1
expo = expo - 8 ' expo reduction of 8 bits per zero byte skipped
Next i ' ubp now points to the MS non-zero byte
uli = *Cast(ULongInt Ptr, ubp - 7) ' normalize bytes, left justify
For i = 63 To 56 Step -1 ' find the MS set bit
If Bit(uli, i) Then Exit For
expo = expo - 1
Next i ' i now points to MSbit
uli = uli Shr (i - 52) ' shift right to put MSbit in bit 52
uli = BitReset(uli, 52) ' kill only the implicit bit now in bit 52
uli = Sign_Bit Or (expo Shl 52) Or uli ' build the double
Return *Cast(Double Ptr, @uli) ' return the bit pattern as a double
End Operator
'----------------------------------------------------------------
' unpack a straight binary string to a decimal ascii string
Operator Bigint.cast() As String
Dim As bigint b
Dim As String d
b = This
d = Chr(0) ' initial decimal output string
Dim As Integer i, j, sign
Dim As UInteger carry
' if negative then negate, append the sign later
If b.s[Len(b.s) - 1] And 128 Then ' negative
sign = - 1
b = -b
End If
' change from base 2 to base 100
For j = Len(b.s) - 1 To 0 Step - 1 ' each byte in string is base 256
carry = b.s[j] ' the next byte to add after multiply
Dim as UByte Ptr k = StrPtr(d) + Len(d) - 1 ' added pointer
For i = Len(d) - 1 To 0 Step - 1
carry += *k * 256
*k = carry Mod 100
carry \= 100
k -= 1
Next i
Do While carry > 0 ' output string overflows
d = Chr(carry Mod 100) & d ' extend output string as needed
carry = carry \ 100
Loop
Next j
' need to split the value of the byte into 2 digits and to convert into ASCII
Dim As String d2 = Chr(d[0] \ 10 + Asc( "0")) + Chr(d[0] Mod 10 + Asc( "0"))
' remove a possible 0 in front of the number
d2 = LTrim(d2, "0")
If Len(d2) = 0 Then Return "+0"
For i = 1 To Len(d) - 1
d2 &= Chr(d[i] \ 10 + Asc( "0")) & Chr(d[i] Mod 10 + Asc( "0"))
Next i
If sign Then d2 = "-" & d2 Else d2 = "+" & d2
Return d2
End Operator
'----------------------------------------------------------------
' Casting to Bigint
Function CBig Overload(a as Byte) as Bigint
Dim As Bigint b=a
Return b
End Function
Function CBig Overload(a as UByte) as Bigint
Dim As Bigint b=a
Return b
End Function
Function CBig Overload(a as Short) as Bigint
Dim As Bigint b=a
Return b
End Function
Function CBig Overload(a as UShort) as Bigint
Dim As Bigint b=a
Return b
End Function
Function CBig Overload(a as Integer) as Bigint
Dim As Bigint b=a
Return b
End Function
Function CBig Overload(a as UInteger) as Bigint
Dim As Bigint b=a
Return b
End Function
Function CBig Overload(a as Long) as Bigint
Dim As Bigint b=a
Return b
End Function
Function CBig Overload(a as ULong) as Bigint
Dim As Bigint b=a
Return b
End Function
Function CBig Overload(a as LongInt) as Bigint
Dim As Bigint b=a
Return b
End Function
Function CBig Overload(a as ULongInt) as Bigint
Dim As Bigint b=a
Return b
End Function
Function CBig Overload(a as Single) as Bigint
Dim As Bigint b=a
Return b
End Function
Function CBig Overload(a as Double) as Bigint
Dim As Bigint b=a
Return b
End Function
Function CBig Overload(a as String) as Bigint
Dim As Bigint b=a
Return b
End Function
'----------------------------------------------------------------
' convert a Bigint to binary (0110000111 etc.)
Function Bin(ByRef s As Bigint) As String
Dim As Long i
Dim As String h ' lsb is string[0] = little endian
For i = Len(s.s)-1 To 0 Step -1
h = h & Bin(s.s[i], 8)
Next i
Return h
End Function
Function Bin (ByRef a As Bigint,ByRef n As ULong) As String
Dim result As String = Right(Bin(a),n)
Return result
End Function
'----------------------------------------------------------------
' convert a Bigint to hexadecimal
Function Hex (ByRef s As Bigint) As String
Dim As Long i
Dim As String h ' lsb is string[0] = little endian
For i = Len(s.s)-1 To 0 Step -1
h = h & Hex(s.s[i], 2)
Next i
h = ltrim(h,any "0")
If Len(h) <> 0 Then
Return h
Else
Return "0"
End If
End Function
Function Hex (ByRef a As Bigint,ByRef n As ULong) As String
Dim result As String = Right(Hex(a),n)
If a < 0 And Len(result) < n Then
result = String(n-Len(result),"F") & result
ElseIf a >= 0 And Len(result) < n Then
result = String(n-Len(result),"0") & result
End If
Return result
End Function
'----------------------------------------------------------------
' convert a Bigint to unsigned hexadecimal (trimmed)
Function Uhex(ByRef s As Bigint) As String
Dim As Long i
Dim As Bigint a = s
If 128 And a.s[Len(a.s)-1] Then ' Bigint is negative
Print "cannot convert negative to uniform"
Sleep : End
End If
Dim As String h ' hex is big-endian
For i = Len(s.s)-1 To 0 Step -1
h = h & Hex(s.s[i], 2)
Next i
If Len(h) <> 8 Then h=LTrim(h,"00000000")
Return h
End Function
Function UhexT(ByRef s As Bigint,ByRef n As ULong) As String
Dim As String h = Uhex(s)
If Len(h) < n Then
h = String(n - Len(h),"0") & h
End If
Return h
End Function
'----------------------------------------------------------------
' convert a Bigint to octal
Function Oct(ByRef a As Bigint) As String
Dim As String b, c
Dim As Bigint s = a
If 128 And a.s[Len(a.s)-1] Then ' extend the sign
s.s = s.s & Chr(255,255,255,255)
Else
s.s = s.s & Bigint_s0
End If
Dim As Long i
Dim As ULongInt u
For i = 1 To Len(a.s) Step 3
b = Mid(s.s, i, 3) ' take three bytes = 24 bits
u = b[0] + 256 * (b[1] + 256 * b[2])
c = Oct(u, 8) + c ' eight symbols per three bytes
Next i
Return c
End Function
Function Oct (ByRef a As Bigint,ByRef n As ULong) As String
Dim result As String = Right(Hex(a),n)
Return result
End Function
'----------------------------------------------------------------
' Bigint to twos compliment binary
Function MkBigint(ByRef a As Bigint) As String
Dim As String s = a.s
Return s
End Function
'----------------------------------------------------------------
' Bigint to unsigned binary
Function MkUBigint(ByRef a As Bigint) As String
Dim As String s
If 128 And a.s[Len(a.s)-1] Then ' Bigint is negative
Print "cannot convert negative to unsigned"
Sleep : End
Else
s = a.s
End If
s = RTrim(s,Bigint_s0)
Return s
End Function
'----------------------------------------------------------------
' Val (as Bigint)
Function ValBigint(ByRef aa As String) As Bigint
Dim c As Bigint = aa ' VAL is integrated in the Constructor(As String)
Return c
End Function
'----------------------------------------------------------------
' Val (as UnSignedBigint)
Function ValUBigint(ByRef aa As String) As Bigint
Dim c As Bigint = aa & "u" ' use the constructor with unsigned suffix
Return c
End Function
'----------------------------------------------------------------
' twos compliment binary to Bigint
Function CVBigint (ByRef a As String) As Bigint
Dim As Bigint b
If Len(a) = 0 Then
b = 0
Return b
End If
If (128 And a[Len(a)-1]) Then ' negative, pad to blocklen with FF
b.s = a & String((4-Len(a) Mod 4),Chr(255))
Else
b.s = a & String((4-Len(a) Mod 4),Chr(0)) ' positive, pad to blocklen with 00
End If
Return b
End Function
'----------------------------------------------------------------
' unsigned binary to Bigint
Function CVUBigint(ByRef a As String) As Bigint
Dim As Bigint b
Dim As Long pad
pad = 4 - (Len(a) Mod 4)
If (pad = 4) And (Len(a) <> 0) Then pad = 0
b.s = a & String(pad,0) ' Pad to blocklen
If (128 And b.s[Len(b.s)-1]) Then b.s &= Bigint_s0 ' make it positive
Return b
End Function
'----------------------------------------------------------------
'===============================================================
' E n d o f B i g I n t e g e r C o d e
'===============================================================
| Zusätzliche Informationen und Funktionen |
|---|
- Das Code-Beispiel wurde am 25.11.2013 von stephanbrunker angelegt.
- Die aktuellste Version wurde am 29.07.2018 von stephanbrunker gespeichert.
|
|
Wer ist online?
Buchempfehlung
FreeBASIC-Chat
FreeBASIC-Nachrichten jetzt auch über Twitter erhalten. Follow us!
Bearbeiten
Versionen