Recent

Author Topic: [Solved] gmp unit problem (arbitrary precision arithmetic)  (Read 6385 times)

HappyLarry

  • Full Member
  • ***
  • Posts: 144
[Solved] gmp unit problem (arbitrary precision arithmetic)
« on: January 18, 2015, 05:42:31 pm »
I have been using the gmp unit and have a reasonable understanding of it.
I have used the floating point 'sqrt' function and it has worked correctly up until now. However when I try to find root of
1247257156019619631959228850474915324796444794921068388944443347008749568000000000000000

I get the answer
35316528085580830944581390125742259962343219.93...

The correct answer is
35316528085580830944581390125749276564383499.14....

I have set the precision as 4*(number of digits in the original number). I have noticed that the error occurs on the 32nd digit and 32*4 = 128 which may be important, but I don't know why nor how to get around it.

Any ideas would be welcome.
« Last Edit: January 18, 2015, 11:57:21 pm by HappyLarry »
Use Lazarus and Free Pascal and stand on the shoulders of giants . . . very generous giants. Thank you.

HappyLarry

  • Full Member
  • ***
  • Posts: 144
Re: gmp unit problem (arbitrary precision arithmetic)
« Reply #1 on: January 18, 2015, 05:56:04 pm »
Here is the procedure which might help).

procedure sqrt2(var strAnswer:string; var Expon:integer; strNum1: string; NumPlaces:integer);
var
    a,b: MPFloat;
    n:integer;
    p:mp_exp_t;

begin
      f_init (a);
      f_init (b);

      f_set_str (a, '1', 10);
      f_set_str (b, strNum1, 10);

      p:=4*NumPlaces;
      f_set_prec (a, p);
      f_sqrt (a,b);

      strAnswer:=f_get_str(n, 10, p, a);
      Expon:=n;
      f_clear (a);
      f_clear (b);

end;
« Last Edit: January 18, 2015, 06:16:07 pm by HappyLarry »
Use Lazarus and Free Pascal and stand on the shoulders of giants . . . very generous giants. Thank you.

marcov

  • Global Moderator
  • Hero Member
  • *****
  • Posts: 7493
Re: gmp unit problem (arbitrary precision arithmetic)
« Reply #2 on: January 18, 2015, 06:05:20 pm »
Precision maybe binary based (number of bits) instead of decimal ?

HappyLarry

  • Full Member
  • ***
  • Posts: 144
Re: gmp unit problem (arbitrary precision arithmetic)
« Reply #3 on: January 18, 2015, 06:13:50 pm »
It is binary-based but I have assumed that the number of bits needed to represent a number in decimal will be 4*(number of decimal digits).

e.g
any 1 digit number can be represented in 4 bits since 24 = 16
any 2 digit number can be represented in 8 bits since 28 = 256
any 3 digit number can be represented in 12 bits since 212 = 4096
any 4 digit number can be represented in 16 bits since 216 = 65536
etc.

That explains why I have set the precision as 4*(number of digits).
Use Lazarus and Free Pascal and stand on the shoulders of giants . . . very generous giants. Thank you.

marcov

  • Global Moderator
  • Hero Member
  • *****
  • Posts: 7493
Re: gmp unit problem (arbitrary precision arithmetic)
« Reply #4 on: January 18, 2015, 06:22:17 pm »
It is binary-based but I have assumed that the number of bits needed to represent a number in decimal will be 4*(number of decimal digits).

Ok. Then I don't know :-)

