Author Topic: Calculating with big numbers  (Read 135 times)

0 Members and 1 Guest are viewing this topic.

Arnold

  • Sr. Member
  • ****
  • Posts: 448
Calculating with big numbers
« on: October 09, 2017, 01:34:55 AM »
Hi Charles,

this is a programming task of Rosettacode
http://www.rosettacode.org/wiki/Arbitrary-precision_integers_(included)

and I am not very happy with my solution.
Calculating up to 1000 digits there is no problem but then the execution time grows very strongly. To calculate 5^19683 will take as much time on my old notebook as for gcc compiled 5^262144. I could not calculate this power with the Oxygenbasic code and I gave up after 2 hours.
Can I improve this code in some way? I know I use too many garbage collections, but I do not know how I can get rid of them.

Roland

Code: [Select]
'http://www.rosettacode.org/wiki/Arbitrary-precision_integers_(included)
'http://www.rosettacode.org/wiki/Arbitrary-precision_integers_(included)#Pure_C

$ filename "Arbitr_Prec.exe"
'include "$/inc/RTL32.inc"

include "$/inc/console.inc"

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

' return = a * b
' Handling of negatives and zeros is not here, since not needed.
function str_mult(string First, string Second) as string

    string a=First, b=Second, r ' result
   
    int ax = 0, bx = 0, rx = 0, al, bl
    int carry, n
    int an, bn, rn, t 'auxiliary

    al = len(a) : bl = len(b)   
    r = space (al+bl)
       
    ' grade-school method of multiplication
    for ax = al to 1 step -1
        carry = 0
        for bx = bl to 1 step -1
           rx = ax + bx
           an=mid(a,ax,1) : bn=mid(b,bx,1) : rn=mid(r,rx,1)                   
           n=an * bn + rn + carry
           mid r,rx, str(mod(n, 10))           
           carry = n \ 10             
        next bx
        rx--       
        int t=mid(r,rx,1) + carry
        mid r,rx,t
    next ax

    ' omit possible first place - zero
    if left(r,1) = "0" then return mid(r,2)   
    return r
end function


function str_exp(int b, int n) as string
  string res, a

  res="1"
  a=str(b)

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

  return res
end function


for z=3 to 3 '4
for x=2 to 3
for y=1 to 2
   z1=str_exp(x,y)
   printl z1
   z2=str_exp(z,z1)
   printl z2

   t1=GetTickCount()
   printl "Calculatating 5^" z "^" x "^" y " which is 5^" z2 & cr

   z3=str_exp(5,z2)

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

printl "Enter ..." : waitkey

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

output:

2
9
Calculatating 5^3^2^1 which is 5^9

Length for result of 5^3^2^1 is 7
First 20 digits: 1953125
Last 20 digits:  1953125
This took 0 seconds

4
81
Calculatating 5^3^2^2 which is 5^81

Length for result of 5^3^2^2 is 57
First 20 digits: 41359030627651383743
Last 20 digits:  82906055450439453125
This took 0.015 seconds

3
27
Calculatating 5^3^3^1 which is 5^27

Length for result of 5^3^3^1 is 19
First 20 digits: 7450580596923828125
Last 20 digits:  7450580596923828125
This took 0 seconds

9
19683
Calculatating 5^3^3^2 which is 5^19683

Length for result of 5^3^3^2 is 13758
First 20 digits: 67080354005587923286
Last 20 digits:  06702899932861328125
This took 277.9 seconds

Enter ...

Charles Pegge

  • Author
  • *****
  • Posts: 3266
    • Oxygen Basic
Re: Calculating with big numbers
« Reply #1 on: October 10, 2017, 02:22:15 AM »

Hi Roland,

Here is your multiplier with most of the stringware replaced with byte overlays. I have verified it against examples/math/mulhuge.o2bas.

It does your large calculation in 1.4 seconds on my notebook, but I have not looked at the rest of the code.

