Thread: Arbitrary precision modulo operation

Arbitrary precision modulo operation

From
Chadwick Boggs
Date:
I need to perform modulo operations on extremely large numbers.  The %
operator is giving me number out of range errors and the mod(x, y)
function simply seems to return the wrong results.  Also, my numerator
is in the format of a quoted string, which the mod function can't take.

Desparately searching for solutions,
Chadwick.


Re: Arbitrary precision modulo operation

From
Chadwick Boggs
Date:
Bruno, perhaps round is an issue.  Thank you.  Here is an example that
should involve no rounding and indeed it works:

Multiply the ten largest integer scale prime numbers:

# select 2147483477::numeric * 2147483489::numeric * 2147483497::numeric
* 2147483543::numeric * 2147483549::numeric * 2147483563::numeric *
2147483579::numeric * 2147483587::numeric * 2147483629::numeric *
2147483647::numeric;
                                            ?column?
------------------------------------------------------------------------------------------------
 2085923946138988916149190605561960475118165298582929035878182900998428077414994652962618167119
(1 row)

Now, modulo this by any of the factors correctly returns 0:

# select
'2085923946138988916149190605561960475118165298582929035878182900998428077414994652962618167119'::numeric%
2147483563;
 ?column?
----------
        0
(1 row)

Also, diviing out all of the factors correctly returns 1:

# select
2085923946138988916149190605561960475118165298582929035878182900998428077414994652962618167119
/ 2147483477 / 2147483489 / 2147483497 / 2147483543 / 2147483549 /
2147483563 / 2147483579 / 2147483587 / 2147483629 / 2147483647;
        ?column?
------------------------
 1.00000000000000000000
(1 row)

This provides the solution to my problem:  I merely needed to cast all
of the numbers to numberic for product and modulo operations.

Thank you,
Chadwick.


Bruno Wolff III wrote:

>On Mon, Apr 26, 2004 at 10:18:52 -0400,
>  Chadwick Boggs <chadwickboggs@yahoo.com> wrote:
>
>
>>I need to perform modulo operations on extremely large numbers.  The %
>>operator is giving me number out of range errors and the mod(x, y)
>>function simply seems to return the wrong results.  Also, my numerator
>>is in the format of a quoted string, which the mod function can't take.
>>
>>
>
>I tried some large values and seem to be getting reasonable results.
>I did get some negative remainders, but I didn't find anything in the
>mod documentation on whether or not these were allowed. For small values
>I always got positive remainders so it isn't consistantly return the
>remainder closest to zero. My guess is that it has to do rounding when
>dividing numerics.
>
>
>


Re: Arbitrary precision modulo operation

From
Bruno Wolff III
Date:
On Mon, Apr 26, 2004 at 13:30:17 -0400,
  Chadwick Boggs <chadwickboggs@yahoo.com> wrote:
> Example of wrong results from modulo operation of arbitrary precision
> numbers:
>
> # select '123456789012345678901234567890'::numeric % 123;
> ?column?
> ----------
>       -6
> (1 row)
>
> # select mod('123456789012345678901234567890'::numeric, 123);
> mod
> -----
>  -6
> (1 row)
>
> The correct result (at least according to another, unnamed, RDBMS):
>
> > select '123456789012345678901234567890' % 123;
> +----------------------------------------+
> | '123456789012345678901234567890' % 123 |
> +----------------------------------------+
> |                                     58 |
> +----------------------------------------+
> 1 row in set (0.00 sec)

I checked this with bc and I got -6 (and 117) as being correct. I would
think the other database was wrong.

It wouldn't happen to be MYSQL would it?

Re: Arbitrary precision modulo operation

From
Bruno Wolff III
Date:
On Mon, Apr 26, 2004 at 10:18:52 -0400,
  Chadwick Boggs <chadwickboggs@yahoo.com> wrote:
> I need to perform modulo operations on extremely large numbers.  The %
> operator is giving me number out of range errors and the mod(x, y)
> function simply seems to return the wrong results.  Also, my numerator
> is in the format of a quoted string, which the mod function can't take.

I tried some large values and seem to be getting reasonable results.
I did get some negative remainders, but I didn't find anything in the
mod documentation on whether or not these were allowed. For small values
I always got positive remainders so it isn't consistantly return the
remainder closest to zero. My guess is that it has to do rounding when
dividing numerics.