(it is 3.322 btw, log 10 / log 2, which comes down to slightly over 106 bits for 32 digits.


HappyLarry

  • Full Member
  • ***
  • Posts: 144
Re: gmp unit problem (arbitrary precision arithmetic)
« Reply #5 on: January 18, 2015, 06:33:29 pm »
I assumed that too little precision was a problem rather than too much!

I have just tried the procedure with
      p:=trunc(3.323 * NumPlaces)+1;
but . . . the answer is the same (wrong).
« Last Edit: January 18, 2015, 06:35:04 pm by HappyLarry »
Use Lazarus and Free Pascal and stand on the shoulders of giants . . . very generous giants. Thank you.

Bart

  • Hero Member
  • *****
  • Posts: 3538
    • Bart en Mariska's Webstek
Re: gmp unit problem (arbitrary precision arithmetic)
« Reply #6 on: January 18, 2015, 07:01:15 pm »
The correct answer is
35316528085580830944581390125749276564383499.14....

35316528085580830944581390125749276564383499.172 is closer to the right answer than the "correct" answer you gave?
35316528085580830944581390125749276564383499.172 ^ 2 = 1247257156019619631959228850474915324796444794943759140060606358369027292979933081630924.685584

35316528085580830944581390125749276564383499.14 ^ 2 = 1247257156019619631959228850474915324796444792683501342583433177915818324931979381510380.7396

Bart

HappyLarry

  • Full Member
  • ***
  • Posts: 144
Re: gmp unit problem (arbitrary precision arithmetic)
« Reply #7 on: January 18, 2015, 08:14:48 pm »
@Bart
You are right. I had checked for the correct answer using bc but copied it out wrongly (bc is easy to use but I can't copy-and-paste off the screen).

The correct answer is (as you say)

35316528085580830944581390125749276564383499.1716787516731376083479

Either way, the answer I get using gmp is not even close. I wondered what I was doing wrong.
Use Lazarus and Free Pascal and stand on the shoulders of giants . . . very generous giants. Thank you.

BeniBela

  • Hero Member
  • *****
  • Posts: 682
    • homepage
Re: gmp unit problem (arbitrary precision arithmetic)
« Reply #8 on: January 18, 2015, 09:56:19 pm »

The correct answer is (as you say)

35316528085580830944581390125749276564383499.1716787516731376083479

Needs more digits:

35316528085580830944581390125749276564383499.171678751673137608347918603494756890875859839579687481657183436994157146340168806365847166873813522

engkin

  • Hero Member
  • *****
  • Posts: 2513
Re: gmp unit problem (arbitrary precision arithmetic)
« Reply #9 on: January 18, 2015, 11:03:44 pm »
I wondered what I was doing wrong.

Two small typos:
1-Perc is number of "bits".
2-In your call to f_get_str you do not pass the number of "digits"
Change your code to:
Code: [Select]
procedure sqrt2(var strAnswer:string; var Expon:integer; strNum1: string; NumPlaces:integer);
var
  a,b: MPFloat;
  n:integer;
  p:mp_exp_t;
  prec: ValUInt;
  digits: ValUInt;
begin
  digits := 200;
  prec := 8*digits;

  f_set_default_prec(prec);
  f_init (a);
  f_init (b);

  f_set_str (a, '1', 10);
  f_set_str (b, strNum1, 10);

  //p:=4*NumPlaces;
  //f_set_prec (a, p);
  f_sqrt (a,b);

  strAnswer:=f_get_str(n, 10, digits, a);
  Expon:=n;

  f_clear (a);
  f_clear (b);
end;

You should get 200 digits:
Quote
35316528085580830944581390125749276564383499.17167875167313760834791860349475689087585983957968748165718343699415714634016880636584716687381352199364066613445449266509253059227517946683412221654516511

HappyLarry

  • Full Member
  • ***
  • Posts: 144
[Solved]:gmp unit problem (arbitrary precision arithmetic)
« Reply #10 on: January 18, 2015, 11:56:18 pm »
This works. Thank you.

Conclusion:
Precision is measured in bits and each digit required needs 8 bits of precision.
A number is stored digit-by-digit like the characters in a string rather than as a binary number like an integer.
Effectively I had half the precision needed.

Thanks to all contributors.

Use Lazarus and Free Pascal and stand on the shoulders of giants . . . very generous giants. Thank you.