help-octave
[Top][All Lists]

## Re: Arbitrary precision integer arithmetic

 From: Etienne Grossmann Subject: Re: Arbitrary precision integer arithmetic Date: Fri, 07 Sep 2018 22:21:09 -0700 User-agent: Internet Messaging Program (IMP) H3 (4.3.10)

Doug Stewart <address@hidden> (Fri, 7 Sep 2018 13:19:25 -0400) wrote:

On Fri, Sep 7, 2018 at 11:56 AM Etienne Grossmann <address@hidden> wrote:

Etienne Grossmann <address@hidden> (Mon, 03 Sep 2018 17:56:44 -0700) wrote:

Doug Stewart <address@hidden> (Mon, 3 Sep 2018 20:39:59 -0400) wrote:

On Mon, Sep 3, 2018 at 8:26 PM Etienne Grossmann <address@hidden> wrote:

Hi all,

I couldn't find arbitrary-precision integer arithmetic for Octave the other day (I didn't look that hard, though), and I thought it would be fun to write a few functions myself.

I made them available (GPL'd) at https://sourceforge.net/p/octave-intinf/code/ci/master/tree/

Cheers,

Etienne

```The supported operations are:

>> s = a = intinf (126)               % Creation of intinf-class object from integer or string
s = 126
>> b = intinf ('-792716404922304611755050042308687')
b = -792716404922304611755050042308687
>> c = a + b                          % Addition
c = -792716404922304611755050042308561
>> d = a - b                          % Subtraction
d = 792716404922304611755050042308813
>> e = c * d                          % Multiplication
e = -628399298632943207296804593622614025887290777301864028724995648093
>> f = (e - a^2) / b                  % Exponentiation, division
f = 792716404922304611755050042308687
>> f == -b                            % Logical operators
ans =  1
>> [q, r] = intinf_divide (b^2, a)    % Divsion and remainder
q = 4987296020896374661085750743036619253073736327792571656547584634
r = 85
>> r < a
ans =  1
>> x = intinf_rand (100)              % Random number (here, with 100 digits)
x = 5957920324041182863967505832677200714276564790135557044108872634353981714576255029666518574621376764
>> tic, a^a, toc                      % Don't be in a big hurry
ans = 4432907660220782149197257457170010056248664733961715006433455717789043517106373872170818953941792055669609014893218047089803712563472169065833738899530142657476809234058293370126853817068631046152741967763913240019546541793769190722594113575550312228000452759781376
Elapsed time is 2.70419 seconds.

Code organization:

* This directory contains non-class m-files, e.g. intinf_div().
* The sub-directory @intinf contains the class m-files, e.g. intinf(), times().
* The sub-sub-directory @intinf/private contains the back-end functions that
actually perform the calculations, e.g. apia_add(), app_mul(). "apia" stands
for "arbitrary-precision integer arithmetic". The integers are represented as
strings. These functions take double or string arguments.

TODO:

* Don't use a gimmick in intinf() to put the apia_X() functions in the Octave
path.
* Better errors.
* Make a package loadable by octave's pkg function.```

--

Did you try the symbolic pkg?

It has all this I believe.
VPA

[snip]

Hi Doug,

caramba! You're right!

Cheers,

Etienne

octave:9> y = vpa('126^126', 400)

y = (sym)

44329076602207821491972574571700100562486647339617150064334557177890435171063738721708189539417920556696090148932180470898037125634721690658337388995301426
57476809234058293370126853817068631046152741967763913240019546541793769190722594113575550312228000452759781376.00000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000

octave:10> a = intinf(126)
a = 126
octave:11> a^a
ans = 4432907660220782149197257457170010056248664733961715006433455717789043517106373872170818953941792055669609014893218047089803712563472169065833738899530142657476809234058293370126853817068631046152741967763913240019546541793769190722594113575550312228000452759781376
octave:12>

----------------------------------------------------------------
This message was sent using IMP, the Internet Messaging Program.

Hi Doug,

would you know if it is possible to not specify the required precision (400, in the above example) in the symbolic package?

Thanks,

Etienne

----------------------------------------------------------------
This message was sent using IMP, the Internet Messaging Program.

I am not sure what you really want!
--
DAS

Hi again,

sorry, I didn't express myself clearly, and I didn't know about the sym() function, only the vpa() function.

Let's say I want to compute 423^567 exactly.

With the vpa() function, I have to specify the length of the result, as in vpa (vpa(423)^567, 1500).

The correct approach is to use the sym() function:  with the sym() function, sym('423^567') gives me the correct answer.

Thanks,

Etienne

----------------------------------------------------------------
This message was sent using IMP, the Internet Messaging Program.