Re: Arbitrary precision modulo operation

From
Chadwick Boggs
Date:
Example of wrong results from modulo operation of arbitrary precision
numbers:

# select '123456789012345678901234567890'::numeric % 123;
 ?column?
----------
       -6
(1 row)

# select mod('123456789012345678901234567890'::numeric, 123);
 mod
-----
  -6
(1 row)

The correct result (at least according to another, unnamed, RDBMS):

 > select '123456789012345678901234567890' % 123;
+----------------------------------------+
| '123456789012345678901234567890' % 123 |
+----------------------------------------+
|                                     58 |
+----------------------------------------+
1 row in set (0.00 sec)


Bruno Wolff III wrote:

>On Mon, Apr 26, 2004 at 10:18:52 -0400,
>  Chadwick Boggs <chadwickboggs@yahoo.com> wrote:
>
>
>>I need to perform modulo operations on extremely large numbers.  The %
>>operator is giving me number out of range errors and the mod(x, y)
>>function simply seems to return the wrong results.  Also, my numerator
>>is in the format of a quoted string, which the mod function can't take.
>>
>>
>
>How large is extremely large?
>You can cast the strings to a numeric type to solve the string problem.
>'numeric' should work for numbers up to about 1000 digits.
>For example:
>area=> select '1234567890'::numeric % '123'::numeric;
> ?column?
>----------
>       39
>(1 row)
>
>If you are getting wrong results you should post a specific example so
>that the developers can figure out what is going wrong.
>
>
>


Re: Arbitrary precision modulo operation

From
Bruno Wolff III
Date:
On Mon, Apr 26, 2004 at 10:18:52 -0400,
  Chadwick Boggs <chadwickboggs@yahoo.com> wrote:
> I need to perform modulo operations on extremely large numbers.  The %
> operator is giving me number out of range errors and the mod(x, y)
> function simply seems to return the wrong results.  Also, my numerator
> is in the format of a quoted string, which the mod function can't take.

How large is extremely large?
You can cast the strings to a numeric type to solve the string problem.
'numeric' should work for numbers up to about 1000 digits.
For example:
area=> select '1234567890'::numeric % '123'::numeric;
 ?column?
----------
       39
(1 row)

If you are getting wrong results you should post a specific example so
that the developers can figure out what is going wrong.

Re: Arbitrary precision modulo operation

From
"Dann Corbit"
Date:
Maple output:
y := 123456789012345678901234567890 mod 123;
                               y := 117

CONNX output (which uses qfloat by S. Moshier):
select mod(123456789012345678901234567890 , 123) {nopassthrough}
117

PariGP output:
? Mod(123456789012345678901234567890, 123)
%4 = Mod(117, 123)
?

> -----Original Message-----
> From: Bruno Wolff III [mailto:bruno@wolff.to]
> Sent: Monday, April 26, 2004 10:42 AM
> To: Chadwick Boggs
> Cc: pgsql-general@postgresql.org
> Subject: Re: [GENERAL] Arbitrary precision modulo operation
>
>
> On Mon, Apr 26, 2004 at 13:30:17 -0400,
>   Chadwick Boggs <chadwickboggs@yahoo.com> wrote:
> > Example of wrong results from modulo operation of arbitrary
> precision
> > numbers:
> >
> > # select '123456789012345678901234567890'::numeric % 123; ?column?
> > ----------
> >       -6
> > (1 row)
> >
> > # select mod('123456789012345678901234567890'::numeric, 123); mod
> > -----
> >  -6
> > (1 row)
> >
> > The correct result (at least according to another, unnamed, RDBMS):
> >
> > > select '123456789012345678901234567890' % 123;
> > +----------------------------------------+
> > | '123456789012345678901234567890' % 123 |
> > +----------------------------------------+
> > |                                     58 |
> > +----------------------------------------+
> > 1 row in set (0.00 sec)
>
> I checked this with bc and I got -6 (and 117) as being
> correct. I would think the other database was wrong.
>
> It wouldn't happen to be MYSQL would it?
>
> ---------------------------(end of
> broadcast)---------------------------
> TIP 1: subscribe and unsubscribe commands go to
> majordomo@postgresql.org
>