Code: [Select]
function str_mult(string a,b) as string

    string r ' result
   
    int ax = 0, bx = 0, rx = 0, al, bl
    int carry, n
    int an, bn, rn, t 'auxiliary

    al = len(a) : bl = len(b)   
    r = stringbytes (al+bl,48) 'FILL 0000..
    byte aa at strptr a
    byte bb at strptr b
    byte rr at strptr r
    indexbase 1
    macro midval(b,i){ b[i]-48 } 'ASCII TO DIGIT VALUE
       
    ' grade-school method of multiplication
    for ax = al to 1 step -1
        carry = 0
        for bx = bl to 1 step -1
           rx = ax + bx
           an=midval(aa,ax)
           bn=midval(bb,bx)
           rn=midval(rr,rx)                   
           n=an * bn + rn + carry
           rr[rx]=mod(n,10)+48
           carry = n \ 10             
        next bx
        rx--       
        'int t=midval(r,rx) + carry
        rr[rx]+=carry 'FINAL CARRY
    next ax

    ' omit possible first place - zero
    if rr[1]=48 then return mid(r,2)   
    return r
end function

/*
a="4225678345671239269412376809064563457830987094214727890460503258124429509352125309123467889090967803"
b="9034098905607687535809534264324877930963452357175890408065643678402109067842097453251206046098340056"

ans1="38175196118078646234796018047078508876075927016257021775208337789644882612626336290568135391711676360295755214293161291134825045621348587392814424177752560047244220841334894056033458232344500841216968"
ans2=str_mult(a,b)
printl ans1
printl ans2
waitkey
*/


Arnold

  • Sr. Member
  • ****
  • Posts: 448
Re: Calculating with big numbers
« Reply #2 on: October 10, 2017, 11:29:53 PM »
Hi Charles,

this is really impressive. Your modifications for str_mult reduced the execution time by almost 99% (!). I also get a correct result for 5**262144. This example demonstrates again how important it is to deal with strings efficiently in time-consuming tasks and I will study the code very carefully. I learned about a new keyword: stringbytes; I already missed this functionality sometimes.

It is also remarkable that the function str_mult runs almost as fast as the multiply function in \examples\math\MulHuge.o2bas which is realized with assembler statements. This is another proof of OxygenBasic's speed.

Thank you again for your help. It will help me to understand some more examples where this approach is applied too.

Roland
« Last Edit: October 10, 2017, 11:52:44 PM by Arnold »

Charles Pegge

  • Author
  • *****
  • Posts: 3266
    • Oxygen Basic
Re: Calculating with big numbers
« Reply #3 on: October 11, 2017, 10:49:20 PM »
Hi Roland,

I'm glad you found it helpful.

There are lots of ways to make it go even faster. For instance, working with base10000 digits instead of base10 digits. In theory, this would give you a 16x boost in speed.

I've also got an algorithm for big number division, something I have always wanted to do. It is based on a series of shifting subtractions, similar to school-teacher punishment long division :)

Will post it later.


Charles Pegge

  • Author
  • *****
  • Posts: 3266
    • Oxygen Basic
Re: Calculating with big numbers
« Reply #4 on: October 15, 2017, 02:09:36 AM »
First, the numbers are converted from ascii to binary by subtracting 48 from each digit. Then a lookup table is prepared with multiples of the divisor x1 x2 x3 .. x9. Then we are ready for the main algorithm.

Comparisons are made using the table values, and the selected table index becomes the value of the current digit. The selected number is then subtracted leaving a possible remainder. We then shift 1 place to the right for the next digit, and repeat the process until the right hand end is reached.

The devil is in the detail :)

