|
|
| |
Pari(3) |
User Contributed Perl Documentation |
Pari(3) |
Math::Pari - Perl interface to PARI.
use Math::Pari;
$a = PARI 2;
print $a**10000;
or
use Math::Pari qw(Mod);
$a = Mod(3,5);
print $a**10000;
This package is a Perl interface to famous library PARI for
numerical/scientific/number-theoretic calculations. It allows use of most PARI
functions as Perl functions, and (almost) seamless merging of PARI and Perl
data. In what follows we suppose prior knowledge of what PARI is (see
<ftp://megrez.math.u-bordeaux.fr/pub/pari>, or Math::libPARI).
- DEFAULT
- By default the package exports functions PARI(), PARIcol(),
PARIvar(), PARImat(), PARImat_tr() and
parse_as_gp() which convert their argument(s) to a PARI object. (In
fact PARI() is just an alias for "new
Math::Pari"). The function PARI() accepts following
data as its arguments
- One integer
- Is converted to a PARI integer.
- One float
- Is converted to a PARI float.
- One string
- Is executed as a PARI expression (so should not contain whitespace).
- PARI object
- Is passed unchanged.
- Reference to a Perl array
- Each element is converted using the same rules, PARI vector-row with these
elements is returned.
- Several of above
- The same as with a reference to array.
- Conflicts of rules in PARI()
- In deciding what rule of the above to apply the preference is given to the
uppermost choice of those available now. If none matches, then the
string rule is used. So PARI(1) returns integer,
"PARI(1.)" returns float,
"PARI("1")" evaluates
1 as a PARI expression (well, the result is the
same as PARI(1), only slower).
Note that for Perl these data are synonymous, since Perl
freely converts between integers, float and strings. However, to
PARI() only what the argument is now is important. If
$v is 1 in the Perl
world, "PARI($v)" may convert it to an
integer, float, or to the result of evaluating the PARI program
1 (all depending on how
$v was created and accessed in Perl).
This is a fundamental limitation of creating an interface
between two systems, both with polymorphic objects, but with subtly
different semantic of the flavors of these objects. In reality, however,
this is rarely a problem.
- PARIcol(), PARImat() and PARImat_tr()
- PARIcol() behaves in the same way as PARI() unless given
several arguments. In the latter case it returns a vector-column instead
of a vector-row.
PARImat() constructs a matrix out of the given
arguments. It will work if PARI() will construct a vector of
vectors given the same arguments. The internal vectors become columns of
the matrix. PARImat_tr() behaves similarly, but the internal
vectors become rows of the matrix.
Since PARI matrices are similar to vector-rows of
vector-columns, PARImat() is quicker, but PARImat_tr()
better corresponds to the PARI input and output forms of matrices:
print PARImat [[1,2], [3,4]]; # prints [1,3;2,4]
print PARImat_tr [[1,2], [3,4]]; # prints [1,2;3,4]
- parse_as_gp()
- Did you notice that when taking a string, PARI() requires that
there is no whitespace there (outside of string constants)? This is
exactly as the "PARI" library parses
strings. However, to simplify human interaction, the
"gp" calculator allows whitespace,
comments, breaking into multiple lines, many independent expressions (such
as function definitions).
We do not include the corresponding C code from the
calculator, but provide a Perl clone. It supports whitespace,
"\\"-comments, and, for multi-line
arguments, it supports trailing "\"
for line-continuation, trailing binary ops, comma, opening
parenthesis/bracket indicate lines with continuation, group of lines in
"{}" joined into one line.
Keep in mind that this is just a convenience function, and no
attempt was performed to make it particularly quick. Moreover, the PARI
user functions (or maybe it is better to call them user macros?) are
currently not automatically importable into Perl, so to access functions
defined in parse_as_gp()' argument may be awkward. (The temporary
fix is to use a temporary convenience function
__wrap_PARI_macro():
parse_as_gp <<EOP;
add2(x) = x + 2
EOP
*add2 = Math::Pari::__wrap_PARI_macro 'add2';
print add2(17);
but keep in mind that the generated this way wrapper is also
not designed to be quick.)
- "use" with arguments
- If arguments are specified in the "use
Math::Pari" directive, the PARI functions appearing as
arguments are exported in the caller context. In this case the function
PARI() and friends is not exported, so if you need them, you should
include them into export list explicitly, or include
":DEFAULT" tag:
use Math::Pari qw(factorint PARI);
use Math::Pari qw(:DEFAULT factorint);
or simply do it in two steps
use Math::Pari;
use Math::Pari 'factorint';
The other tags recognized are
":PARI",
":all",
"prec=NUMBER", number tags (e.g.,
":4"), overloaded constants tags
(":int",
":float",
":hex") and section names tags. The
number tags export functions from the PARI library from the given class
(except for ":PARI", which exports all
of the classes). Tag ":all" exports
all of the exportable symbols and
":PARI".
Giving "?" command to
"gp" (PARI calculator) lists
the following classes:
1: Standard monadic or dyadic OPERATORS
2: CONVERSIONS and similar elementary functions
3: TRANSCENDENTAL functions
4: NUMBER THEORETICAL functions
5: Functions related to ELLIPTIC CURVES
6: Functions related to general NUMBER FIELDS
7: POLYNOMIALS and power series
8: Vectors, matrices, LINEAR ALGEBRA and sets
9: SUMS, products, integrals and similar functions
10: GRAPHIC functions
11: PROGRAMMING under GP
One can use section names instead of number tags. Recognized
names are
:standard :conversions :transcendental :number :elliptic
:fields :polynomials :vectors :sums :graphic :programming
One can get the list of all of the functions accessible by
"Math::Pari", or the accessible
functions from the given section using listPari() function.
Starting from version 5.005 of Perl, three constant-overload
tags are supported: ":int",
":float",
":hex". If used, all the
integer/float/hex-or-octal-or-binary literals in Perl will be
automatically converted to became PARI objects. For example,
use Math::Pari ':int';
print 2**1000;
is equivalent to
print PARI(2)**PARI(1000);
(The support for this Perl feature is buggy before the Perl
version 5.005_57 - unless Perl uses mymalloc options; you can check for
this with "perl
-V:usemymalloc".) Note also that (at least
with some versions of Perl) one should enable
':float' for conversion of long integer literals
(Perl may consider them as floats, since they won't fit into Perl
integers); note that it is PARI which determines which PARI subtype is
assigned to each such literal:
use Math::Pari ':float', 'type_name';
print type_name 22222222222222222222222;
prints "t_INT".
This package supports all the functions from the PARI library with a
signature which can be recognized by Math::Pari. This means that when
you update the PARI library, the newly added functions will we available
without any change to this package; only a recompile is needed. In fact no
recompile will be needed if you link libPARI dynamically (you need to modify
the Makefile manually to do this).
You can "reach" unsupported functions via going directly
to PARI parser using the string flavor of PARI() function, as in
3 + PARI('O(x^17)');
For some "unreachable" functions there is a special
wrapper functions, such as
"O(variable,power)").
The following functions are specific to GP calculator, thus are
not available to Math::Pari in any way:
default error extern input print print1 printp printp1
printtex quit read system whatnow write write1 writetex
whatnow() function is useless, since Math::Pari does not
support the "compatibility" mode (with older PARI library). The
functionality of print(), write() and variants is available
via automatic string translation, and pari_print() function and its
variants (see "Printout functions").
default() is the only important function with functionality
not supported by the current interface. Note however, that four most
important default() actions are supported by allocatemem(),
setprimelimit(), setprecision() and
setseriesprecision() functions. (When called without arguments, these
functions return the current values.)
allocatemem($bytes) should not be called from inside Math::Pari
functions (such as forprimes()).
Arguments to PARI functions are automatically converted to
"long" or a PARI object depending on the
signature of the actual library function. The arguments are forced into
the given type, so even if "gp" rejects your
code similar to
func(2.5); # func() takes a long in C
arguing that a particular argument should be of
"type T_INT" (i.e., a Pari integer), the
corresponding code will work in
"Math::Pari", since 2.5 is silently
converted to "long", per the function
signature.
PARI functions return a PARI object or a Perl's integer depending on what the
actual library function returns.
Some PARI functions are available in "gp"
(i.e., in "PARI" calculator) via infix
notation only. In "Math::Pari" these
functions are available in functional notations too. Some other convenience
functions are also made available.
- Infix, prefix and postfix operations
- are available under names
gneg, gadd, gsub, gmul, gdiv, gdivent, gmod, gpui,
gle, gge, glt, ggt, geq, gne, gegal, gor, gand,
gcmp, gcmp0, gcmp1, gcmp_1.
"gdivent" means euclidean
quotient, "gpui" is power,
"gegal" checks whether two objects are
equal, "gcmp" is applicable to two
real numbers only, "gcmp0",
"gcmp1",
"gcmp_1" compare with 0, 1 and -1
correspondingly (see PARI user manual for details, or Math::libPARI).
Note that all these functions are more readily available via operator
overloading, so instead of
gadd(gneg($x), $y)
one can write
-$x+$y
(as far as overloading may be triggered, see overload, so we
assume that at least one of $x or
$y is a PARI object).
- Conversion functions
-
pari2iv, pari2nv, pari2num, pari2pv, pari2bool
convert a PARI object to an integer, float, integer/float
(whatever is better), string, and a boolean value correspondingly. Most
the time you do not need these functions due to automatic
conversions.
- Printout functions
-
pari_print, pari_pprint, pari_texprint
perform the same conversions to strings as their PARI
counterparts, but do not print the result. The difference of
pari_print() with pari2pv() is the number of significant
digits they output, and whitespace in the output. pari2pv(),
which is intended for "computer-readable strings", outputs as
many digits as is supported by the current precision of the number;
while pari_print(), which targets human-readable strings, takes
into account the currently specified output precision too.
- Constant functions
- Some mathematical constants appear as function without arguments in PARI.
These functions are available in Math::Pari too. If you export them as in
use Math::Pari qw(:DEFAULT Pi I Euler);
they can be used as barewords in your program:
$x = Pi ** Euler;
- Low-level functions
- For convenience of low-level PARI programmers some low-level functions are
made available as well (all except type_name() and
changevalue() are not exportable):
typ($x)
lg($x)
lgef($x)
lgefint($x)
longword($x, $n)
type_name($x)
changevalue($name,$newvalue)
Here longword($x,$n) returns $n-th
word in the memory representation of $x
(including non-code words). type_name() differs from the PARI
function type(): type() returns a PARI object, while
type_name() returns a Perl string. (PARI objects of string type
behave very non-intuitive w.r.t. string comparison functions; remember
that they are compared using lex() to the results of
evaluation of other argument of comparison!)
The function listPari($number) outputs a list of names of PARI
functions in the section $number. Use
listPari(-1) to get the list across all of the sections.
- Uncompatible functions
-
O
Since implementing "O(7**6)"
would be very tedious, we provide a two-argument form
"O(7,6)" instead (meaning the same as
"O(7^6)" in PARI). Note that with
polynomials there is no problem like this one, both
"O($x,6)" and
"O($x**6)" work.
ifact(n)
integer factorial functions, available from
"gp" as
"n!".
PARI has a big collection of functions which loops over some set. Such a
function takes two special arguments: loop variable, and the code to
execute in the loop.
The code can be either a string (which contains PARI code to
execute - thus should not contain whitespace), or a Perl code reference. The
loop variable can be a string giving the name of PARI variable (as in
fordiv(28, 'j', 'a=a+j+j^2');
or
$j= 'j';
fordiv(28, $j, 'a=a+j+j^2');
), a PARI monomial (as in
$j = PARI 'j';
fordiv(28, $j, sub { $a += $j + $j**2 });
), or a "delayed Math::Pari variable" (as in
$j = PARIvar 'j';
fordiv(28, $j, 'a=a+j+j^2');
). If none of these applies, as in
my $j; # Have this in a separate statement
fordiv(28, $j, sub { $a += $j + $j**2 });
then during the execution of the
"sub", Math::Pari would autogenerate a
PARI variable, and would put its value in $j; this
value of $j is temporary only, the old contents of
$j is restored when fordiv() returns.
Note that since you have no control over this name, you will not
be able to use this variable from your PARI code; e.g.,
$j = 7.8;
fordiv(28, $j, 'a=a+j+j^2');
will not make "j" mirror
$j (unless you explicitly set up
"j" to be a no-argument PARI function
mirroring $j, see "Accessing Perl functions
from PARI code").
Caveats. There are 2 flavors of the "code"
arguments (string/"sub"), and 4 types of
the "variable" arguments
(string/monomial/"PARIvar"/other).
However, not all 8 combinations make sense. As we already explained, an
"other" variable cannot work with a "string" code.
Useless musing alert! Do not read the rest of this section!
Do not use "string" variables with
"sub" code, and do not ask why!
Additionally, the following code will not do what you expect
$x = 0;
$j = PARI 'j';
fordiv(28, 'j', sub { $x += $j } ); # Use $j as a loop variable!
since the PARI function "fordiv"
localizes the PARI variable "j"
inside the loop, but $j will still reference the old
value; the old value is a monomial, not the index of the loop (which is an
integer each time "sub" is called). The
simplest workaround is not to use the above syntax (i.e., not mixing literal
loop variable with Perl loop code, just using $j as
the second argument to "fordiv" is
enough):
$x = 0;
$j = PARI 'j';
fordiv(28, $j, sub { $x += $j } );
Alternately, one can make a delayed variable
$j which will always reference the same thing
"j" references in PARI now by using
"PARIvar" constructor
$x = 0;
$j = PARIvar 'j';
fordiv(28, 'j', sub { $x += $j } );
(This problem is similar to
$ref = \$_; # $$ref is going to be old value even after
# localizing $_ in Perl's grep/map
not accessing localized values of $_ in
the plain Perl.)
Another possible quirk is that
fordiv(28, my $j, sub { $a += $j + $j**2 });
will not work too - by a different reason.
"my" declarations change the
meaning of $j only after the end of
the current statement; thus $j inside
"sub" will access a different
variable $j (typically a non-lexical, global
variable $j) than one you declared on this line.
Use the same name inside PARI code:
sub counter { $i += shift; }
$i = 145;
PARI 'k=5' ;
fordiv(28, 'j', 'k=k+counter(j)');
print PARI('k'), "\n";
prints
984
Due to a difference in the semantic of
variable-number-of-parameters-functions between PARI and Perl, if the Perl
subroutine takes a variable number of arguments (via
"@" in the prototype or a missing
prototype), up to 6 arguments are supported when this function is called
from PARI. If called from PARI with fewer arguments, the rest of arguments
will be set to be integers "PARI 0".
Note also that no direct import of Perl variables is available yet
(but you can write a function wrapper for this):
sub getv () {$v}
There is an unsupported (and undocumented ;-) function for
explicitly importing Perl functions into PARI, possibly with a different
name, and possibly with explicitly specifying number of arguments.
Functions from PARI library may take as arguments and/or return values the
objects of C type "GEN". In Perl these data
are encapsulated into special kind of Perl variables: PARI objects. You can
check for a variable $obj to be a PARI object using
ref $obj and $obj->isa('Math::Pari');
Most the time you do not need this due to automatic conversions
and overloading.
If very lazy, one can code in Perl the same way one does it in PARI. Variables
in PARI are denoted by barewords, as in "x",
and in the default configuration (no warnings, no strict) Perl allows the same
- up to some extent. Do not do this, since there are many surprising problems.
Some bareletters denote Perl operators, like
"q",
"x",
"y",
"s". This can lead to errors in Perl
parsing your expression. E.g.,
print sin(tan(t))-tan(sin(t))-asin(atan(t))+atan(asin(t));
may parse OK after "use Math::Pari qw(sin
tan asin atan)". Why?
After importing, the word "sin"
will denote the PARI function sin(), not Perl operator sin().
The difference is subtle: the PARI function implicitly forces its
arguments to be converted PARI objects; it gets 't'
as the argument, which is a string, thus is converted to what
"t" denotes in PARI - a monomial. While
the Perl operator sin() grants overloading (i.e., it will call PARI
function sin() if the argument is a PARI object), it does not
force its argument; given 't' as argument, it
converts it to what sin() understands, a float (producing
0.), so will give 0. as the
answer.
However
print sin(tan(y))-tan(sin(y))-asin(atan(y))+atan(asin(y));
would not compile. You should avoid lower-case barewords used as
PARI variables, e.g., do
$y = PARI 'y';
print sin(tan($y))-tan(sin($y))-asin(atan($y))+atan(asin($y));
to get
-1/18*y^9+26/4725*y^11-41/1296*y^13+328721/16372125*y^15+O(y^16)
(BTW, it is a very good exercise to get the leading term by
hand).
Well, the same advice again: do not use barewords anywhere in your
program!
Whenever an arithmetic operation includes at least one PARI object, the other
arguments are converted to a PARI object and the corresponding PARI library
functions is used to implement the operation. Currently the following
arithmetic operations are overloaded:
unary -
+ - * / % ** abs cos sin exp log sqrt
<< >>
<= == => < > != <=>
le eq ge lt gt ne cmp
| & ^ ~
Numeric comparison operations are converted to
"gcmp" and friends, string comparisons
compare in lexicographical order using
"lex".
Additionally, whenever a PARI object appears in a situation that
requires integer, numeric, boolean or string data, it is converted to the
corresponding type. Boolean conversion is subject to usual PARI pitfalls
related to imprecise zeros (see documentation of
"gcmp0" in PARI reference).
For details on overloading, see overload.
Note that a check for equality is subject to same pitfalls as in
PARI due to imprecise values. PARI may also refuse to compare data of
different types for equality if it thinks this may lead to counterintuitive
results.
Note also that in PARI the numeric ordering is not defined for
some types of PARI objects. For string comparison operations we use
PARI-lexicographical ordering.
In the versions of perl earlier than 5.003 overloading used a different
interface, so you may need to convert "use
overload" line to %OVERLOAD, or, better,
upgrade.
Starting from version 2.0, this module comes without a PARI library included.
For the source of PARI library see
<ftp://megrez.math.u-bordeaux.fr/pub/pari>.
Note that the PARI notations should be used in the string arguments to
PARI() function, while the Perl notations should be used otherwise.
- "^"
- Power is denoted by "**" in Perl.
- "\" and "\/"
- There are no such operators in Perl, use the word forms
"gdivent(x,y)" and
"gdivround(x,y)" instead.
- "~"
- There is no postfix "~" Perl operator.
Use mattranspose() instead.
- "'"
- There is no postfix "'" Perl operator.
Use deriv() instead.
- "!"
- There is no postfix "!" Perl operator.
Use factorial()/ifact() instead (returning a real or an
integer correspondingly).
- big integers
- Perl converts big literal integers to doubles if they could not be
put into C integers (the particular flavor can be found in the
output of "perl -V" in newer version of
Perl, look for
"ivtype"/"ivsize").
If you want to input such an integer, use
while ($x < PARI('12345678901234567890')) ...
instead of
while ($x < 12345678901234567890) ...
Why? Because conversion to double leads to precision loss
(typically above 1e15, see perlnumber), and you will get something like
12345678901234567168 otherwise.
Starting from version 5.005 of Perl, if the tag
":int" is used on the 'use Math::Pari'
line, all of the integer literals in Perl will be automatically
converted to became PARI objects. E.g.,
use Math::Pari ':int';
print 2**1000;
is equivalent to
print PARI(2)**PARI(1000);
Similarly, large integer literals do not lose precision.
This directive is lexically scoped. There is a similar tag
":hex" which affects hexadecimal,
octal and binary constants. One may also need to use tag
":float" for auto-conversion of large
integer literals which Perl considers as floating point literals (see
""use" with arguments" for
details).
- doubles
- Doubles in Perl are typically of precision approximately 15 digits (see
perlnumber). When you use them as arguments to PARI functions, they are
converted to PARI real variables, and due to intermediate
15-digits-to-binary conversion of Perl variables the result may be
different than with the PARI many-digits-to-binary conversion. E.g.,
"PARI(0.01)" and
"PARI('0.01')" differ at 19-th place, as
setprecision(38);
print pari_print(0.01), "\n",
pari_print('0.01'), "\n";
shows.
Note that setprecision() changes the output format of
pari_print() and friends, as well as the default internal
precision. The generic PARI===>string conversion does not take into
account the output format, thus
setprecision(38);
print PARI(0.01), "\n",
PARI('0.01'), "\n",
pari_print(0.01), "\n";
will print all the lines with different number of digits after
the point: the first one with 22, since the double 0.01 was converted to
a low-precision PARI object, the second one with 41, since internal form
for precision 38 requires that many digits for representation, and the
last one with 39 to have 38 significant digits.
Starting from version 5.005 of Perl, if the tag
":float" is used on the
"use Math::Pari" line, all the float
literals in Perl will be automatically converted to became PARI objects.
E.g.,
use Math::Pari ':float';
print atan(1.);
is equivalent to
print atan(PARI('1.'));
Similarly, large float literals do not lose precision.
This directive is lexically scoped.
- array base
- Arrays are 1-based in PARI, are 0-based in Perl. So while array access is
possible in Perl, you need to use different indices:
$nf = PARI 'nf'; # assume that PARI variable nf contains a number field
$a = PARI('nf[7]');
$b = $nf->[6];
Now $a and $b
contain the same value.
- matrices
- Note that "PARImat([[...],...,[...])"
constructor creates a matrix with specified columns, while in PARI the
command "[1,2,3;4,5,6]" creates a matrix
with specified rows. Use a convenience function PARImat_tr() which
will transpose a matrix created by PARImat() to use the same order
of elements as in PARI.
- builtin perl functions
- Some PARI functions, like "length" and
"eval", are Perl (semi-)reserved words.
To reach these functions, one should either import them:
use Math::Pari qw(length eval);
or call them with prefix (like
&length) or the full name (like
"Math::Pari::length").
If you have Term::Gnuplot Perl module installed, you may use high-resolution
graphic primitives of PARI. Before the usage you need to establish a
link between Math::Pari and Term::Gnuplot by calling link_gnuplot().
You can change the output filehandle by calling set_plot_fh(), and
output terminal by calling plotterm(), as in
use Math::Pari qw(:graphic asin);
open FH, '>out.tex' or die;
link_gnuplot(); # automatically loads Term::Gnuplot
set_plot_fh(\*FH);
plotterm('emtex');
ploth($x, .5, .999, sub {asin $x});
close FH or die;
libPARI documentation is included, see Math::libPARI. It is converted from
Chapter 3 of PARI/GP documentation by the gphelp script of
GP/PARI.
No environment variables are used.
- A few of PARI functions are available indirectly only.
- Using overloading constants with the Perl versions below 5.005_57 could
lead to segfaults (at least without "-D
usemymalloc"), as in:
use Math::Pari ':int';
for ( $i = 0; $i < 10 ; $i++ ) { print "$i\n" }
- It may be possible that conversion of a Perl value which has both the
integer slot and the floating slot set may create a PARI integer, even if
the actual value is not an integer.
- problems with refcounting of array elements and Mod().
Workaround: make the modulus live longer than the result of
Mod(). Until Perl version 5.6.1, one
should exercise a special care so that the modulus goes out of scope on
a different statement than the result:
{ my $modulus = 125;
{ my $res = Mod(34, $modulus);
print $res;
}
$fake = 1; # A (fake) statement here is required
}
Here $res is destructed before the
"$fake = 1" statement,
$modulus is destructed before the first
statement after the provided block. However, if you remove the
"$fake = 1" statement, both these
variables are destructed on the first statement after the provided block
(and in a wrong order!).
In 5.6.1 declaring
$modulus before $res is
all that is needed to circumvent the same problem:
{ my $modulus = 125;
my $res = Mod(34, $modulus);
print $res;
} # destruction will happen in a correct order.
Access to array elements may result in similar problems. Hard
to fix since in PARI the data is not refcounted.
- Legacy implementations of dynalinking require the code of DLL to be
compiled to be "position independent" code (PIC). This slows
down the execution, while allowing sharing the loaded copy of the DLL
between different processes. [On contemporary architectures the same
effect is allowed without the position-independent hack.]
Currently, PARI assembler files are not position-independent.
When compiled for the dynamic linking on legacy systems, this creates a
DLL which cannot be shared between processes. Some legacy systems are
reported to recognize this situation, and load the DLL as a non-shared
module. However, there may be systems (are there?) on which this can
cause some "problems".
Summary: if the dynaloading on your system requires some kind
of "-fPIC" flag, using
"assembler" compiles (anything but
"machine=none") *may* force you to do
a static build (i.e., creation of a custom Perl executable with
perl Makefile.PL static
make perl
make test_static
).
- isprime() is a misnomer before PARI version 2.3
In older versions of PARI, the one-argument variant of the
function isprime() is actually checking for probable primes.
Moreover, it has certain problems.
POSSIBLE WORKAROUND (not needed for newer PARI): before
version 2.3 of PARI, to get probability of misdetecting a prime below
1e-12, call isprime() twice; below 1e-18, call it 3 times; etc.
(The algorithm is probabilistic, and the implementation is such that the
result depends on which calls to isprime() were performed
ealier.)
The problems: first, while the default algorithm (before
version 2.3) gives practically acceptable results in non-adversarial
situations, the worst-case behaviour is significantly worse than the
average behaviour. The algorithm is looking for so-called
"witnesses" (with up to 10 tries) among random integers;
usually, witnesses are abundant. However, there are non-prime numbers
for which the fraction of witnesses is close to the theoretical mininum,
0.75; with 10 random tries, the probability of missing a witness for
such numbers is close to 1e-6. (The known worst-case numbers M have
phi(M)/4 non-witnesses, with M=P(2P-1), prime P, 2P-1 and 4|P+1; the
proportion of such numbers near K is expected to be
const/sqrt(K)log(K)^2. Note that numbers which have more than about 5%
non-witnesses may also be candidates for false positives. Conjecturally,
they are of the form (aD+1)(bD+1) with a<b, ab <= const, prime
aD+1, and bD+1, and D not divisible by high power of 2 (above a=1, b=2
and D is odd); the proportion of such numbers may have a similar
asymptotic const/sqrt(K)log(K)^2.)
Second, the random number generator is "reset to known
state" when PARI library is initialized. That means that the
behaviour is actually predictable if one knows which calls to
isprime() are performed; an adversary can find non-primes M which
will trigger a false positive exactly on the Nth call to isprime(M) (for
particular values of N). With enough computing resources, one can find
non-primes M for which N is relatively small (with M about 1e9, one can
achieve N as low as 1000). Compare with similar (but less abundant)
examples for simpler algorithm, Carmichael numbers
<http://en.wikipedia.org/wiki/Carmichael_numbers>; see also
numbers with big proportion of non-witnesses
<http://oeis.org/A090659> and numbers with many non-witnesses
<http://oeis.org/A141768>, and the conjecture about proportion
<http://web.archive.org/web/*/http://www.ma.iup.edu/MAA/proceedings/vol1/higgins.pdf>.
See the discussion of isprime()
<https://rt.cpan.org/Public/Bug/Display.html?id=93652>.
When Math::Pari is loaded, it examines variables
$Math::Pari::initmem and
$Math::Pari::initprimes. They specify up to which
number the initial list of primes should be precalculated, and how large
should be the arena for PARI calculations (in bytes). (These values have safe
defaults.)
Since setting these values before loading requires either a
"BEGIN" block, or postponing the loading
("use" vs.
"require"), it may be more convenient to
set them via Math::PariInit:
use Math::PariInit qw( primes=12000000 stack=1e8 );
"use Math::PariInit" also
accepts arbitrary Math::Pari import directives, see Math::PariInit.
These values may be changed at runtime too, via
allocatemem() and setprimelimit(), with performance penalties
for recalculation/reallocation.
Ilya Zakharevich, ilyaz@cpan.org
Hey! The above document had some coding errors, which are explained
below:
- Around line 849:
- Expected '=item *'
Visit the GSP FreeBSD Man Page Interface. Output converted with ManDoc. |