Re: Arbitrary precision modulo operation

From
Paul Tillotson
Date:
I see there are a few misconceptions about numeric and modulus on here:

(1) A modulus operation on a numeric type should NOT have rounding
errors.  The whole point of numeric is that it is an arbitrary precision
BASE 10 representation of your number.  The modulus returns the (whole
number) remainder as a result of a division.

(2) the modulus operator/function is, AFAIK, supposed to return the
modulus with the SAME SIGN as the divisor, so I think this is a bug.
That's what every other modulus operator that I have ever seen does.
Would you mind doing

foodb=> SELECT version();

(3) MySQL just rounds large numbers to the highest value that the type
will support, and apparently, no arbitrary precision types are listed on
this page:

http://dev.mysql.com/doc/mysql/en/Numberic_type_overview.html

You can't expect to get a modulus from an out-of-range number.

Paul Tillotson

Chadwick Boggs wrote:

> Example of wrong results from modulo operation of arbitrary precision
> numbers:
>
> # select '123456789012345678901234567890'::numeric % 123;
> ?column?
> ----------
>       -6
> (1 row)
>
> # select mod('123456789012345678901234567890'::numeric, 123);
> mod
> -----
>  -6
> (1 row)
>
> The correct result (at least according to another, unnamed, RDBMS):
>
> > select '123456789012345678901234567890' % 123;
> +----------------------------------------+
> | '123456789012345678901234567890' % 123 |
> +----------------------------------------+
> |                                     58 |
> +----------------------------------------+
> 1 row in set (0.00 sec)
>
>
> Bruno Wolff III wrote:
>
>> On Mon, Apr 26, 2004 at 10:18:52 -0400,
>>  Chadwick Boggs <chadwickboggs@yahoo.com> wrote:
>>
>>
>>> I need to perform modulo operations on extremely large numbers.  The
>>> % operator is giving me number out of range errors and the mod(x, y)
>>> function simply seems to return the wrong results.  Also, my
>>> numerator is in the format of a quoted string, which the mod
>>> function can't take.
>>>
>>
>>
>> How large is extremely large?
>> You can cast the strings to a numeric type to solve the string problem.
>> 'numeric' should work for numbers up to about 1000 digits.
>> For example:
>> area=> select '1234567890'::numeric % '123'::numeric;
>> ?column?
>> ----------
>>       39
>> (1 row)
>>
>> If you are getting wrong results you should post a specific example so
>> that the developers can figure out what is going wrong.
>>
>>
>>
>
>
> ---------------------------(end of broadcast)---------------------------
> TIP 2: you can get off all lists at once with the unregister command
>    (send "unregister YourEmailAddressHere" to majordomo@postgresql.org)
>
>



Re: Arbitrary precision modulo operation

From
"Dann Corbit"
Date:
> -----Original Message-----
> From: Paul Tillotson [mailto:pntil@shentel.net]
> Sent: Monday, April 26, 2004 4:41 PM
> To: pgsql-general@postgresql.org
> Subject: Re: [GENERAL] Arbitrary precision modulo operation
>
>
> I see there are a few misconceptions about numeric and
> modulus on here:
>
> (1) A modulus operation on a numeric type should NOT have
> rounding errors.  The whole point of numeric is that it is an
> arbitrary precision BASE 10 representation of your number.

This is true

> The modulus returns the (whole
> number) remainder as a result of a division.

This is true if the numeric values are integers.

When the values are not integral, some non-integral results can be
returned.

2.50 % 2.50 is 0
But 13.89 modulo 3.50 is 3.39
If you work it out on paper, you will see that 3.39 is the correct
remainder.

> (2) the modulus operator/function is, AFAIK, supposed to
> return the modulus with the SAME SIGN as the divisor, so I
> think this is a bug. That's what every other modulus operator
> that I have ever seen does. Would you mind doing

I would agree that a positive modulus is preferable.  However, the
negative result is also mathematically correct.

> foodb=> SELECT version();
>
> (3) MySQL just rounds large numbers to the highest value that
> the type will support, and apparently, no arbitrary precision
> types are listed on this page:
>
> http://dev.mysql.com/doc/mysql/en/Numberic_type_overview.html
>
> You can't expect to get a modulus from an out-of-range number.