Code: [Select]
  function div(string sa,sb) as string
  ====================================
  string a="0"+zerotrim(sa)
  int la=len a
  string b="0"+zerotrim(sb)
  int lb=len b
  tobin a : tobin b
  int lo=la-lb+1
  string o=nuls lo
  '
  'CREATE LOOKUP TABLE 1xB 2xB ..9xB
  indexbase 1
  string tbl[9]
  int i,j,k
  sbyte *ti,*tj
  tbl[1]=b : j=1
  for i=2 to 9
    tbl[i]=b
    @ti = strptr(tbl[i])+lb-1
    @tj = strptr(tbl[j])+lb-1
    cy=0
    for k=1 to lb 'add digits
      ti=ti+tj+cy
      if ti>=10 then ti-=10 : cy=1 else cy=0
      @ti-- : @tj--
    next
    j++
  next
  '
  sbyte aa at strptr(a)
  sbyte bb at strptr(b)
  sbyte oo at strptr(o)
  sbyte d,carry
  int i,j,k,c,cm
  sys pa=@aa
  for i = 1 to lo 'setting dividend digits
    c=9 'default digit
    for k=1 to 9
      cm=CompareMid(i,a,tbl[k]) 'lookup multiple of b
      if cm<0 then c=k-1 : exit for 'set digit
    next
    if c then 'perform subtraction
      @aa=pa+i+lb-2
      @bb=strptr(tbl[c])+lb-1
      cb=lb
      carry=0
      for j=1 to lb 'subtraction
        d=aa-bb-carry
        if d<0 then d+=10 : carry=1 else carry=0
        aa=d
        @bb-- : @aa--
      next 'next subtraction digit
    end if 'subtraction
    oo=c 'set dividend digit as subraction count
    @oo++
  next  'next dividend digit
  remainder=nullTrim a
  o=nulltrim o
  if not o then o=chr 0
  toasc o
  toasc remainder
  return o
  end function
  '

supporting functions
Code: [Select]
  macro tobin(s,  bb,le) {
  ========================
  int le=len s
  sbyte bb at strptr s
  while le
    bb-=48
    @bb++ : le--
  wend
  }
  macro toasc(s,  bb,le) {
  ========================
  int le=len s
  sbyte bb at strptr s
  while le
    bb+=48
    @bb++ : le--
  wend
  }

  function compareMid(int i, string a,b) as int
  =============================================
  i-- 'convert to base 0
  int la=len(a)-i
  int lb=len b
  byte aa at i+strptr a
  byte bb at strptr b
  if lb>la
    return -1
  else
    while lb
      if aa>bb then return 1
      if aa<bb then return -1
      lb-- : @aa++ : @bb++
    wend
  end if
  end function
  '
  function nullTrim(string a) as string
  =====================================
  byte aa at strptr a
  int la=len a
  int c=1
  while la
    if aa>0 then return mid a,c
    la-- : c++ : @aa++
  wend
  end function
  '
  function zeroTrim(string a) as string
  =====================================
  byte aa at strptr a
  int la=len a
  int c=1
  while la
    if aa>48 then return mid a,c
    la-- : c++ : @aa++
  wend
  end function


Arnold

  • Sr. Member
  • ****
  • Posts: 448
Re: Calculating with big numbers
« Reply #5 on: October 15, 2017, 09:32:16 AM »
Hey Charles,

you have already finished what I was about to do. But since I will never achieve your accuracy I will focus on your code. Will you work out the functions as a library? There are already many things possible, e.g.:

The idea to explore big numbers arose from Robbek's message:
http://www.oxygenbasic.org/forum/index.php?topic=1231.msg11994#msg11994

and this can be done easily with your functions:

Code: [Select]
...
print "fibonacci" cr
sa=0
sb=1
time1=GetTickCount()
for i=1 to 4000
  sr=add sa,sb
'  print i tab sr cr
  sa=sb
  sb=sr
next
i-=1 ' adjust i
time=(GetTickCount()-time1) / 1000
print "fibo " i " = " sr cr
printl "This took " time " Seconds" cr
...

With a modern pc this will take 0.0078 seconds or less and this is really fast if one takes into account that the code does not play with changing the bases and/or grouping the digits to integers.

I have not yet explored the functions intensively but I noticed you also added some code for negative numbers. Is it already possible to calculate 1-1000, or -100*1 or -2*-3 ? I was about to implement Karatsuba's method for multiplication, but I noticed that I will need to calculate negative interim results in some cirumstances.

Roland
« Last Edit: October 15, 2017, 11:27:13 PM by Arnold »