Author Topic: fibonacci(4784969)  (Read 31634 times)

Offline AIR

  • BASIC Developer
  • Posts: 932
  • Coder
Re: fibonacci(4784969)
« Reply #30 on: April 22, 2019, 08:40:18 AM »
riveraa@dpi:~/Projects/Python/gmp $ time ./fibogmp.py
1000000


real    0m0.362s
user    0m0.362s
sys     0m0.000s

Is this on PC or RPI 3 B+ ?


RPI 3B+


riveraa@dpi: ~/Projects/Python/gmp $ uname -a
Linux dpi 4.14.98-v7+ #1200 SMP Tue Feb 12 20:27:48 GMT 2019 armv7l GNU/Linux


AIR.



Offline jalih

  • Advocate
  • Posts: 109
Re: fibonacci(4784969)
« Reply #31 on: April 22, 2019, 08:47:31 AM »
RPI 3B+

Very good time! How long does it take to output the whole result?

Offline AIR

  • BASIC Developer
  • Posts: 932
  • Coder
Re: fibonacci(4784969)
« Reply #32 on: April 22, 2019, 08:58:29 AM »
Printing is always expensive....



real    0m1.895s
user    0m1.817s
sys     0m0.050s

« Last Edit: April 22, 2019, 09:26:46 AM by AIR »

Offline John

  • Forum Support / SB Dev
  • Posts: 3510
    • ScriptBasic Open Source Project
Re: fibonacci(4784969)
« Reply #33 on: April 22, 2019, 09:52:04 AM »
GMP is the only way to fly. (with unreal numbers)
« Last Edit: April 22, 2019, 06:37:38 PM by John »

Offline jalih

  • Advocate
  • Posts: 109
Re: fibonacci(4784969)
« Reply #34 on: April 23, 2019, 07:13:40 AM »
Printing is always expensive....

real    0m1.895s
user    0m1.817s
sys     0m0.050s

I made a small modification that gave about 26 percent more speed to fibonacci calculation. My current 8th version on my PC calculates and prints 1000000 digit number, if output is redirected into a file in 243 msec. Not sure, how RPI 3 B+ would perform.
« Last Edit: April 23, 2019, 07:24:02 AM by jalih »

Offline AIR

  • BASIC Developer
  • Posts: 932
  • Coder
Re: fibonacci(4784969)
« Reply #35 on: April 23, 2019, 12:45:42 PM »

I made a small modification that gave about 26 percent more speed to fibonacci calculation. My current 8th version on my PC calculates and prints 1000000 digit number, if output is redirected into a file in 243 msec. Not sure, how RPI 3 B+ would perform.

Are you saving as text or data?

Because if I convert the number to a byte array in GO, and save that, these are the timings I get on my MBP:

real    0m0.240s
user    0m0.232s
sys    0m0.009s



Offline jalih

  • Advocate
  • Posts: 109
Re: fibonacci(4784969)
« Reply #36 on: April 23, 2019, 01:21:45 PM »
Are you saving as text or data?

I format the result into string as it's a big-float and I want to display just the integer part. After that I print how many digits, followed by a carriage return and the string containing actual million digit number. Output is directed into file from command-line.

Offline AIR

  • BASIC Developer
  • Posts: 932
  • Coder
Re: fibonacci(4784969)
« Reply #37 on: April 23, 2019, 01:41:57 PM »
Jali, you say you "format" the number into a string rather than say "convert".  Interesting.

Offline John

  • Forum Support / SB Dev
  • Posts: 3510
    • ScriptBasic Open Source Project
Re: fibonacci(4784969)
« Reply #38 on: April 23, 2019, 03:39:16 PM »
It's getting to the point you can't even say fibo before the results are returned.

Offline AIR

  • BASIC Developer
  • Posts: 932
  • Coder
Re: fibonacci(4784969)
« Reply #39 on: April 23, 2019, 06:55:01 PM »
Final GO version (I'm fibonacci'd out already!):