Should it not (therefore) throw an error of some sort?

Re: Arbitrary precision modulo operation

From
Alvaro Herrera
Date:
On Mon, Apr 26, 2004 at 12:48:45PM -0700, Dann Corbit wrote:
> Maple output:
> y := 123456789012345678901234567890 mod 123;
>                                y := 117

PgSQL 7.3.6 gives the right answer (117), 7.4 gets it wrong (-6). Most
likely a bug was introduced when NUMERIC was rewritten.  Strange it
hasn't been noticed before.

--
Alvaro Herrera (<alvherre[a]dcc.uchile.cl>)
"Vivir y dejar de vivir son soluciones imaginarias.
La existencia está en otra parte" (Andre Breton)

Re: Arbitrary precision modulo operation

From
Harald Fuchs
Date:
In article <408D4729.303@yahoo.com>,
Chadwick Boggs <chadwickboggs@yahoo.com> writes:

> Example of wrong results from modulo operation of arbitrary precision
> numbers:

> # select '123456789012345678901234567890'::numeric % 123;
>  ?column?
> ----------
>        -6
> (1 row)

> # select mod('123456789012345678901234567890'::numeric, 123);
>  mod
> -----
>   -6
> (1 row)

> The correct result (at least according to another, unnamed, RDBMS):

>> select '123456789012345678901234567890' % 123;
> +----------------------------------------+
> | '123456789012345678901234567890' % 123 |
> +----------------------------------------+
> |                                     58 |
> +----------------------------------------+
> 1 row in set (0.00 sec)

Is the name of the other, unnamed RDBMS by chance starting with an
'M'?  ;-)

Anyway, both are wrong.  According to GNU MP, the correct answer is
117.  Thus PostgreSQL is at least almost right (117 - 123 = -6).

Re: Arbitrary precision modulo operation

From
Paul Tillotson
Date:
Alvaro Herrera wrote:

>On Mon, Apr 26, 2004 at 12:48:45PM -0700, Dann Corbit wrote:
>
>
>>Maple output:
>>y := 123456789012345678901234567890 mod 123;
>>                               y := 117
>>
>>
>
>PgSQL 7.3.6 gives the right answer (117), 7.4 gets it wrong (-6). Most
>likely a bug was introduced when NUMERIC was rewritten.  Strange it
>hasn't been noticed before.
>
>
>
mod(x, y) is computed as x - trunc(x / y) * y in the mod_var() function
(I think).

However, it appears that the division operator itself is rounding up,
such that the trunc() function  (which ought to round down) does no good
as a round up has already occurred.

Thus, the value of (x / y) is 1 too large, and so  x % y is actually
giving you (x % y) - y, a negative number.  I tried looking at how the
division actually works, but it is over my head at least for the 30
minute perusal.

Regards,
Paul Tillotson

-----------------------------------------------------------------------

