$ filename "KaratsubaMultiply.exe"

'uses rtl32
'uses rtl64

% karatsuba

uses BigNumbers.o2bas



function multiply(string a,b,int d=0) as string
'http://www.rosettacode.org/wiki/Arbitrary-precision_integers_(included)#Pure_C
'https://en.wikipedia.org/wiki/Multiplication_algorithm#Karatsuba_multiplication
'https://en.wikipedia.org/wiki/Karatsuba_algorithm#Pseudocode
   string high1,high2,low1,low2
   string z0,z1,z2, t1,t2,t1a, t1b, r1,r2,r3, result
   int le,le1,le2,t
   d++
   'print d cr 'nesting level 1..6..
   '
   le1=len(a) : le2=len(b)
   if le1 < 256 or le2 < 256 then
     return bignum::mul1(a,b)
   end if
   if le1 < 4000 or le2 < 4000 then
     return bignum::mul2(a,b)
   end if
   '
   'calculates the size of the numbers
   le=le1
   if le2>le1 then le=le2
   int m=le
   int m2 = m>>1 

   'split the digit sequences about the middle
   'high1, low1 = split_at(num1, m2)
   'high2, low2 = split_at(num2, m2)

   high1=left(a,le1-m2) : low1=mid(a,-m2)
   high2=left(b,le2-m2) : low2=mid(b,-m2)

    ' 3 calls made to numbers approximately half the size
    ' z0 = karatsuba(low1,low2)
    ' z1 = karatsuba((low1+high1),(low2+high2))
    ' z2 = karatsuba(high1,high2)

    z2 = multiply(high1,high2,d)
    z0 = multiply(low1,low2,d)

      t1a=bignum::add(high1,low1)
      t1b=bignum::add(high2,low2)
    z1 = multiply(t1a,t1b,d)

    ' return (z2*10^(2*m2))+((z1-z2-z0)*10^(m2))+(z0)           

      t1 = bignum::sub(z1,z2)
      t2 = bignum::sub(t1,z0)

      r1=z2 & stringbytes(m2+m2,48)     
      r2=t2 & stringbytes(m2,48)
      r3 = bignum::add(r1, r2)
    result = bignum::add(r3, z0)
 
   return bignum::zerotrim result
end function


'TEST

'Use Windows functions
! GetTickCount lib "kernel32.dll" () as sys

function str_exp(int b, int n) as string
========================================
  string result, a

  result="1"
  a=str(b)

  if n=0 then return result
 
  while n > 0
    if n and 1 then 'mod(n,2)=1 then
      result = multiply(result, a)     
    end if
    n = n >> 1 '\= 2
    if n=0 then exit while   
    a = multiply(a, a)
  wend

  return result
end function


int x,y,z
string z1,z2,z3
int t1,t2
'
for z=3 to 4'4
for x=2 to 3
for y=1 to 2
   z1=str_exp(x,y)
   printl z1
   z2=str_exp(z,z1)
   print z2 cr
   
   t1=GetTickCount()
   print "Calculatating 5^" z "^" x "^" y " which is 5^" z2 & cr   cr
   print "Please wait ..." & cr cr

   z3=str_exp(5,z2)
   t2=(GetTickCount()-t1)/1000
   printl "Length for result of 5^" z "^" x "^" y " is " len(z3) cr
   printl "First 20 digits: " left(z3,20) cr
   printl "Last 20 digits:  " mid(z3,-20) cr
   printl "This took " str(t2,3) " seconds" cr
   print cr
next y
next x
next z

print "Enter ... " cr : waitkey

/*
GCC compiled:
Length of 5^4^3^2 is 183231
First 20 digits:  62060698786608744707
Last 20 digits:   92256259918212890625
*/                                                                                    
             