Code: Go
  1. package main
  2.  
  3. import (
  4.     "fmt"
  5.     "io/ioutil"
  6.  
  7.     big "github.com/ncw/gmp"
  8. )
  9.  
  10. func Mul(x, y *big.Int) *big.Int {
  11.     return big.NewInt(0).Mul(x, y)
  12. }
  13.  
  14. func Add(x, y *big.Int) *big.Int {
  15.     return big.NewInt(0).Add(x, y)
  16. }
  17.  
  18. func fib(n int) *big.Int {
  19.     v1, v2, v3 := big.NewInt(1), big.NewInt(1), big.NewInt(0)
  20.     s := fmt.Sprintf("%b", n)
  21.     for i := 1; i < len(s); i++ {
  22.         calc := Mul(v2, v2)
  23.         v1, v2, v3 = Add(Mul(v1, v1), calc), Mul(Add(v1, v3), v2), Add(Mul(v3, v3), calc)
  24.         if s[i] == 49 {
  25.             v1, v2, v3 = Add(v1, v2), v1, v2
  26.         }
  27.     }
  28.     return v2
  29. }
  30.  
  31. func main() {
  32.  
  33.     result := fib(4784969).String()
  34.     fmt.Println(len(result))
  35.  
  36.     ioutil.WriteFile("output.txt", []byte(result), 0644)
  37. }

[riveraa@Kiara ~/Projects/go/Fib] $ time ./Fib
1000000

real    0m0.130s
user    0m0.122s
sys    0m0.007s


    Hardware Overview:

      Model Name: MacBook Pro
      Model Identifier: MacBookPro14,1
      Processor Name: Intel Core i5
      Processor Speed: 2.3 GHz
      Number of Processors: 1
      Total Number of Cores: 2
      L2 Cache (per Core): 256 KB
      L3 Cache: 4 MB
      Hyper-Threading Technology: Enabled
      Memory: 16 GB


On my RasPi, the same program:

real    0m2.246s
user    0m2.218s
sys    0m0.040s


AIR.

Offline John

  • Forum Support / SB Dev
  • Posts: 3510
    • ScriptBasic Open Source Project
Re: fibonacci(4784969)
« Reply #40 on: April 23, 2019, 07:13:52 PM »
Nice finish!

Offline AIR

  • BASIC Developer
  • Posts: 932
  • Coder
Re: fibonacci(4784969)
« Reply #41 on: April 23, 2019, 08:52:30 PM »
Admittedly, that laptop is pretty fast considering it's only dual core, but it's only about a year old.

Here's a more realistic type of run on my Linux server, which is running on an old Optiplex 7010 (7yr old HW, 16GB ram, 8TB storage):

riveraa@nas:~/Projects/go/Fib$ time ./Fib
1000000

real    0m0.235s
user    0m0.216s
sys    0m0.016s


Hardware:

Architecture:          x86_64
CPU op-mode(s):        32-bit, 64-bit
Byte Order:            Little Endian
CPU(s):                4
On-line CPU(s) list:   0-3
Thread(s) per core:    1
Core(s) per socket:    4
Socket(s):             1
NUMA node(s):          1
Vendor ID:             GenuineIntel
CPU family:            6
Model:                 58
Model name:            Intel(R) Core(TM) i5-3470 CPU @ 3.20GHz

Offline petelomax

  • BASIC Developer
  • Posts: 9
  • Author of Phix
    • The Phix Programming Language
Re: fibonacci(4784969)
« Reply #42 on: May 18, 2019, 02:13:29 PM »
Code: Text
  1. include mpfr.e
  2.  
  3. mpz res = NULL, prev, next
  4. integer lastn
  5. atom t0 = time()
  6.  
  7. function fibonampz(integer n) -- resumable, works for -ve numbers, yields mpz
  8. integer absn = abs(n)
  9.     if res=NULL or absn!=abs(lastn)+1 then
  10.         if res=NULL then
  11.             prev = mpz_init(0)
  12.             res = mpz_init(1)
  13.             next = mpz_init()
  14.         else
  15.             if n==lastn then return res end if
  16.         end if
  17.         mpz_fib2_ui(res,prev,absn)
  18.     else
  19.         if lastn<0 and remainder(lastn,2)=0 then
  20.             mpz_mul_si(res,res,-1)
  21.         end if
  22.         mpz_add(next,res,prev)
  23.         {prev,res,next} = {res,next,prev}
  24.     end if
  25.     if n<0 and remainder(n,2)=0 then
  26.         mpz_mul_si(res,res,-1)
  27.     end if
  28.     lastn = n
  29.     return res
  30. end function
  31.  
  32. string s = mpz_get_str(fibonampz(4784969))
  33. integer l = length(s)
  34. s[40..-40] = "..."
  35. ?{l,s}
  36. ?elapsed(time()-t0)
  37.  

Offline John

  • Forum Support / SB Dev
  • Posts: 3510
    • ScriptBasic Open Source Project
Re: fibonacci(4784969)
« Reply #43 on: May 18, 2019, 06:29:19 PM »
Looks good Pete!