[paul@pjt4 paul]$ bc
bc 1.06
Copyright 1991-1994, 1997, 1998, 2000 Free Software Foundation, Inc.
This is free software with ABSOLUTELY NO WARRANTY.
For details type `warranty'.
111111111111111111 / 6
18518518518518518

[paul@pjt4 bin]$ ./psql -U postgres template1
Welcome to psql 7.4.2, the PostgreSQL interactive terminal.

Type:  \copyright for distribution terms
       \h for help with SQL commands
       \? for help on internal slash commands
       \g or terminate with semicolon to execute query
       \q to quit

template1=# select 111111111111111111::numeric / 6;
     ?column?
-------------------
 18518518518518519
(1 row)

template1=# select 111111111111111111 / 6;
     ?column?
-------------------
 18518518518518518
(1 row)

template1=# select version();

version

---------------------------------------------------------------------------------------------------------
 PostgreSQL 7.4.2 on i686-pc-linux-gnu, compiled by GCC gcc (GCC) 3.3.2
20031022 (Red Hat Linux 3.3.2-1)
(1 row)






Re: Arbitrary precision modulo operation

From
Tom Lane
Date:
Paul Tillotson <pntil@shentel.net> writes:
> mod(x, y) is computed as x - trunc(x / y) * y in the mod_var() function

> However, it appears that the division operator itself is rounding up,

This is because div_var rounds its output to the number of fractional
digits specified by select_div_scale, which saith

    /*
     * The result scale of a division isn't specified in any SQL standard.
     * For PostgreSQL we select a result scale that will give at least
     * NUMERIC_MIN_SIG_DIGITS significant digits, so that numeric gives a
     * result no less accurate than float8; but use a scale not less than
     * either input's display scale.
     */

You get the "correct" answer if you force a couple more digits to be
carried by increasing either input's dscale, eg

regression=# select 123456789012345678901234567890.00 % 123;
 ?column?
----------
   117.00
(1 row)

It could be that select_div_scale needs to allow some more slop, or that
that routine is okay but mod_var shouldn't depend on it to select the
intermediate result scale for its internal division.  Comments?

            regards, tom lane

Re: Arbitrary precision modulo operation

From
Bruno Wolff III
Date:
On Wed, Apr 28, 2004 at 02:01:00 -0400,
  Tom Lane <tgl@sss.pgh.pa.us> wrote:
> Paul Tillotson <pntil@shentel.net> writes:
> > mod(x, y) is computed as x - trunc(x / y) * y in the mod_var() function
>
> It could be that select_div_scale needs to allow some more slop, or that
> that routine is okay but mod_var shouldn't depend on it to select the
> intermediate result scale for its internal division.  Comments?

One option would be to define a separate division operator that always
returns an integral value and that is truncated toward 0 and use that
for the mod function. People might find this operator useful in itself.

Another option would be for mod to check if the remainder doesn't have the
same sign as the divisor and if so keeping adding the divisor to it until
it does.

Re: Arbitrary precision modulo operation

From
"Dann Corbit"
Date:

> -----Original Message-----
> From: Bruno Wolff III [mailto:bruno@wolff.to]
> Sent: Wednesday, April 28, 2004 12:05 AM
> To: Tom Lane
> Cc: Paul Tillotson; pgsql-general@postgresql.org
> Subject: Re: [GENERAL] Arbitrary precision modulo operation
>
>
> On Wed, Apr 28, 2004 at 02:01:00 -0400,
>   Tom Lane <tgl@sss.pgh.pa.us> wrote:
> > Paul Tillotson <pntil@shentel.net> writes:
> > > mod(x, y) is computed as x - trunc(x / y) * y in the mod_var()
> > > function
> >
> > It could be that select_div_scale needs to allow some more slop, or
> > that that routine is okay but mod_var shouldn't depend on
> it to select
> > the intermediate result scale for its internal division.  Comments?
>
> One option would be to define a separate division operator
> that always returns an integral value and that is truncated
> toward 0 and use that for the mod function. People might find
> this operator useful in itself.

It will give some wrong results.  The result of mod should be the
remainder after division, which is not always integral in the case of
numeric fixed point or floating point.
Consider the output of this program:

#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include <time.h>
int             main(void)
{
    int             i;
    double          den,
                    other;
    srand((unsigned) time(NULL));
    for (i = 0; i < 10; i++) {
        den = rand() + 1.0;
        other = 1.0 + rand() / (rand() + 1.0);
        printf("%f mod %f is %f\n", other, den, fmod(other, den));
        printf("%f mod %f is %f\n", den, other, fmod(den, other));
    }
    return 0;
}

> Another option would be for mod to check if the remainder
> doesn't have the same sign as the divisor and if so keeping
> adding the divisor to it until it does.

I would suggest computation in sufficient digits of accuracy to get a
correct answer.  If any precision is lost in numeric calculations then
an error of PLOSS should be returned (unless somehow there is a total
loss of precision).

Re: Arbitrary precision modulo operation

From
"Dann Corbit"
Date:
I do not know how division is performed in numeric quantities for
PostgreSQL.

However, if the algorithm is being revised, here is a possible outline:

1.  Convert to double and perform the division to get an estimate
2.  Use Newton's method for division, starting at 30 digits and doubling
the digits for each iteration.
3.  Stop when you have reached 50% excess digits or so
4.  Round to the correct value.

Since the computations use smaller significant digits until the last
iteration, there is typically a large savings compared to performing all
calculations in full precision.

So, if the target is 70 digits you will have:

~15 digits precision after the double precision floating point divide
30 digits after one iteration
60 digits after two iterations
120 digits after three iterations

Then round to the correct length.

See:
http://www.azillionmonkeys.com/qed/adiv.html

For a short discussion of Newton's method for division.

Re: Arbitrary precision modulo operation

From
"Dann Corbit"
Date:
/*
** Hopefully, I am not annoying anyone with all of this.
** This is a brief demonstration of how Newton's division algorithm
works:
*/
#include <stdio.h>
double          divide(double x, double y)
{
    unsigned __int64 est;
    double          y1 = 1.0 / y;
    double          y2;
    puts("\nForming reciprocal:");
    est = y1 * 10000;
    y1 = est / 10000.0;         /* make an estimate good to 4 decimal
places... */
    printf("x=%.7f, y=%.7f, y1=%.7f\n", x, y, y1);
    y1 *= (2 - y * y1); /* should be 8 places */
    printf("x=%.14f, y=%.14f, y1=%.14f\n", x, y, y1);
    y1 *= (2 - y * y1); /* In theory ~16 digits */
    printf("x=%.20f, y=%.20f, y1=%.20f\n", x, y, y1);
    y1 *= (2 - y * y1); /* One more is really needed here... */
    printf("x=%.20f, y=%.20f, y1=%.20f\n", x, y, y1);
    puts("\ndividing:");
    printf("result of division is:%f\n", x * y1);
    return 0;
}

int             main(void)
{
    divide(1.0, 2.0);
    divide(1.0, 17.0);
    divide(3.0, 19.0);
    divide(19.0, 3.0);
    return 0;
}

Re: Arbitrary precision modulo operation

From
Alvaro Herrera
Date:
On Wed, Apr 28, 2004 at 03:06:25PM -0700, Dann Corbit wrote:
> /*
> ** Hopefully, I am not annoying anyone with all of this.
> ** This is a brief demonstration of how Newton's division algorithm
> works:
> */

Well, have you looked at Postgres implementation?  I don't know if it's
Newton or not, but it certainly is not simple (I didn't take the time to
understand it).

Have a look at src/backend/utils/adt/numeric.c function div_var()

 Cc: list trimmed ...

--
Alvaro Herrera (<alvherre[a]dcc.uchile.cl>)
"¿Cómo puedes confiar en algo que pagas y que no ves,
y no confiar en algo que te dan y te lo muestran?" (Germán Poo)

Re: Arbitrary precision modulo operation

From
Bruno Wolff III
Date:
On Wed, Apr 28, 2004 at 14:02:57 -0700,
  Dann Corbit <DCorbit@connx.com> wrote:
>
>
> > -----Original Message-----
> > From: Bruno Wolff III [mailto:bruno@wolff.to]
> >
> > One option would be to define a separate division operator
> > that always returns an integral value and that is truncated
> > toward 0 and use that for the mod function. People might find
> > this operator useful in itself.
>
> It will give some wrong results.  The result of mod should be the
> remainder after division, which is not always integral in the case of
> numeric fixed point or floating point.
> Consider the output of this program:

The remainder may not be integral but the quotient should be. However the
idea of getting the quotient, multiplying by the divisor and subtracting
from the dividend is not very good from an efficiancy point of view.

Re: Arbitrary precision modulo operation

From
Tom Lane
Date:
"Dann Corbit" <DCorbit@connx.com> writes:
> I would suggest computation in sufficient digits of accuracy to get a
> correct answer.

Fine.  How many is that, exactly?

            regards, tom lane

Re: Arbitrary precision modulo operation

From
"Dann Corbit"
Date:
> -----Original Message-----
> From: Tom Lane [mailto:tgl@sss.pgh.pa.us]
> Sent: Wednesday, April 28, 2004 9:14 PM
> To: Dann Corbit
> Cc: Bruno Wolff III; Paul Tillotson; pgsql-general@postgresql.org
> Subject: Re: [GENERAL] Arbitrary precision modulo operation
>
>
> "Dann Corbit" <DCorbit@connx.com> writes:
> > I would suggest computation in sufficient digits of
> accuracy to get a
> > correct answer.
>
> Fine.  How many is that, exactly?

Here is what I would suggest:

Using the outline I proposed before (starting with a floating point
divide of DBL_DIG digits of precision), keep doubling the precision
until the precision is 5 digits larger than either operand.  If the last
doubling makes the precision larger (quite likely) simply reduce it to
the smaller margin.

Something like this (pseudocode):

numeric          divide(numeric x, numeric y)
{
    /* double starting estimate of quotient*/
    numeric y1 = 1.0 / y::double;
    numeric two = 2:numeric;
    /* Need to collect the maximal precision of either operand */
    total_precision = get_max_precision(x,y);

    y1 *= (two - y * y1); /* use numeric 30 math or total precision+5
whichever is less */
    if (total_precision <= 25) return x*y; /* correctly rounded -->
banker's rounding? */
    y1 *= (two - y * y1); /* use numeric 60 math or total precision+5
whichever is less */
    if (total_precision <= 50) return x*y; /* correctly rounded -->
banker's rounding? */
    y1 *= (two - y * y1); /* use numeric 120 math or total precision+5
whichever is less */
    if (total_precision <= 110) return x*y; /* correctly rounded -->
banker's rounding? */
    y1 *= (two - y * y1); /* use numeric 240 math or total precision+5
whichever is less */
    if (total_precision <= 230) return x*y; /* correctly rounded -->
banker's rounding? */
    y1 *= (two - y * y1); /* use numeric 480 math or total precision+5
whichever is less */
    if (total_precision <= 470) return x*y; /* correctly rounded -->
banker's rounding? */
    y1 *= (two - y * y1); /* use numeric 960 math or total precision+5
whichever is less */
    if (total_precision <= 950) return x*y; /* correctly rounded -->
banker's rounding? */
    y1 *= (two - y * y1); /* use maximum precision math or total
precision+5 whichever is less */
    return x*y;
}

Re: Arbitrary precision modulo operation

From
"Dann Corbit"
Date:
> -----Original Message-----
> From: Tom Lane [mailto:tgl@sss.pgh.pa.us]
> Sent: Thursday, April 29, 2004 2:31 PM
> To: Dann Corbit
> Cc: Bruno Wolff III; Paul Tillotson; pgsql-general@postgresql.org
> Subject: Re: [GENERAL] Arbitrary precision modulo operation
>
>
> "Dann Corbit" <DCorbit@connx.com> writes:
> > From: Tom Lane [mailto:tgl@sss.pgh.pa.us]
> >> Fine.  How many is that, exactly?
>
> > Here is what I would suggest:
> > Using the outline I proposed before (starting with a floating point
> > divide of DBL_DIG digits of precision), keep doubling the precision
> > until the precision is 5 digits larger than either operand.
>  If the last
> > doubling makes the precision larger (quite likely) simply
> reduce it to
> > the smaller margin.
>
> And this guarantees a correct answer why?
>
> AFAIK div_var is already correct per its spec, which is that it
> generates an answer rounded to the requested number of digits.
> The question at hand is what number of digits to request.
>
> After thinking about it I don't see any reason that DBL_DIG has
> anything to do with a non-surprising answer ... much less
> "DBL_DIG + 5"
> which seems picked out of the air ...

DBL_DIG is there (from float.h) to define the precision of the first
operation.  After that, the precision doubles each time.

5 was pretty much picked out of the air (it would be good enough to
ensure correctness but is overkill).  I assumed that you could lose up
to one ULP per operation.

Do it any way you like.  It was just a suggestion.  It is clear that mod
is not working correctly now, and I was trying to be helpful.

Re: Arbitrary precision modulo operation

From
Tom Lane
Date:
"Dann Corbit" <DCorbit@connx.com> writes:
> From: Tom Lane [mailto:tgl@sss.pgh.pa.us]
>> Fine.  How many is that, exactly?

> Here is what I would suggest:
> Using the outline I proposed before (starting with a floating point
> divide of DBL_DIG digits of precision), keep doubling the precision
> until the precision is 5 digits larger than either operand.  If the last
> doubling makes the precision larger (quite likely) simply reduce it to
> the smaller margin.

And this guarantees a correct answer why?

AFAIK div_var is already correct per its spec, which is that it
generates an answer rounded to the requested number of digits.
The question at hand is what number of digits to request.

After thinking about it I don't see any reason that DBL_DIG has
anything to do with a non-surprising answer ... much less "DBL_DIG + 5"
which seems picked out of the air ...

            regards, tom lane