Is this written in your BASIC?

Is the BIGINT library like GMP?

How long does it take for the 1 mil digit Fibo?
« Last Edit: May 18, 2019, 08:23:14 PM by John »

Offline John

  • Forum Support / SB Dev
  • Posts: 3510
    • ScriptBasic Open Source Project
Re: fibonacci(4784969)
« Reply #44 on: June 11, 2019, 04:02:12 PM »
AIR,

I'm trying to get the GMP extension module that I added Heater's routines to allow me to do the 1 mil fibo using the BIGINT math. I'm getting a seg fault with the return of the string. Your FIBO function code works great. Can you give me a hand adapting that same method to the other functions I added?

What has me confused is in your fibo() function you create a buffer the size of n but the resulting length is 1 million bytes. Should  I use the _ui versions of the function calls? I'm still a bit confused about initial buffer size.

Heater's C version
Code: C
  1. //
  2. // An experiment in doing integer arithmetic using GMP with all numbers represented by strings.
  3. //
  4. // By heater.
  5. // Modified June 11, 2019 to use base 32 strings for intermediate results.
  6. //
  7. #include <stdio.h>
  8. #include <gmp.h>
  9. #include <stdlib.h>
  10. #include <memory.h>
  11. #include <string.h>
  12.  
  13. // Number base used for internal calculations by GMP.
  14. int BASE = 32;
  15.  
  16. mpz_t op1;
  17. mpz_t op2;
  18. mpz_t res;
  19.  
  20. // Functions letis, addis, subis and mulis do large integer arithmetic on integers represented by strings.
  21.  
  22. void writeis(const char *s) {
  23.     mpz_set_str (op1, s, BASE);
  24.     char *buf=mpz_get_str (0, 10, op1);
  25.     puts(buf);
  26.     free(buf);
  27. }
  28.  
  29. char* letis(const char* s) {
  30.     return strdup(s);
  31. }
  32.  
  33. char* addis(const char* s1, const char* s2) {
  34.     mpz_set_str (op1, s1, BASE);
  35.     mpz_set_str (op2, s2, BASE);
  36.     mpz_add (res, op1, op2);  // result = x * y
  37.     char* res_string = mpz_get_str (0, BASE, res);
  38.     return res_string;
  39. }
  40.  
  41. char* subis(const char* s1, const char* s2) {
  42.     mpz_set_str (op1, s1, BASE);
  43.     mpz_set_str (op2, s2, BASE);
  44.     mpz_sub (res, op1, op2);  // result = x * y
  45.     char* res_string = mpz_get_str (0, BASE, res);
  46.     return res_string;
  47. }
  48.  
  49. char* mulis(const char* s1, const char* s2) {
  50.     mpz_set_str (op1, s1, BASE);
  51.     mpz_set_str (op2, s2, BASE);
  52.     mpz_mul (res, op1, op2);  // result = x * y
  53.     char* res_string = mpz_get_str (0, BASE, res);
  54.     return res_string;
  55. }
  56.  
  57. char* memo[3];
  58.  
  59. void init_memo() {
  60.     memo[0] = letis("0");
  61.     memo[1] = letis("1");
  62.     memo[2] = letis("1");
  63. }
  64.  
  65. void clean_memo() {
  66.     free(memo[0]);
  67.     free(memo[1]);
  68.     free(memo[2]);
  69. }
  70.  
  71.  
  72. // Return the n'th Fibonacci number as a decimal string for integer n
  73. char* fibois (int n) {
  74.     char* res;
  75.     if (n <= 2) {
  76.         return letis(memo[n]);
  77.     }
  78.  
  79.     int k = (n / 2);
  80.     char* fk = fibois(k);
  81.     char* fk1 = fibois(k + 1);
  82.     char* a;
  83.     char* b;
  84.     if ((n % 2) == 0) {
  85.         a = addis(fk1, fk1);
  86.         b = subis(a, fk);
  87.         res = mulis(fk, b);
  88.     } else {
  89.         a = mulis(fk, fk);
  90.         b = mulis(fk1, fk1);
  91.         res = addis(a, b);
  92.     }
  93.     free(a);
  94.     free(b);
  95.     free(fk);
  96.     free(fk1);
  97.     return res;
  98. }
  99.  
  100. int main(int argc, char* argv[]) {
  101.     char* f;
  102.     int n = 4784969;               // The first Fibonacci number with a million digits
  103.  
  104.     if (argc >= 2) {
  105.         n = atoi(argv[1]);
  106.     }
  107.  
  108.     mpz_init(op1);
  109.     mpz_init(op2);
  110.     mpz_init(res);
  111.  
  112.     init_memo();
  113.  
  114.     f = fibois(n);
  115.     writeis(f);
  116.     free(f);
  117.  
  118.     clean_memo();
  119.  
  120.     mpz_clear(op1);
  121.     mpz_clear(op2);
  122.     mpz_clear(res);
  123.  
  124.     return (0);
  125. }
  126.  



interface.c
Code: C
  1. /* GMP Extension Module
  2. UXLIBS: -lgmp
  3. */
  4.  
  5. #include <stdio.h>
  6. #include <stdlib.h>
  7. #include <string.h>
  8. #include <unistd.h>
  9. #include <memory.h>
  10. #include <gmp.h>
  11. #include "../../basext.h"
  12.  
  13.  
  14. int BASE = 32;
  15. mpz_t op1;
  16. mpz_t op2;
  17. mpz_t res;
  18. char* memo[3];
  19.  
  20.  
  21. static char* let(const char* s) {
  22.     return strdup(s);
  23. }
  24.  
  25.  
  26. /**************************
  27.  Extension Module Functions
  28. **************************/
  29.  
  30. typedef struct _ModuleObject {
  31.   void *HandleArray;
  32. }ModuleObject,*pModuleObject;
  33.  
  34.  
  35. besVERSION_NEGOTIATE
  36.   return (int)INTERFACE_VERSION;
  37. besEND
  38.  
  39.  
  40. besSUB_START
  41.   pModuleObject p;
  42.  
  43.   besMODULEPOINTER = besALLOC(sizeof(ModuleObject));
  44.   if( besMODULEPOINTER == NULL )return 0;
  45.  
  46.   p = (pModuleObject)besMODULEPOINTER;
  47.   return 0;
  48. besEND
  49.  
  50.  
  51. besSUB_FINISH
  52.   pModuleObject p;
  53.  
  54.   p = (pModuleObject)besMODULEPOINTER;
  55.   if( p == NULL )return 0;
  56.   return 0;
  57. besEND
  58.  
  59.  
  60. /*************
  61.  GMP Functions
  62. *************/
  63.  
  64.  
  65. besFUNCTION(fibo)
  66.   int fval;
  67.  
  68.   besARGUMENTS("i")
  69.     &fval
  70.   besARGEND
  71.  
  72.   char buf[fval+1];
  73.   memset(buf,0,1);
  74.   mpz_t res;
  75.   mpz_init(res);
  76.  
  77.   mpz_fib_ui(res, fval);
  78.  
  79.   gmp_snprintf( buf,sizeof(buf),"%Zd", res );
  80.  
  81.   besRETURN_STRING(buf);
  82. besEND
  83.  
  84.  
  85. besFUNCTION(writeis)
  86.   const char *s;
  87.  
  88.   besARGUMENTS("z")
  89.     &s
  90.   besARGEND
  91.  
  92.   mpz_set_str (op1, s, BASE);
  93.   char *buf=mpz_get_str (0, 10, op1);
  94.   puts(buf);
  95.   free(buf);
  96.   besRETURNVALUE = NULL;
  97.  
  98. besEND
  99.  
  100.  
  101. besFUNCTION(letis)
  102.   const char* s;
  103.  
  104.   besARGUMENTS("z")
  105.     &s
  106.   besARGEND
  107.  
  108.   besRETURN_STRING(strdup(s));
  109.  
  110. besEND
  111.  
  112.  
  113. besFUNCTION(addis)
  114.   const char* s1;
  115.   const char* s2;
  116.  
  117.   besARGUMENTS("zz")
  118.     &s1, &s2
  119.   besARGEND
  120.  
  121.   mpz_set_str (op1, s1, BASE);
  122.   mpz_set_str (op2, s2, BASE);
  123.   mpz_add (res, op1, op2);
  124.   char* res_string = mpz_get_str (0, BASE, res);
  125.   besRETURN_STRING(res_string);
  126.  
  127. besEND
  128.  
  129.  
  130. besFUNCTION(subis)
  131.   const char* s1;
  132.   const char* s2;
  133.  
  134.   besARGUMENTS("zz")
  135.     &s1, &s2
  136.   besARGEND
  137.  
  138.   mpz_set_str (op1, s1, BASE);
  139.   mpz_set_str (op2, s2, BASE);
  140.   mpz_sub (res, op1, op2);
  141.   char* res_string = mpz_get_str (0, BASE, res);
  142.   besRETURN_STRING(res_string);
  143.  
  144. besEND
  145.  
  146.  
  147. besFUNCTION(mulis)
  148.   const char* s1;
  149.   const char* s2;
  150.  
  151.   besARGUMENTS("zz")
  152.     &s1, &s2
  153.   besARGEND
  154.  
  155.   mpz_set_str (op1, s1, BASE);
  156.   mpz_set_str (op2, s2, BASE);
  157.   mpz_mul (res, op1, op2);
  158.   char* res_string = mpz_get_str (0, BASE, res);
  159.   besRETURN_STRING(res_string);
  160.  
  161. besEND
  162.  
  163.  
  164. besFUNCTION(init_memo)
  165.   memo[0] = let("0");
  166.   memo[1] = let("1");
  167.   memo[2] = let("1");
  168.   besRETURNVALUE = NULL;
  169.  
  170. besEND
  171.  
  172. besFUNCTION(clean_memo)
  173.   free(memo[0]);
  174.   free(memo[1]);
  175.   free(memo[2]);
  176.   besRETURNVALUE = NULL;
  177.  
  178. besEND
  179.  
  180.  
  181. besFUNCTION(init_globals)
  182.   mpz_init(op1);
  183.   mpz_init(op2);
  184.   mpz_init(res);
  185.   besRETURNVALUE = NULL;
  186.  
  187. besEND
  188.  
  189.  
  190. besFUNCTION(clear_globals)
  191.   mpz_clear(op1);
  192.   mpz_clear(op2);
  193.   mpz_clear(res);
  194.   besRETURNVALUE = NULL;
  195.  
  196. besEND
  197.  

gmp.bas
Code: ScriptBasic
  1.  
  2. MODULE GMP
  3.  
  4. DECLARE SUB ::FIBO              ALIAS "fibo"            LIB "gmp"
  5. DECLARE SUB ::BI_WRITE          ALIAS "writeis"         LIB "gmp"
  6. DECLARE SUB ::BI_LET            ALIAS "letis"           LIB "gmp"
  7. DECLARE SUB ::BI_ADD            ALIAS "addis"           LIB "gmp"
  8. DECLARE SUB ::BI_SUB            ALIAS "subis"           LIB "gmp"
  9. DECLARE SUB ::BI_MUL            ALIAS "mulis"           LIB "gmp"
  10. DECLARE SUB ::BI_INIT_MEMO      ALIAS "init_memo"       LIB "gmp"
  11. DECLARE SUB ::BI_CLEAN_MEMO     ALIAS "clean_memo"      LIB "gmp"
  12. DECLARE SUB ::BI_INIT_GLOBALS   ALIAS "init_globals"    LIB "gmp"
  13. DECLARE SUB ::BI_CLEAR_GLOBALS  ALIAS "clear_globals"   LIB "gmp"
  14.  
  15. END MODULE
  16.  

gmpfibo.sb
Code: ScriptBasic
  1. ' Return the n'th Fibonacci number as a decimal string for integer n
  2.  
  3. IMPORT gmp.bas
  4.  
  5. SPLITA STRING(4,"0") BY "" To memo
  6.  
  7. FUNCTION fibois(n)
  8.   IF n <= 2 THEN
  9.     fibis = GMP::BI_LET(memo[n])
  10.   END IF
  11.  
  12.   k = (n / 2)
  13.   fk = fibois(k)
  14.   fk1 = fibois(k + 1)
  15.   a = ""
  16.   b = ""
  17.   IF n % 2 = 0 THEN
  18.     a = GMP::BI_ADD(fk1, fk1)
  19.     b = GMP::BI_SUB(a, fk)
  20.     res = GMP::BI_MUL(fk, b)
  21.   ELSE
  22.     a = GMP::BI_MUL(fk, fk)
  23.     b = GMP::BI_MUL(fk1, fk1)
  24.     res = GMP::BI_ADD(a, b)
  25.   END IF
  26.   fibois = res
  27. END FUNCTION
  28.  
  29. ' MAIN
  30.  
  31. n = 4784969
  32. GMP::BI_INIT_GLOBALS()
  33. GMP::BI_INIT_MEMO()
  34. f = fibois(n)
  35. GMP::BI_WRITE(f)
  36. GMP::BI_CLEAN_MEMO()
  37. GMP::BI_CLEAR_GLOBALS()
  38.  
  39. PRINT LEN(f),"\n"
  40.  
« Last Edit: June 11, 2019, 09:46:01 PM by John »