NAME

libPARI - Technical Reference Guide for Low-Level Functions


DESCRIPTION

In this chapter, we describe all public low-level functions of the PARI library. These essentially include functions for handling all the PARI types. Higher level functions, such as arithmetic or transcendental functions, are described in Chapter 3 of the GP user's manual. A general introduction to the major concepts of PARI programming can be found in Chapter 4.

Many other undocumented functions can be found throughout the source code. These private functions are more efficient than the library wrappers, but sloppier on argument checking and damage control. Use them at your own risk!

Important advice: generic routines eventually call lower level functions. Optimize your algorithms first, not overhead and conversion costs between PARI routines. For generic operations, use generic routines first, don't waste time looking for the most specialized one available unless you identify a genuine bottleneck. The PARI source code is part of the documentation; look for inspiration there.

We let BIL abbreviate BITS_IN_LONG. The type long denotes a BIL-bit signed long integer. The type ulong is defined as unsigned long. The word stack always refer to the PARI stack, allocated through an initial pari_init call. Refer to Chapters 1--2 and 4 for general background.


Initializing the library

The following functions enable you to start using the PARI functions in a program, and cleanup without exiting the whole program.

General purpose

void pari_init(size_t size, ulong maxprime) initialize the library, with a stack of size bytes and a prime table up to the maximum of maxprime and 2^{16}. Unless otherwise mentionned, no PARI function will function properly before such an initialization.

void pari_close(void) stop using the library (assuming it was initialized with pari_init) and frees all allocated objects.

Technical functions

void pari_init_opts(size_t size, ulong maxprime, ulong opts) as pari_init, more flexible. opts is a mask of flags among the following:

INIT_JMPm: install pari error handler. When an exception is raised, the program is terminated with exit(1).

INIT_SIGm: install pari signal handler.

INIT_DFTm: initialize the GP_DATA environment structure. This one must be enabled once. If you close pari, then restart it, you need not reinitialize GP_DATA; if you do not, then old values are restored.

void pari_close_opts(ulong init_opts) as pari_close, for a library initialized with a mask of options using pari_init_opts. opts is a mask of flags among

INIT_SIGm: restore SIG_DFL default action for signals tampered with by pari signal handler.

INIT_DFTm: frees the GP_DATA environment structure.

void pari_sig_init(void (*f)(int)) install the signal handler f (see signal(2)): the signals SIGBUS, SIGFPE, SIGINT, SIGBREAK, SIGPIPE and SIGSEGV are concerned.

Notions specific to the GP interpreter

An entree is the generic object associated to an identifier (a name) in GP's interpreter, be it a built-in or user function, or a variable. For a function, it has at least the following fields:

char *name : the name under which the interpreter knows us.

ulong valence : obsolete, set it to 1.

void *value : a pointer to the C function to call.

long menu : an integer from 1 to 11 (to which group of function help do we belong).

char *code : the prototype code.

char *help : the help text for the function.

A routine in GP is described to the analyzer by an entree structure. Built-in pari routines are grouped in modules, which are arrays of entree structs, the last of which satisfy name = NULL (sentinel).

There are currently six modules in GP: general functions (functions_basic), gp-specific functions (functions_fp), gp-specific highlevel functions (functions_highlevel), member functions, and two modules of obsolete functions. The function pari_init initializes the interpreter and declares all symbols in functions_basic. You may declare further functions on a case by case basis or as a whole module using

void pari_add_function(entree *ep) adds a single routine to the table of symbols in the interpreter. It assumes pari_init has been called.

void pari_add_module(entree *mod) adds all the routines in module mod to the table of symbols in the interpreter. It assumes pari_init has been called.

For instance, gp implements a number of private routines, which it adds to the default set via the call

    pari_add_module(functions_gp);
    pari_add_module(functions_highlevel);


Handling GENs

Almost all these functions are either macros or inlined. Unless mentioned otherwise, they do not evaluate their arguments twice. Most of them are specific to a set of types, although no consistency checks are made: e.g. one may access the sign of a t_PADIC, but the result is meaningless.

Length conversions

long ndec2nlong(long x) converts a number of decimal digits to a number of words. Returns 1 + floor(x x BIL log _2 10).

long ndec2prec(long x) converts a number of decimal digits to a number of codewords. This is equal to 2 + ndec2nlong(x).

long prec2ndec(long x) converts a number of of codewords to a number of decimal digits.

long nbits2nlong(long x) converts a number of bits to a number of words. Returns the smallest word count containing x bits, i.e ceil(x / BIL).

long nbits2prec(long x) converts a number of bits to a number of codewords. This is equal to 2 + nbits2nlong(x).

long nchar2nlong(long x) converts a number of bytes to number of words. Returns the smallest word count containing x bytes, i.e ceil(x / sizeof(long)).

long bit_accuracy(long x) converts a t_REAL length into a number of significant bits. Returns (x - 2)BIL. The macro bit_accuracy_mul(x,y) computes the same thing multiplied by y.

Read type-dependent information

long typ(GEN x) returns the type number of x. The header files included through pari.h define symbolic constants for the GEN types: t_INT etc. Never use their actual numerical values. E.g to determine whether x is a t_INT, simply check

    if (typ(x) == t_INT) { }

The types are internally ordered and this simplifies the implementation of commutative binary operations (e.g addition, gcd). Avoid using the ordering directly, as it may change in the future; use type grouping macros instead (Label se:typegroup).

long lg(GEN x) returns the length of x in BIL-bit words.

long lgefint(GEN x) returns the effective length of the t_INT x in BIL-bit words.

long signe(GEN x) returns the sign (-1, 0 or 1) of x. Can be used for t_INT, t_REAL, t_POL and t_SER (for the last two types, only 0 or 1 are possible).

long gsigne(GEN x) same as signe, but also valid for t_FRAC (and marginally less efficient for the other types). Raise a type error if typ(x) is not among those three.

long expi(GEN x) returns the binary exponent of the real number equal to the t_INT x. This is a special case of gexpo.

long expo(GEN x) returns the binary exponent of the t_REAL x.

long gexpo(GEN x) same as expo, but also valid when x is not a t_REAL (returns the largest exponent found among the components of x). When x is an exact 0, this returns -HIGHEXPOBIT, which is lower than any valid exponent.

long valp(GEN x) returns the p-adic valuation (for a t_PADIC) or X-adic valuation (for a t_SER, taken with respect to the main variable) of x.

long precp(GEN x) returns the precision of the t_PADIC x.

long varn(GEN x) returns the variable number of the t_POL or t_SER x (between 0 and MAXVARN).

long gvar(GEN x) returns the main variable number when any variable at all occurs in the composite object x (the smallest variable number which occurs), and BIGINT otherwise.

long degpol(GEN x) returns the degree of t_POL x, assuming its leading coefficient is non-zero (an exact 0 is impossible, but an inexact 0 is allowed). By convention the degree of an exact 0 polynomial is -1. If the leading coefficient of x is 0, the result is undefined.

int precision(GEN x) If x is of type t_REAL, returns the precision of x (the length of x in BIL-bit words if x is not zero, and a reasonable quantity obtained from the exponent of x if x is numerically equal to zero). If x is of type t_COMPLEX, returns the minimum of the precisions of the real and imaginary part. Otherwise, returns 0 (which stands in fact for infinite precision).

int gprecision(GEN x) as precision for scalars; returns the lowest precision encountered among the components otherwise.

long sizedigit(GEN x) returns 0 if x is exactly 0. Otherwise, returns gexpo(x) multiplied by log _{10}(2). This gives a crude estimate for the maximal number of decimal digits of the components of x.

Eval type-dependent information.

These routines convert type-dependant information to bitmask to fill the codewords of GEN objects (see Label se:impl). E.g for a t_REAL z:

    z[1] = evalsigne(-1) | evalexpo(2)

Compatible components of a codeword for a given type can be OR-ed as above.

ulong evaltyp(long x) convert type x to bitmask (first codeword of all GENs)

long evallg(long x) convert length x to bitmask (first codeword of all GENs). Raise overflow error if x is so large that the corresponding length cannot be represented

long _evallg(long x) as evallg without the overflow check.

ulong evalvarn(long x) convert variable number x to bitmask (second codeword of t_POL and t_SER)

long evalsigne(long x) convert sign x (in -1,0,1) to bitmask (second codeword of t_INT, t_REAL, t_POL, t_SER)

long evalprecp(long x) convert p-adic (X-adic) precision x to bitmask (second codeword of t_PADIC, t_SER)

long evalvalp(long x) convert p-adic (X-adic) valuation x to bitmask (second codeword of t_PADIC, t_SER). Raise overflow error if x is so large that the corresponding valuation cannot be represented

long _evalvalp(long x) same as evalvalp without the overflow check.

long evalexpo(long x) convert exponent x to bitmask (second codeword of t_REAL). Raise overflow error if x is so large that the corresponding exponent cannot be represented

long _evalexpo(long x) same as evalexpo without the overflow check.

long evallgefint(long x) convert effective length x to bitmask (second codeword t_INT). This should be less or equal than the length of the t_INT, hence there is no overflow check for the effective length.

long evallgeflist(long x) convert effective length x to bitmask (second codeword t_LIST). This should be less or equal than the length of the t_LIST, hence there is no overflow check for the effective length.

Set type-dependent information.

Use these macros with extreme care since usually the corresponding information is set otherwise, and the components and further codeword fields (which are left unchanged) may not be compatible with the new information.

void settyp(GEN x, long s) sets the type number of x to s.

void setlg(GEN x, long s) sets the length of x to s. This is an efficient way of truncating vectors, matrices or polynomials.

void setlgefint(GEN x, long s) sets the effective length of the t_INT x to s. The number s must be less than or equal to the length of x.

void setsigne(GEN x, long s) sets the sign of x to s. If x is a t_INT or t_REAL, s must be equal to -1, 0 or 1, and if x is a t_POL or t_SER, s must be equal to 0 or 1.

void setexpo(GEN x, long s) sets the binary exponent of the t_REAL x to s. The value s must be a 24-bit signed number.

void setvalp(GEN x, long s) sets the p-adic or X-adic valuation of x to s, if x is a t_PADIC or a t_SER, respectively.

void setprecp(GEN x, long s) sets the p-adic precision of the t_PADIC x to s.

void setvarn(GEN x, long s) sets the variable number of the t_POL or t_SER x to s (where 0 <= s <= MAXVARN).

Type groups

. In the following macros, t denotes the type of a GEN. Some of these macros may evaluate their argument twice. Always use them as in

    long tx = typ(x);
    if (is_intreal_t(tx)) { }

int is_recursive_t(long t) true iff t is a recursive type (the recursive types are t_INT, t_REAL, t_STR or t_VECSMALL).

int is_intreal_t(long t) true iff t is t_INT or t_REAL.

int is_rational_t(long t) true iff t is t_INT or t_FRAC.

int is_vec_t(long t) true iff t is t_VEC or t_COL.

int is_matvec_t(long t) true iff t is t_MAT, t_VEC or t_COL.

int is_scalar_t(long t) true iff t is a scalar, i.e a t_INT, t_REAL, t_INTMOD, t_FRAC, t_COMPLEX, t_PADIC, t_QUAD, or t_POLMOD.

int is_extscalar_t(long t) true iff t is a scalar (see is_scalar_t) or t is t_POL.

int is_const_t(long t) true iff t is a scalar which is not t_POLMOD.

Accessors and components

The first two functions return GEN components as copies on the stack:

GEN compo(GEN x, long n) creates a copy of the n-th true component (i.e. not counting the codewords) of the object x.

GEN truecoeff(GEN x, long n) creates a copy of the coefficient of degree n of x if x is a scalar, t_POL or t_SER, and otherwise of the n-th component of x.

On the contrary, the following routines return the address of a GEN component. No copy is made on the stack:

GEN constant_term(GEN x) returns the address the constant term of t_POL x. By convention, a 0 polynomial (whose sign is 0) has gen_0 constant term.

GEN leading_term(GEN x) returns the address the leading term of t_POL x. This may be an inexact 0.

GEN gel(GEN x, long i) returns the address of the x[i] entry of x. (el stands for element.)

GEN gcoeff(GEN x, long i, long j) returns the address of the x[i,j] entry of t_MAT x, i.e. the coefficient at row i and column j.

GEN gmael(GEN x, long i, long j) returns the address of the x[i][j] entry of x. (mael stands for multidimensional array element.)

GEN gmael2(GEN A, long x1, long x2) is an alias for gmael. Similar macros gmael3, gmael4, gmael5 are available.


Handling the PARI stack

Allocating memory on the stack

GEN cgetg(long n, long t) allocates memory on the stack for an object of length n and type t, and initializes its first codeword.

GEN cgeti(long n) allocates memory on the stack for a t_INT of length n, and initializes its first codeword. Identical to cgetg(n,t_INT).

GEN cgetr(long n) allocates memory on the stack for a t_REAL of length n, and initializes its first codeword. Identical to cgetg(n,t_REAL).

GEN cgetc(long n) allocates memory on the stack for a t_COMPLEX, whose real and imaginary parts are t_REALs of length n.

GEN cgetp(GEN x) creates space sufficient to hold the t_PADIC x, and sets the prime p and the p-adic precision to those of x, but does not copy (the p-adic unit or zero representative and the modulus of) x.

GEN new_chunk(size_t n) allocates a GEN with n components, without filling the required code words. This is the low-level constructor underlying cgetg, which calls new_chunk then sets the first code word. It works by simply returning the address ((GEN)avma) - n, after checking that it is larger than (GEN)bot.

char* stackmalloc(size_t n) allocates memory on the stack for n chars (not n GENs). This is faster than using malloc, and easier to use in most situations when temporary storage is needed. In particular there is no need to free individually all variables thus allocated: a simple avma = oldavma might be enough. On the other hand, beware that this is not permanent independant storage, but part of the stack.

Objects allocated through these last two functions cannot be gerepile'd. They are not valid GENs since they have no PARI type.

Garbage collection.

See Label se:garbage for a detailed explanation and many examples.

void cgiv(GEN x) frees object x if it is the last created on the stack (otherwise nothing happens).

GEN gerepile(pari_sp p, pari_sp q, GEN x) general garbage collector for the stack.

void gerepileall(pari_sp av, int n, ...) cleans up the stack from av on (i.e from avma to av), preserving the n objects which follow in the argument list (of type GEN*). E.g: gerepileall(av, 2, &x, &y) preserves x and y.

void gerepileallsp(pari_sp av, pari_sp ltop, int n, ...) cleans up the stack between av and ltop, updating the n elements which follow n in the argument list (of type GEN*). Check that the elements of g have no component between av and ltop, and assumes that no garbage is present between avma and ltop. Analogous to (but faster than) gerepileall otherwise.

GEN gerepilecopy(pari_sp av, GEN x) cleans up the stack from av on, preserving the object x. Special case of gerepileall (case n = 1), except that the routine returns the preserved GEN instead of updating its adress through a pointer.

void gerepilemany(pari_sp av, GEN* g[], int n) alternative interface to gerepileall

void gerepilemanysp(pari_sp av, pari_sp ltop, GEN* g[], int n) alternative interface to gerepileallsp.

void gerepilecoeffs(pari_sp av, GEN x, int n) cleans up the stack from av on, preserving x[0],..., x[n-1] (which are GENs).

void gerepilecoeffssp(pari_sp av, pari_sp ltop, GEN x, int n) cleans up the stack from av to ltop, preserving x[0], ..., x[n-1] (which are GENs). Same assumptions as in gerepilemanysp, of which this is a variant. For instance

    z = cgetg(3, t_COMPLEX);
    av = avma; garbage(); ltop = avma;
    z[1] = fun1();
    z[2] = fun2();
    gerepilecoeffssp(av, ltop, z + 1, 2);
    return z;

cleans up the garbage between av and ltop, and connects z and its two components. This is marginally more efficient than the standard

    av = avma; garbage(); ltop = avma;
    z = cgetg(3, t_COMPLEX);
    z[1] = fun1();
    z[2] = fun2(); return gerepile(av, ltop, z);

GEN gerepileupto(pari_sp av, GEN q) analogous to (but faster than) gerepilecopy. Assumes that q is connected and that its root was created before any component.

GEN gerepileuptoint(pari_sp av, GEN q) analogous to (but faster than) gerepileupto. Assumes further that q is a t_INT. The length and effective length of the resulting t_INT are equal.

GEN gerepileuptoleaf(pari_sp av, GEN q) analogous to (but faster than) gerepileupto. Assumes further that q is a leaf, i.e a non-recursive type (is_recursive_t(typ(q)) is non-zero). Contrary to gerepileuptoint, gerepileuptoleaf leaves length and effective length of a t_INT unchanged.

void stackdummy(pari_sp av, pari_sp ltop) inhibits the memory area between av included and ltop excluded with respect to gerepile, in order to avoid a call to gerepile(av, ltop,...). The stack space is not reclaimed though.

More precisely, this routine assumes that av is recorded earlier than ltop, then marks the specified stack segment as a non-recursive type of the correct length. Thus gerepile will not inspect the zone, at most copy it. To be used in the following situation:

    av0 = avma; z = cgetg(t_VEC, 3);
    gel(z,1) = HUGE(); av = avma; garbage(); ltop = avma;
    gel(z,2) = HUGE(); stackdummy(av, ltop);

Compared to the orthodox

    gel(z,2) = gerepile(av, ltop, gel(z,2));

or even more wasteful

    z = gerepilecopy(av0, z);

we temporarily lose (av - ltop) words but save a costly gerepile. In principle, a garbage collection higher up the call chain should reclaim this later anyway.

Without the stackdummy, if the [av, ltop] zone is arbitrary (not even valid GENs as could happen after direct truncation via setlg), we would leave dangerous data in the middle of z, which would be a problem for a later

    gerepile(..., ... , z);

And even if it were made of valid GENs, inhibiting the area makes sure gerepile will not inspect their components, saving time.

Another natural use in low-level routines is to ``shorten'' an existing GEN z to its first l-1 components:

    setlg(z, l);
    stackdummy((pari_sp)(z + lg(z)), (pari_sp)(z + l));

or to its last l components:

    long L = lg(z) - l;
    stackdummy((pari_sp)(z + L), (pari_sp)z);
    z += L; setlg(z, L);

Copies and clones

GEN gclone(GEN x) creates a new permanent copy of the object x on the heap. The clone bit of the result is set.

void gunclone(GEN x) delete the clone x (created by gclone). Fatal error if x not a clone.

GEN gcopy(GEN x) creates a new copy of the object x on the stack.

int isonstack(GEN x) true iff x belongs to the stack. This is a macro whose argument is evaluated several times.

void copyifstack(GEN x, GEN y) sets y = gcopy(x) if x belongs to the stack, and y = x otherwise. This macro evaluates its arguments once, contrary to

    y = isonstack(x)? gcopy(x): x;

void icopyifstack(GEN x, GEN y) as copyifstack assuming x is a t_INT.

long taille(GEN x) returns the total number of BIL-bit words occupied by the tree representing x.

void traverseheap(void(*f)(GEN, void *), void *data) this applies f(x, data) to each object x on the PARI heap, most recent first. Mostly for debugging purposes.

GEN getheap() a simple wrapper around traverseheap. Returns a two-component row vector giving the number of objects on the heap and the amount of memory they occupy in long words.


Level 0 kernel (operations on ulongs)

Micro-kernel.

Level 0 operations simulate basic operations of the 68020 processor on which PARI was originally implemented. They need ``global'' ulong variables overflow (which will contain only 0 or 1) and hiremainder to function properly. However, for certain architectures these are replaced with local variables for efficiency; and the `functions' mentioned below are really chunks of inlined assembler code. So, a routine using one of these lowest-level functions where the description mentions either hiremainder or overflow must declare the corresponding

    LOCAL_HIREMAINDER;
    LOCAL_OVERFLOW;

in a declaration block. Variables hiremainder and overflow then become available in the enclosing block. For instance a loop over the powers of an ulong p protected from overflows could read

   while (pk < lim)
   {
     LOCAL_HIREMAINDER;
     ...
     pk = mulll(pk, p); if (hiremainder) break;
   }

ulong addll(ulong x, ulong y) adds x and y, returns the lower BIL bits and puts the carry bit into overflow.

ulong addllx(ulong x, ulong y) adds overflow to the sum of the x and y, returns the lower BIL bits and puts the carry bit into overflow.

ulong subll(ulong x, ulong y) subtracts x and y, returns the lower BIL bits and put the carry (borrow) bit into overflow.

ulong subllx(ulong x, ulong y) subtracts overflow from the difference of x and y, returns the lower BIL bits and puts the carry (borrow) bit into overflow.

int bfffo(ulong x) returns the number of leading zero bits in x. That is, the number of bit positions by which it would have to be shifted left until its leftmost bit first becomes equal to 1, which can be between 0 and BIL-1 for nonzero x. When x is 0, the result is undefined.

ulong mulll(ulong x, ulong y) multiplies x by y, returns the lower BIL bits and stores the high-order BIL bits into hiremainder.

ulong addmul(ulong x, ulong y) adds hiremainder to the product of x and y, returns the lower BIL bits and stores the high-order BIL bits into hiremainder.

ulong divll(ulong x, ulong y) returns the Euclidean quotient of (hiremainder << BIL)+x by y and stores the remainder into hiremainder. An error occurs if the quotient cannot be represented by an ulong, i.e. if initially hiremainder >= y.

Modular kernel.

The following routines are not part of the level 0 kernel per se, but implement modular operations on words in terms of the above. They are written so that no overflow may occur. Let m >= 1 be the modulus; all operands representing classes modulo m are assumed to belong to [0,m-1[. The result may be wrong for a number of reasons otherwise: it may not be reduced, overflow can occur, etc.

ulong Fl_add(ulong x, ulong y, ulong m) returns the smallest positive representative of x + y modulo m.

ulong Fl_neg(ulong x, ulong m) returns the smallest positive representative of -x modulo m.

ulong Fl_sub(ulong x, ulong y, ulong m) returns the smallest positive representative of x - y modulo m.

long Fl_center(ulong x, ulong m, ulong mo2) returns the representative in ]-m/2,m/2] of x modulo m. Assume 0 <= x < m and mo2 = m >> 1.

ulong Fl_mul(ulong x, ulong y, ulong m) returns the smallest positive representative of x y modulo m.

ulong Fl_inv(ulong x, ulong m) returns the smallest positive representative of x^{-1} modulo m. If x is not invertible mod m, raise an exception.

ulong Fl_div(ulong x, ulong y, ulong m) returns the smallest positive representative of x y^{-1} modulo m. If y is not invertible mod m, raise an exception.

ulong Fl_pow(ulong x, ulong n, ulong m) returns the smallest positive representative of x^n modulo m.

ulong Fl_sqrt(ulong x, ulong p) returns the square root of x modulo p (smallest positive representative). Assumes p to be prime, and x to be a square modulo p.

ulong gener_Fl(ulong p) returns a primitive root modulo p, assuming p is prime.

ulong gener_Fl_local(ulong p, GEN L), see gener_Fp_local, L is an Flv.

long krouu(ulong x, ulong y) returns the Kronecker symbol (x|y), i.e.-1, 0 or 1. Assumes y is non-zero. If y is an odd prime, this is the Legendre symbol.


Level 1 kernel (operations on longs, integers and reals)

Note: Many functions consist of an elementary operation, immediately followed by an assignment statement. They will be introduced as in the following example:

GEN gadd[z](GEN x, GEN y[, GEN z]) followed by the explicit description of the function

GEN gadd(GEN x, GEN y)

which creates its result on the stack, returning a GEN pointer to it, and the parts in brackets indicate that there exists also a function

void gaddz(GEN x, GEN y, GEN z)

which assigns its result to the pre-existing object z, leaving the stack unchanged. All such functions are obtained using macros (see the file paricom.h), hence you can easily extend the list. These assignment variants are inefficient; don't use them.

Creation

GEN cgeti(long n) allocates memory on the PARI stack for a t_INT of length n, and initializes its first codeword. Identical to cgetg(n,t_INT).

GEN cgetr(long n) allocates memory on the PARI stack for a t_REAL of length n, and initializes its first codeword. Identical to cgetg(n,t_REAL).

GEN cgetc(long n) allocates memory on the PARI stack for a t_COMPLEX, whose real and imaginary parts are t_REALs of length n.

GEN real_1(long prec) create a t_REAL equal to 1 to prec words of accuracy.

GEN real_m1(long prec) create a t_REAL equal to -1 to prec words of accuracy.

GEN real_0_bit(long bit) create a t_REAL equal to 0 with exponent -bit.

GEN real_0(long prec) is a shorthand for

    real_0_bit( -bit_accuracy(prec) )

GEN int2n(long n) creates a t_INT equal to 1 << n (i.e 2^n if n >= 0, and 0 otherwise).

GEN int2u(ulong n) creates a t_INT equal to 2^n.

GEN real2n(long n, long prec) create a t_REAL equal to 2^n to prec words of accuracy.

GEN stroi(char *s) convert the character string s to a t_INT.

GEN stror(char *s, long prec) convert the character string s to a t_REAL of precision prec.

Assignment.

In this section, the z argument in the z-functions must be of type t_INT or t_REAL.

void mpaff(GEN x, GEN z) assigns x into z (where x and z are t_INT or t_REAL). Assumes that lg(z) > 2.

void affii(GEN x, GEN z) assigns the t_INT x into the t_INT z.

void affir(GEN x, GEN z) assigns the t_INT x into the t_REAL z. Assumes that lg(z) > 2.

void affiz(GEN x, GEN z) assigns t_INT x into t_INT or t_REAL z. Assumes that lg(z) > 2.

void affsi(long s, GEN z) assigns the long s into the t_INT z. Assumes that lg(z) > 2.

void affsr(long s, GEN z) assigns the long s into the t_REAL z. Assumes that lg(z) > 2.

void affsz(long s, GEN z) assigns the long s into the t_INT or t_REAL z. Assumes that lg(z) > 2.

void affui(ulong u, GEN z) assigns the ulong u into the t_INT z. Assumes that lg(z) > 2.

void affur(ulong u, GEN z) assigns the ulong u into the t_REAL z. Assumes that lg(z) > 2.

void affrr(GEN x, GEN z) assigns the t_REAL x into the t_REAL z.

The function affrs and affri do not exist. So don't use them.

Copy

GEN icopy(GEN x) copy relevant words of the t_INT x on the stack: the length and effective length of the copy are equal.

GEN rcopy(GEN x) copy the t_REAL x on the stack.

GEN mpcopy(GEN x) copy the t_INT or t_REAL x on the stack. Contrary to icopy, mpcopy preserves the original length of a t_INT.

Conversions

GEN itor(GEN x, long prec) converts the t_INT x to a t_REAL of length prec and return the latter. Assumes that prec > 2.

long itos(GEN x) converts the t_INT x to a long if possible, otherwise raise an exception.

long itos_or_0(GEN x) converts the t_INT x to a long if possible, otherwise return 0.

ulong itou(GEN x) converts the t_INT |x| to an ulong if possible, otherwise raise an exception.

long itou_or_0(GEN x) converts the t_INT |x| to an ulong if possible, otherwise return 0.

GEN stoi(long s) creates the t_INT corresponding to the long s.

GEN stor(long s, long prec) converts the long s into a t_REAL of length prec and return the latter. Assumes that prec > 2.

GEN utoi(ulong s) converts the ulong s into a t_INT and return the latter.

GEN utoipos(ulong s) converts the non-zero ulong s into a t_INT and return the latter.

GEN utoineg(ulong s) converts the non-zero ulong s into a t_INT and return the latter.

GEN utor(ulong s, long prec) converts the ulong s into a t_REAL of length prec and return the latter. Assumes that prec > 2.

GEN rtor(GEN x, long prec) converts the t_REAL x to a t_REAL of length prec and return the latter. If prec < lg(x), round properly. If prec > lg(x), padd with zeroes. Assumes that prec > 2.

The following function is also available as a special case of mkintn:

GEN u2toi(ulong a, ulong b)

Returns the GEN equal to 2^{32} a + b, assuming that a,b < 2^{32}. This does not depend on sizeof(long): the behaviour is as above on both 32 and 64-bit machines.

Integer parts

GEN ceilr(GEN x) smallest integer larger or equal to the t_REAL x (i.e. the ceil function).

GEN floorr(GEN x) largest integer smaller or equal to the t_REAL x (i.e. the floor function).

GEN roundr(GEN x) rounds the t_REAL x to the nearest integer (towards + oo ).

GEN truncr(GEN x) truncates the t_REAL x (not the same as floorr if x is and negative).

GEN mpceil[z](GEN x[, GEN z]) as ceilr except that x may be a t_INT.

GEN ceil_safe(GEN x), x being a real number (not necessarily a t_REAL) returns an integer which is larger than any possible incarnation of x. (Recall that a t_REAL represents an interval of possible values.)

GEN mpfloor[z](GEN x[, GEN z]) as floorr except that x may be a t_INT.

GEN mpround[z](GEN x[, GEN z]) as roundr except that x may be a t_INT.

GEN mptrunc[z](GEN x[, GEN z]) as truncr except that x may be a t_INT.

GEN diviiround(GEN x, GEN y) if x and y are t_INTs, returns the quotient x/y of x and y, rounded to the nearest integer. If x/y falls exactly halfway between two consecutive integers, then it is rounded towards + oo (as for roundr).

Valuation and shift

long vals(long s) 2-adic valuation of the long s. Returns -1 if s is equal to 0.

long vali(GEN x) 2-adic valuation of the t_INT x. Returns -1 if x is equal to 0.

GEN mpshift[z](GEN x, long n[, GEN z]) shifts the t_INT or t_REAL x by n. If n is positive, this is a left shift, i.e. multiplication by 2^{n}. If n is negative, it is a right shift by -n, which amounts to the truncation of the quotient of x by 2^{-n}.

GEN shifti(GEN x, long n) shifts the t_INT x by n.

GEN shiftr(GEN x, long n) shifts the t_REAL x by n.

long Z_pvalrem(GEN x, GEN p, GEN *r) applied to t_INTs x != 0 and p, |p| > 1, returns the highest exponent e such that p^{e} divides x. The quotient x/p^{e} is returned in *r. In particular, if p is a prime, this returns the valuation at p of x, and *r is the prime-to-p part of x.

long Z_pval(GEN x, GEN p) as Z_pvalrem but only returns the ``valuation''.

long Z_lvalrem(GEN x, ulong p, GEN *r) as Z_pvalrem, except that p is an ulong (p > 1).

long Z_lval(GEN x, ulong p) as Z_pval, except that p is an ulong (p > 1).

long u_lvalrem(ulong x, ulong p, ulong *r) as Z_pvalrem, except the inputs/outputs are now ulongs.

long u_pvalrem(ulong x, GEN p, ulong *r) as Z_pvalrem, except x and r are now ulongs.

long u_lval(ulong x, ulong p) as Z_pval, except the inputs/outputs are now ulongs.

Factorization

GEN Z_factor(GEN n) factors the t_INT n. The ``primes'' in the factorization are actually strong pseudoprimes.

long Z_issquarefree(GEN x) returns 1 if the t_INT n is square-free, and 0 otherwise.

long Z_issquare(GEN n) returns 1 if t_INT n is a square, and 0 otherwise. This is tested first modulo small prime powers, then sqrtremi is called.

long Z_issquarerem(GEN n, GEN *sqrtn) as Z_issquare. If n is indeed a square, set sqrtn to its integer square root.

int isprime(GEN n), returns 1 if the t_INT n is a (fully proven) prime number and 0 otherwise.

int uisprime(ulong p), returns 1 if p is a prime number and 0 otherwise.

long Z_issquarerem(GEN n, GEN *sqrtn) as Z_issquare. If n is indeed a square, set sqrtn to its integer square root.

long uissquarerem(ulong n, ulong *sqrtn) as Z_issquarerem, for an ulong operand n.

Generic unary operators.

Let ``op'' be a unary operation among

op = neg: negation (-x).

op = abs: absolute value (|x|).

The names and prototypes of the low-level functions corresponding to op are as follows. The result is of the same type as x.

GEN mpop(GEN x) creates the result of op applied to the t_INT or t_REAL x.

GEN opi(GEN x) creates the result of op applied to the t_INT x.

GEN opr(GEN x) creates the result of op applied to the t_REAL x.

GEN mpopz(GEN x, GEN z) assigns the result of applying op to the t_INT or t_REAL x into the t_INT or t_REAL z.

Remark: it has not been considered useful to include functions void opsz(long,GEN), void opiz(GEN,GEN) and void oprz(GEN, GEN).

Comparison operators

int mpcmp(GEN x, GEN y) compares the t_INT or t_REAL x to the t_INT or t_REAL y. The result is the sign of x-y.

int cmpii(GEN x, GEN y) compares the t_INT x to the t_INT y.

int cmpir(GEN x, GEN y) compares the t_INT x to the t_REAL y.

int cmpis(GEN x, long s) compares the t_INT x to the long s.

int cmpsi(long s, GEN x) compares the long s to the t_INT x.

int cmpsr(long s, GEN x) compares the long s to the t_REAL x.

int cmpri(GEN x, GEN y) compares the t_REAL x to the t_INT y.

int cmprr(GEN x, GEN y) compares the t_REAL x to the t_REAL y.

int cmprs(GEN x, long s) compares the t_REAL x to the long s.

int equalii(GEN x, GEN y) compares the t_INTs x and y. The result is 1 if x = y, 0 otherwise.

int equalsi(long s, GEN x)

int equalis(GEN x, long s) compare the t_INT x and the long s. The result is 1 if x = y, 0 otherwise.

int equalui(ulong s, GEN x)

int equaliu(GEN x, ulong s) compare the t_INT x and the ulong s. The result is 1 if |x |= y, 0 otherwise.

int absi_cmp(GEN x, GEN y) compares the t_INTs x and y. The result is the sign of |x| - |y|.

int absi_equal(GEN x, GEN y) compares the t_INTs x and y. The result is 1 if |x |= |y|, 0 otherwise.

int absr_cmp(GEN x, GEN y) compares the t_REALs x and y. The result is the sign of |x| - |y|.

Generic binary operators.

Let ``op'' be a binary operation among

op = add: addition (x + y). The result is a t_REAL unless both x and y are t_INTs (or longs).

op = sub: subtraction (x - y). The result is a t_REAL unless both x and y are t_INT (or longs).

op = mul: multiplication (x * y). The result is a t_REAL unless both x and y are t_INTs (or longs), or if x or y is an exact 0.

op = div: division (x / y). In the case where x and y are both t_INTs or longs, the result is the Euclidean quotient, where the remainder has the same sign as the dividend x. It is the ordinary division otherwise. If one of x or y is a t_REAL, the result is a t_REAL unless x is an exact 0. A division-by-0 error occurs if y is equal to 0.

op = rem: remainder (``x % y''). This operation is defined only when x and y are longs or t_INT. The result is the Euclidean remainder corresponding to div, i.e. its sign is that of the dividend x. The result is always a t_INT.

op = mod: true remainder (x % y). This operation is defined only when x and y are longs or t_INTs. The result is the true Euclidean remainder, i.e. non-negative and less than the absolute value of y.

The names and prototypes of the low-level functions corresponding to op are as follows. In this section, the z argument in the z-functions must be of type t_INT or t_REAL. t_INT is only allowed when no `r' appears in the argument code (no t_REAL operand is involved).

GEN mpop[z](GEN x, GEN y[, GEN z]) applies op to the t_INT or t_REAL x and y.

GEN opsi[z](long s, GEN x[, GEN z]) applies op to the long s and the t_INT x.

GEN opsr[z](long s, GEN x[, GEN z]) applies op to the long s and the t_REAL x.

GEN opss[z](long s, long t[, GEN z]) applies op to the longs s and t.

GEN opii[z](GEN x, GEN y[, GEN z]) applies op to the t_INTs x and y.

GEN opir[z](GEN x, GEN y[, GEN z]) applies op to the t_INT x and the t_REAL y.

GEN opis[z](GEN x, long s[, GEN z]) applies op to the t_INT x and the long s.

GEN opri[z](GEN x, GEN y[, GEN z]) applies op to the t_REAL x and the t_INT y.

GEN oprr[z](GEN x, GEN y[, GEN z]) applies op to the t_REALx and y.

GEN oprs[z](GEN x, long s[, GEN z]) applies op to the t_REAL x and the long s.

Some miscellaneous routines whose meaning should be clear from their names:

GEN muluu(ulong x, ulong y)

GEN mului(ulong x, GEN y)

GEN muliu(GEN x, ulong y)

GEN sqri(GEN x) squares the t_INT x

GEN truedivii(GEN x, GEN y) returns the true Euclidean quotient (with non-negative remainder less than |y|).

GEN truedivis(GEN x, long y) returns the true Euclidean quotient (with non-negative remainder less than |y|).

GEN centermodii(GEN x, GEN y, GEN y2), given t_INTs x, y, returns z congruent to x modulo y, such that -y/2 <= z < y/2. Assumes that y2 = shifti(y, -1). the representative of ssquares the t_INT x

Modulo to longs.

The following variants of modii do not clutter the stack:

long smodis(GEN x, long y) computes the true Euclidean remainder of the t_INT x by the long y. This is the non-negative remainder, not the one whose sign is the sign of x as in the div functions.

long smodsi(long x, GEN y) computes the true Euclidean remainder of the long x by a t_INT y.

long smodss(long x, long y) computes the true Euclidean remainder of the long x by a t_long y.

ulong umodiu(GEN x, ulong y) computes the true Euclidean remainder of the t_INT x by the ulong y.

ulong umodui(ulong x, GEN y) computes the true Euclidean remainder of the ulong x by the t_INT |y|.

The routine smodsi does not exist, since it would not always be defined: for a negative x, its result x + |y| would in general not fit into a long. Use either umodui or modsi.

Exact division and divisibility

void diviiexact(GEN x, GEN y) returns the Euclidean quotient x / y, assuming y divides x. Uses Jebelean algorithm (Jebelean-Krandick bidirectional exact division is not implemented).

void diviuexact(GEN x, ulong y) returns the Euclidean quotient |x| / y (note the absolue value!), assuming y divides x and y is non-zero.

int dvdii(GEN x, GEN y) if the t_INT y divides the t_INT x, returns 1 (true), otherwise returns 0 (false).

int dvdiiz(GEN x, GEN y, GEN z) if the t_INT y divides the t_INT x, assigns the quotient to the t_INT z and returns 1 (true), otherwise returns 0 (false).

int dvdisz(GEN x, long y, GEN z) if the t_long y divides the t_INT x, assigns the quotient to the t_INT z and returns 1 (true), otherwise returns 0 (false).

int dvdiuz(GEN x, ulong y, GEN z) if the t_ulong y divides the t_INT x, assigns the quotient |x|/y to the t_INT z and returns 1 (true), otherwise returns 0 (false).

Division with remainder.

The following functions return two objects, unless specifically asked for only one of them --- a quotient and a remainder. The quotient is returned and the remainder is returned through the variable whose address is passed as the r argument. The term true Euclidean remainder refers to the non-negative one (mod), and Euclidean remainder by itself to the one with the same sign as the dividend (rem). All GENs, whether returned directly or through a pointer, are created on the stack.

GEN dvmdii(GEN x, GEN y, GEN *r) returns the Euclidean quotient of the t_INT x by a t_INT y and puts the remainder into *r. If r is equal to NULL, the remainder is not created, and if r is equal to ONLY_REM, only the remainder is created and returned. In the generic case, the remainder is created after the quotient and can be disposed of individually with a cgiv(r). The remainder is always of the sign of the dividend x. If the remainder is 0 set r = gen_0.

void dvmdiiz(GEN x, GEN y, GEN z, GEN t) assigns the Euclidean quotient of the t_INTs x and y into the t_INT or t_REAL z, and the Euclidean remainder into the t_INT or t_REAL t.

Analogous routines dvmdis[z], dvmdsi[z], dvmdss[z] are available, where s denotes a long argument. But the following routines are in general more flexible:

long sdivss_rem(long s, long t, long *r) computes the Euclidean quotient and remainder of the longs s and t. Puts the remainder into *r, and returns the quotient. The remainder is of the sign of the dividend s, and has strictly smaller absolute value than t.

long sdivsi_rem(long s, GEN x, long *r) computes the Euclidean quotient and remainder of the long s by the t_INT x. As sdivss_rem otherwise.

long sdivsi(long s, GEN x) as sdivsi_rem, without remainder.

GEN divis_rem(GEN x, long s, long *r) computes the Euclidean quotient and remainder of the t_INT x by the long s. As sdivss_rem otherwise.

GEN diviu_rem(GEN x, ulong s, long *r) computes the Euclidean quotient and remainder of the t_INT x by the ulong s. As sdivss_rem otherwise.

GEN divsi_rem(long s, GEN y, long *r) computes the Euclidean quotient and remainder of the t_long s by the GEN y. As sdivss_rem otherwise.

GEN divss_rem(long x, long y, long *r) computes the Euclidean quotient and remainder of the t_long x by the long y. As sdivss_rem otherwise.

GEN truedvmdii(GEN x, GEN y, GEN *r), as dvmdii but with a non-negative remainder.

Square root and remainder

GEN sqrtremi(GEN N, GEN *r), returns the integer square root S of the non-negative t_INT N (rounded towards 0) and puts the remainder R into *r. Precisely, N = S^2 + R with 0 <= R <= 2S. If r is equal to NULL, the remainder is not created. In the generic case, the remainder is created after the quotient and can be disposed of individually with cgiv(R). If the remainder is 0 set R = gen_0.

Uses a divide and conquer algorithm (discrete variant of Newton iteration) due to Paul Zimmermann (``Karatsuba Square Root'', INRIA Research Report 3805 (1999)).

GEN sqrti(GEN N), returns the integer square root S of the non-negative t_INT N (rounded towards 0). This is identical to sqrtremi(N, NULL).

Pseudo-random integers

long random_bits(long k) returns a random 0 <= x < 2^k. Assumes that 0 <= k < 31.

long pari_rand31(long k) as random_bits with k = 31.

GEN randomi(GEN n) returns a random t_INT between 0 and n - 1. The result is pasted from successive calls to pari_rand31.

Modular operations.

In this subsection, all GENs are t_INT.

ulong Fp_powu(GEN x, ulong n, GEN m) raises x to the n-th power modulo p (smallest non-negative residue).

GEN Fp_pow(GEN x, GEN n, GEN m) returns x^n modulo p (smallest non-negative residue).

GEN Fp_inv(GEN a, GEN m) returns an inverse of a modulo m (smallest non-negative residue). Raise an error if a is not invertible.

GEN Fp_invsafe(GEN a, GEN m) as Fp_inv, but return NULL if a is not invertible.

int invmod(GEN a, GEN m, GEN *g), return 1 if a modulo m is invertible, else return 0 and set g = gcd(a,m).

GEN Fp_sqrt(GEN x, GEN p) returns a square root of x modulo p (the smallest non-negative residue), where x, p are t_INTs, and p is assumed to be prime. Return NULL if x is not a quadratic residue modulo p.

GEN Fp_sqrtn(GEN x, GEN n, GEN p, GEN *zn) returns an n-th root of x modulo p (smallest non-negative residue), where x, n, p are t_INTs, and p is assumed to be prime. Return NULL if x is not an n-th power residue. Otherwise, if zn is non-NULL set it to a primitive n-th root of 1.

long kross(long x, long y) returns the Kronecker symbol (x|y), i.e.-1, 0 or 1. If y is an odd prime, this is the Legendre symbol. (Contrary to krouu, kross also supports y = 0)

long krois(GEN x, long y) returns the Kronecker symbol (x|y) of t_INT x and long y. As kross otherwise.

long krosi(long x, GEN y) returns the Kronecker symbol (x|y) of long x and t_INT y. As kross otherwise.

long kronecker(GEN x, GEN y) returns the Kronecker symbol (x|y) of t_INTs x and y. As kross otherwise.

GEN gener_Fp(GEN p) returns a primitive root modulo p, assuming p is prime.

GEN gener_Fp_local(GEN p, GEN L), L being a vector of primes dividing p - 1, returns an integer x which is a generator of the \ell-Sylow of F_p^* for every \ell in L. In other words, x^{(p-1)/\ell} != 1 for all such \ell. In particular, returns Fp_gener(p) if L contains all primes dividing p - 1.

Miscellaneous functions

void addumului(ulong a, ulong b, GEN x) return a + b|X|.

long cgcd(long x, long y), returns the GCD of the t_longs x and y.

long cbezout(long a,long b, long *u,long *v), returns the GCD d of a and b and sets u, v to the Bezout coefficients such that au + bv = d.

GEN bezout(GEN a,GEN b, GEN *u,GEN *v), returns the GCD d of t_INTs a and b and sets u, v to the Bezout coefficients such that au + bv = d.

GEN factoru(ulong n), returns the factorization of n. The result is a 2-component vector [P,E], where P and E are t_VECSMALL containing the prime divisors of n, and the v_p(n).

GEN factoru_pow(ulong n), returns the factorization of n. The result is a 3-component vector [P,E,C], where P, E and C are t_VECSMALL containing the prime divisors of n, the v_p(n) and the p^{v_p(n)}.

GEN gcdii(GEN x, GEN y), returns the GCD of the t_INTs x and y.

GEN lcmii(GEN x, GEN y), returns the LCM of the t_INTs x and y.

long maxss(long x, long y), return the largest of x and y.

long minss(long x, long y), return the smallest of x and y.

GEN powuu(ulong n, ulong k), returns n^k.

GEN powiu(GEN n, ulong k), assumes n is a t_INT and returns n^k.

ulong upowuu(ulong n, ulong k), returns n^k modulo 2^BIL. This is meant to be used for tiny k, where in fact n^k fits into an ulong.

void rdivii(GEN x, GEN y, long prec), assuming x and y are both of type t_INT, return the quotient x/y as a t_REAL of precision prec.

void rdivis(GEN x, long y, long prec), assuming x is of type t_INT, return the quotient x/y as a t_REAL of precision prec.

void rdivsi(long x, GEN y, long prec), assuming y is of type t_INT, return the quotient x/y as a t_REAL of precision prec.

void rdivss(long x, long y, long prec), return the quotient x/y as a t_REAL of precision prec.


Level 2 kernel (modular arithmetic)

These routines implement univariate polynomial arithmetic and linear algebra over finite fields, in fact over finite rings of the form (Z/pZ)[X]/(T), where p is not necessarily prime and T belongs to (Z/pZ)[X] is possibly reducible; and finite extensions thereof. All this can be emulated with t_INTMOD and t_POLMOD coefficients and using generic routines, at a considerable loss of efficiency. Also, some specialized routines are available that have no obvious generic equivalent.

Naming scheme. A function name is built in the following way

A_1_..._A_nfun for an operation fun with n arguments of class A_1,..., A_n. A class name is given by a base ring followed by a number of code letters. Base rings are among

Fl: Z/lZ where l < 2^{BIL} is not necessarily prime. Implemented using ulongs

Fp: Z/pZ where p is a t_INT, not necessarily prime. Implemented as t_INTs z, preferably satisfying 0 <= z < p. More precisely, any t_INT can be used as an Fp, but reduced inputs are treated more efficiently. Outputs from Fpxxx routines are reduced.

Fq: Z[X]/(p,T(X)), p a t_INT, T a t_POL with Fp coefficients or NULL (in which case no reduction modulo T is performed). Implemented as t_POLs z with Fp coefficients, deg (z) < deg T.

Z: the integers Z, implemented as t_INTs.

z: the integers Z, implemented using (signed) longs.

Q: the rational numbers Q, implemented as t_INTs and t_FRACs.

Rg: a commutative ring, whose elements can be gadd-ed, gmul-ed, etc.

Possible letters are:

X: polynomial in X (t_POL in a fixed variable), e.g. FpX means Z/pZ[X]

Y: polynomial in Y != X. E.g. FpXY means ((Z/pZ)[Y])[X]

V: vector (t_VEC or t_COL), treated as a line vector (independantly of the actual type). E.g. ZV means Z^k for some k.

C: vector (t_VEC or t_COL), treated as a column vector (independantly of the actual type). The difference with V is purely semantic.

M: matrix (t_MAT). E.g. QM means a matrix with rational entries

Q: representative (t_POL) of a class in a polynomial quotient ring. E.g. an FpXQ belongs to (Z/pZ)[X]/(T(X)), FpXQV means a vector of such elements, etc.

x, m, v, c, q: as their uppercase counterpart, but coefficient arrays are implemented using t_VECSMALLs, which coefficient understood as ulongs.

x (and q) are implemented by a t_VECSMALL whose first coefficient is used as a code-word and the following are the coefficients , similarly to a t_POL. This is known as a 'POLSMALL'.

  C<m> are implemented by a C<t_MAT> whose components (columns) are
C<t_VECSMALL>s. This is know as a 'MATSMALL'.

v and c are regular t_VECSMALLs. Difference between the two is purely semantic.

Omitting the letter means the argument is a scalar in the base ring. Standard functions fun are

add: add

sub: subtract

mul: multiply

sqr: square

div: divide (Euclidean quotient)

rem: Euclidean remainder

divrem: return Euclidean quotient, store remainder in a pointer argument.

gcd: GCD

extgcd: return GCD, store Bezout coefficients in pointer arguments

pow: exponentiate

compo: composition

ZX, ZV, ZM.

A ZV (resp. a ZM, resp. a ZX) is a t_VEC or t_COL (resp. t_MAT, resp. t_POL) with t_INT coefficients.

ZV
GEN ZV_add(GEN x, GEN y) adds x and y.

GEN ZV_sub(GEN x, GEN y) subtracts x and y.

ZX
GEN ZX_renormalize(GEN x, long l), as normalizepol, where l = lg(x), in place.

GEN ZX_add(GEN x,GEN y) adds x and y.

GEN ZX_sub(GEN x,GEN y) subtracts x and y.

GEN ZX_neg(GEN x,GEN p) returns -x.

GEN ZX_Z_add(GEN x,GEN y) adds the integer y to the polynomial x.

GEN ZX_Z_mul(GEN x,GEN y) multiplies the polynomial x by the integer y.

GEN ZX_mul(GEN x,GEN y) multiplies x and y.

GEN ZX_sqr(GEN x,GEN p) returns x^2.

GEN ZX_caract(GEN T, GEN A, long v) returns the characteristic polynomial of Mod(A, T), where T is a ZX, A is a ZX. More generally, A is allowed to be a QX, hence possibly has rational coefficients, assuming the result is a ZX, i.e. the algebraic number Mod(A,T) is integral over Z.

GEN ZX_disc(GEN T) returns the discriminant of the ZX T.

int ZX_is_squarefree(GEN T) returns 1 if the ZX T is squarefree, 0 otherwise.

GEN ZX_resultant(GEN A, GEN B) returns the resultant of the ZX A and B.

GEN ZX_QX_resultant(GEN A, GEN B) returns the resultant of the ZX A and the QX B, assuming the result is an integer.

GEN ZY_ZXY_resultant(GEN A, GEN B) under the assumption that A in Z[Y], B in Q[Y][X], and R = {Res}_Y(A, B) belongs to Z[X], returns the resultant R.

GEN ZY_ZXY_rnfequation(GEN A, GEN B, long *lambda), assume A in Z[Y], B in Q[Y][X], and R = {Res}_Y(A, B) belongs to Z[X]. If lambda = NULL, returns R as in ZY_ZXY_resultant. Otherwise, lambda must point to some integer, e.g. 0 which is used as a seed. The function then finds a small lambda belongs to Z (starting from *lambda) such that R_lambda(X) := {Res}_Y(A, B(X + lambda Y)) is squarefree, resets *lambda to the chosen value and returns R_{lambda}.

GEN ZM_inv(GEN M, GEN d) if M is a ZM and d is a t_INT such that M' := dM^{-1} is integral, return M'. It is allowed to set d = NULL, in which case, the determinant of M is computed and used instead.

GEN QM_inv(GEN M, GEN d) as above, with M a QM. We still assume that M' has integer coefficients.

FpX.

Let p an understood t_INT, to be given in the function arguments; in practice p is not assumed to be prime, but be wary. An Fp object is a t_INT belonging to [0, p-1], an FpX is a t_POL in a fixed variable whose coefficients are Fp objects. Unless mentionned otherwise, all outputs in this section are FpXs. All operations are understood to take place in (Z/pZ)[X].

Basic operations. In what follows p is always a t_INT,
not necessarily prime.

GEN Rg_to_Fp(GEN z, GEN p), z a scalar which can be mapped to Z/pZ: a t_INT, a t_INTMOD whose modulus is divisible by p, a t_FRAC whose denominator is coprime to p, or a t_PADIC with underlying prime p. Returns lift(z * Mod(1,p)), normalized.

GEN RgX_to_FpX(GEN z, GEN p), z a t_POL, returns the FpX obtained by applying Rg_to_Fp coefficientwise.

GEN RgV_to_FpV(GEN z, GEN p), z a t_VEC or t_COL, returns the FpV (as a t_VEC) obtained by applying Rg_to_Fp coefficientwise.

GEN RgC_to_FpC(GEN z, GEN p), z a t_VEC or t_COL, returns the FpC (as a t_COL) obtained by applying Rg_to_Fp coefficientwise.

GEN FpX_to_mod(GEN z, GEN p), z a ZX. Returns z * Mod(1,p), normalized. Hence the returned value has t_INTMOD coefficients.

GEN FpX_red(GEN z, GEN p), z a ZX, returns lift(z * Mod(1,p)), normalized.

GEN FpXV_red(GEN z, GEN p), z a t_VEC of ZX. Applies FpX_red componentwise and returns the result (and we obtain a vector of FpXs).

Now, except for p, the operands and outputs are all FpX objects. Results are undefined on other inputs.

GEN FpX_add(GEN x,GEN y, GEN p) adds x and y.

GEN FpX_neg(GEN x,GEN p) returns -x.

GEN FpX_renormalize(GEN x, long l), as normalizepol, where l = lg(x), in place.

GEN FpX_sub(GEN x,GEN y,GEN p) subtracts y from x.

GEN FpX_mul(GEN x,GEN y,GEN p) multiplies x and y.

GEN FpX_sqr(GEN x,GEN p) returns x^2.

GEN FpX_divrem(GEN x, GEN y, GEN p, GEN *pr) returns the quotient of x by y, and sets pr to the remainder.

GEN FpX_div(GEN x, GEN y, GEN p) returns the quotient of x by y.

GEN FpX_div_by_X_x(GEN A, GEN a, GEN p, GEN *r) returns the quotient of the FpX A by (X - a), and sets r to the remainder A(a).

GEN FpX_rem(GEN x, GEN y, GEN p) returns the remainder x mod y

GEN FpX_gcd(GEN x, GEN y, GEN p) returns a (not necessarily monic) greatest common divisor of x and y.

GEN FpX_extgcd(GEN x, GEN y, GEN p, GEN *u, GEN *v) returns d = {GCD}(x,y), and sets *u, *v to the Bezout coefficients such that *ux + *vy = d.

GEN FpX_center(GEN z, GEN p) returns the polynomial whose coefficient belong to the symmetric residue system (clean version of centermod, which assumes the coefficients already belong to [0,p-1]).

Miscellaneous operations
GEN FpX_normalize(GEN z, GEN p) divides the FpX z by its leading coefficient. If the latter is 1, z itself is returned, not a copy. If not, the inverse remains uncollected on the stack.

GEN FpX_Fp_add(GEN y, GEN x, GEN p) add the Fp x to the FpX y, possibly modifying the argument y (thus the operation uses constant time instead of linear linear). This function is not suitable for gerepileupto nor for gerepile.

GEN FpX_Fp_mul(GEN y, GEN x, GEN p) multiplies the FpX y by the Fp x.

GEN FpX_rescale(GEN P, GEN h, GEN p) returns h^{ deg (P)} P(x/h). P is an FpX and h is a non-zero Fp (the routine would work with any non-zero t_INT but is not efficient in this case).

GEN FpX_eval(GEN x, GEN y, GEN p) evaluates the FpX x at the Fp y. The result is an Fp.

GEN FpXV_FpC_mul(GEN V, GEN W, GEN p) multiplies a line vector ofFpX by a column vector of Fp (scalar product). The result is an FpX.

GEN FpXV_prod(GEN V, GEN p), V being a vector of FpX, returns their product.

GEN FpV_roots_to_pol(GEN V, GEN p, long v), V being a vector of INTs, returns the monic FpX prod_i (pol_x[v] - V[i]).

GEN FpX_chinese_coprime(GEN x,GEN y, GEN Tx,GEN Ty, GEN Tz, GEN p) returns an FpX, congruent to x mod Tx and to y mod Ty. Assumes Tx and Ty are coprime, and Tz = Tx * Ty or NULL (in which case it is computed within).

GEN FpV_polint(GEN x, GEN y, GEN p) returns the FpX interpolation polynomial with value y[i] at x[i]. Assumes lengths are the same, components are t_INTs, and the x[i] are distinct modulo p.

long FpX_is_squarefree(GEN f, GEN p) returns 1 if the FpX f is squarefree, 0 otherwise.

long FpX_is_irred(GEN f, GEN p) returns 1 if the FpX f is irreducible, 0 otherwise. Assumes that p is prime. If f has few factors, FpX_nbfact(f,p) == 1 is much faster.

long FpX_is_totally_split(GEN f, GEN p) returns 1 if the FpX f splits into a product of distinct linear factors, 0 otherwise. Assumes that p is prime.

GEN FpX_factor(GEN f, GEN p), factors the FpX f. Assumes that p is prime. The returned value v has two components: v[1] is a vector of distinct irreducible (FpX) factors, and v[2] is a t_VECSMALL of corresponding exponents. The order of the factors is deterministic (the computation is not).

long FpX_nbfact(GEN f, GEN p), assuming the FpX f is squarefree, returns the number of its irreducible factors. Assumes that p is prime.

long FpX_degfact(GEN f, GEN p), as FpX_factor, but the degrees of the irreducible factors are returned instead of the factors themselves (as a t_VECSMALL). Assumes that p is prime.

long FpX_nbroots(GEN f, GEN p) returns the number of distinct roots in Z/pZ of the FpX f. Assumes that p is prime.

GEN FpX_roots(GEN f, GEN p) returns the roots in Z/pZ of the FpX f (without multiplicity, as a vector of Fps). Assumes that p is prime.

GEN FpX_rand(long d, long v, GEN p) returns a random FpX in variable v, of degree less than d.

GEN FpX_resultant(GEN x, GEN y, GEN p) returns the resultant of x and y, both FpX. The result is a t_INT belonging to [0,p-1].

GEN FpY_FpXY_resultant(GEN a, GEN b, GEN p), a a t_POL of t_INTs (say in variable Y), b a t_POL (say in variable X) whose coefficients are either t_POLs in Z[Y] or t_INTs. Returns {Res}_Y(a, b), which is an FpX. The function assumes that Y has lower priority than X.

FpXQ, Fq.

Let p a t_INT and T an FpX for p, both to be given in the function arguments; an FpXQ object is an FpX whose degree is strictly less than the degree of T. An Fq is either an FpXQ or an Fp. Both represent a class in (Z/pZ)[X] / (T), in which all operations below take place. In addition, Fq routines also allow T = NULL, in which case no reduction mod T is performed on the result.

For efficiency, the routines in this section may leave small unused objects behind on the stack (their output is still suitable for gerepileupto). Besides T and p, arguments are either FpXQ or Fq depending on the function name. (All Fq routines accept FpXQs by definition, not the other way round.)

GEN Rg_to_FpXQ(GEN z, GEN T, GEN p), z a GEN which can be mapped to F_p[X]/(T): anything Rg_to_Fp can be applied to, a t_POL to which RgX_to_FpX can be applied to, a t_POLMOD whose modulus is divisible by T (once mapped to a FpX), a suitable t_RFRAC. Returns z as an FpXQ, normalized.

GEN RgX_to_FpXQX(GEN z, GEN T, GEN p), z a t_POL, returns the FpXQ obtained by applying Rg_to_FpXQ coefficientwise.

GEN RgX_to_FqX(GEN z, GEN T, GEN p), z a t_POL, returns the FpXQ obtained by applying Rg_to_FpXQ coefficientwise and simplifying scalars to t_INTs.

GEN Fq_red(GEN x, GEN T, GEN p), x a ZX or t_INT, reduce it to an Fq (T = NULL is allowed iff x is a t_INT).

GEN FqX_red(GEN x, GEN T, GEN p), x a t_POL whose coefficients are ZXs or t_INTs, reduce them to Fqs. (If T = NULL, as FpXX_red(x, p).)

GEN FqV_red(GEN x, GEN T, GEN p), x a vector of ZXs or t_INTs, reduce them to Fqs. (If T = NULL, only reduce components mod p to FpXs or Fps.)

GEN FpXQ_mul(GEN y, GEN x, GEN T,GEN p)

GEN FpXQ_sqr(GEN y, GEN T, GEN p)

GEN FpXQ_div(GEN x, GEN y, GEN T,GEN p)

GEN FpXQ_inv(GEN x, GEN T, GEN p) computes the inverse of x

GEN FpXQ_invsafe(GEN x,GEN T,GEN p), as FpXQ_inv, returning NULL if x is not invertible.

GEN FpXQ_pow(GEN x, GEN n, GEN T, GEN p) computes x^n.

GEN Fq_add(GEN x, GEN y, GEN T/*unused*/, GEN p)

GEN Fq_sub(GEN x, GEN y, GEN T/*unused*/, GEN p)

GEN Fq_mul(GEN x, GEN y, GEN T, GEN p)

GEN Fq_neg(GEN x, GEN T, GEN p)

GEN Fq_neg_inv(GEN x, GEN T, GEN p) computes -x^{-1}

GEN Fq_inv(GEN x, GEN pol, GEN p) computes x^{-1}, raising an error if x is not invertible.

GEN Fq_invsafe(GEN x, GEN pol, GEN p) as Fq_inv, but returns NULL if x is not invertible.

GEN Fq_pow(GEN x, GEN n, GEN pol, GEN p) returns x^n.

GEN FpXQ_charpoly(GEN x, GEN T, GEN p) returns the characteristic polynomial of x

GEN FpXQ_minpoly(GEN x, GEN T, GEN p) returns the minimal polynomial of x

GEN FpXQ_powers(GEN x, long n, GEN T, GEN p) returns [x^0, ..., x^n] as a t_VEC of FpXQs.

GEN FpX_FpXQ_compo(GEN f,GEN x,GEN T,GEN p) returns f(x).

GEN FpX_FpXQV_compo(GEN f,GEN V,GEN T,GEN p) returns f(x), assuming that V was computed by FpXQ_powers(x, n, T, p).

FpXX.

Contrary to what the name implies, an FpXQX is a t_POL whose coefficients are either t_INTs or t_FpXs. This reduce memory overhead at the expense of consistency.

GEN FpXX_add(GEN x, GEN y, GEN p) adds x and y.

GEN FpXX_red(GEN z, GEN p), z a t_POL whose coefficients are either ZXs or t_INTs. Returns the t_POL equal to z with all components reduced modulo p.

GEN FpXX_renormalize(GEN x, long l), as normalizepol, where l = lg(x), in place.

FpXQX, FqX.

Contrary to what the name implies, an FpXQX is a t_POL whose coefficients are Fqs. So the only difference between FqX and FpXQX routines is that T = NULL is not allowed in the latter. (It was thought more useful to allow t_INT components than to enforce strict consistency, which would not imply any efficiency gain.)

Basic operations
GEN FqX_mul(GEN x, GEN y, GEN T, GEN p)

GEN FqX_Fq_mul(GEN P, GEN U, GEN T, GEN p) multiplies the FqX y by the Fq x.

GEN FqX_normalize(GEN z, GEN T, GEN p) divides the FqX z by its leading term.

GEN FqX_sqr(GEN x, GEN T, GEN p)

GEN FqX_divrem(GEN x, GEN y, GEN T, GEN p, GEN *z)

GEN FqX_div(GEN x, GEN y, GEN T, GEN p)

GEN FqX_rem(GEN x, GEN y, GEN T, GEN p)

GEN FqX_gcd(GEN P, GEN Q, GEN T, GEN p)

GEN FpXQX_red(GEN z, GEN T, GEN p) z a t_POL whose coefficients are ZXs or t_INTs, reduce them to FpXQs.

GEN FpXQX_mul(GEN x, GEN y, GEN T, GEN p)

GEN FpXQX_sqr(GEN x, GEN T, GEN p)

GEN FpXQX_divrem(GEN x, GEN y, GEN T, GEN p, GEN *pr)

GEN FpXQX_gcd(GEN x, GEN y, GEN T, GEN p)

GEN FpXQX_extgcd(GEN x, GEN y, GEN T, GEN p, GEN *ptu, GEN *ptv)

GEN FpXQYQ_pow(GEN x, GEN n, GEN S, GEN T, GEN p), x and T being FpXQXs, returns x^n modulo S.

GEN FpXQXV_prod(GEN V, GEN T, GEN p), V being a vector of FpXQX, returns their product.

GEN FqV_roots_to_pol(GEN V, GEN T, GEN p, long v), V being a vector of Fqs, returns the monic FqX prod_i (pol_x[v] - V[i]).

Miscellaneous operations
GEN init_Fq(GEN p, long n, long v) returns an irreducible polynomial of degree n over F_p, in variable v.

long FqX_is_squarefree(GEN P, GEN T, GEN p)

GEN FqX_factor(GEN x, GEN T, GEN p) same output convention as FpX_factor. Assumes p is prime and T irreducible in F_p[X].

GEN FpX_factorff_irred(GEN P, GEN T, GEN p). Assumes p prime and T irreducible in F_p[X]. P being an irreducible FpX, factors it over the finite field F_p[Y]/(T(Y)) and returns the vector of irreducible FqXs factors (the exponents, being all equal to 1, are not included).

GEN FpX_ffisom(GEN P, GEN Q, GEN p). Assumes p prime, P, Q are ZXs, both irreducible mod p, and deg (P) | deg Q. Outputs a monomorphism between F_p[X]/(P) and F_p[X]/(Q), as a polynomial R such that Q | P(R) in F_p[X]. If P and Q have the same degree, it is of course an isomorphism.

void FpX_ffintersect(GEN P, GEN Q, long n, GEN p, GEN *SP,GEN *SQ, GEN MA,GEN MB)\hfil

Assumes p is prime, P, Q are ZXs, both irreducible mod p, and n divides both the degree of P and Q. Compute SP and SQ such that the subfield of F_p[X]/(P) generated by SP and the subfield of F_p[X]/(Q) generated by SQ are isomorphic of degree n. The polynomials P and Q do not need to be of the same variable. If MA (resp. MB) is not NULL, it must be the matrix of the Frobenius map in F_p[X]/(P) (resp. F_p[X]/(Q)).

GEN FpXQ_ffisom_inv(GEN S, GEN T, GEN p). Assumes p is prime, T a ZX, which is irreducible modulo p, S a ZX representing an automorphism of F_q := F_p[X]/(T). (S(X) is the image of X by the automorphism.) Returns the inverse automorphism of S, in the same format, i.e. an FpX H such that H(S) = X modulo (T, p).

long FqX_nbfact(GEN u, GEN T, GEN p). Assumes p is prime and T irreducible in F_p[X].

long FqX_nbroots(GEN f, GEN T, GEN p) Assumes p is prime and T irreducible in F_p[X].

FpV, FpM, FqM.

A ZV (resp. a ZM) is a t_VEC or t_COL (resp. t_MAT) with t_INT coefficients. An FpV or FpM, with respect to a given t_INT p, is the same with Fp coordinates; operations are understood over Z/pZ. An FqM is a matrix with Fq coefficients (with respect to given T, p), not necessarily reduced (i.e arbitrary t_INTs and ZXs in the same variable as T).

Basic operations
GEN FpC_to_mod(GEN z, GEN p), z a ZC. Returns Col(z) * Mod(1,p), hence a t_COL with t_INTMOD coefficients.

GEN FpV_to_mod(GEN z, GEN p), z a ZV. Returns Vec(z) * Mod(1,p), hence a t_VEC with t_INTMOD coefficients.

GEN FpM_to_mod(GEN z, GEN p), z a ZM. Returns z * Mod(1,p), hence with t_INTMOD coefficients.

GEN FpC_red(GEN z, GEN p), z a ZC. Returns lift(Col(z) * Mod(1,p)), hence a t_COL.

GEN FpV_red(GEN z, GEN p), z a ZV. Returns lift(Vec(z) * Mod(1,p)), hence a t_VEC

GEN FpM_red(GEN z, GEN p), z a ZM. Returns lift(z * Mod(1,p)), which is an FpM.

GEN FpC_Fp_mul(GEN x, GEN y, GEN p) multiplies the ZC x (seen as a column vector) by the t_INT y and reduce modulo p to obtain an FpC.

GEN FpC_FpV_mul(GEN x, GEN y, GEN p) multiplies the ZC x (seen as a column vector) by the ZV y (seen as a row vector, assumed to have compatible dimensions), and reduce modulo p to obtain an FpM.

GEN FpM_mul(GEN x, GEN y, GEN p) multiplies the two ZMx and y (assumed to have compatible dimensions), and reduce modulo p to obtain an FpM.

GEN FpM_FpC_mul(GEN x, GEN y, GEN p) multiplies the ZM x by the ZC y (seen as a column vector, assumed to have compatible dimensions), and reduce modulo p to obtain an FpC.

GEN FpV_FpC_mul(GEN x, GEN y, GEN p) multiplies the ZV x (seen as a row vector) by the ZC y (seen as a column vector, assumed to have compatible dimensions), and reduce modulo p to obtain an Fp.

Fp-linear algebra. The implementations are not
asymptotically efficient (O(n^3) standard algorithms).

GEN FpM_deplin(GEN x, GEN p) returns a non-trivial kernel vector, or NULL if none exist.

GEN FpM_gauss(GEN a, GEN b, GEN p) as gauss

GEN FpM_image(GEN x, GEN p) as image

GEN FpM_intersect(GEN x, GEN y, GEN p) as intersect

GEN FpM_inv(GEN x, GEN p) returns the inverse of x, or NULL if x is not invertible.

GEN FpM_invimage(GEN m, GEN v, GEN p) as inverseimage

GEN FpM_ker(GEN x, GEN p) as ker

long FpM_rank(GEN x, GEN p) as rank

GEN FpM_indexrank(GEN x, GEN p) as indexrank but returns a t_VECSMALL

GEN FpM_suppl(GEN x, GEN p) as suppl

Fq-linear algebra
GEN FqM_gauss(GEN a, GEN b, GEN T, GEN p) as gauss

GEN FqM_ker(GEN x, GEN T, GEN p) as ker

GEN FqM_suppl(GEN x, GEN T, GEN p) as suppl

Flx

Let p an understood ulong, assumed to be prime, to be given the the function arguments; an Fl is an ulong belonging to [0,p-1], an Flx z is a t_VECSMALL representing a polynomial with small integer coefficients. Specifically z[0] is the usual codeword, z[1] = evalvarn(v) for some variable v, then the coefficients by increasing degree. An FlxX is a t_POL whose coefficients are Flxs.

In the following, an argument called sv is of the form evalvarn(v) for some variable number v.

Basic operations
ulong Rg_to_Fl(GEN z, ulong p), z which can be mapped to Z/pZ: a t_INT, a t_INTMOD whose modulus is divisible by p, a t_FRAC whose denominator is coprime to p, or a t_PADIC with underlying prime p. Returns lift(z * Mod(1,p)), normalized, as an Fl.

GEN Flx_red(GEN z, ulong p) converts from zx with non-negative coefficients to Flx (by reducing them mod p).

GEN Flx_add(GEN x, GEN y, ulong p)

GEN Flx_neg(GEN x, ulong p)

GEN Flx_neg_inplace(GEN x, ulong p), same as Flx_neg, in place (x is destroyed).

GEN Flx_sub(GEN x, GEN y, ulong p)

GEN Flx_mul(GEN x, GEN y, ulong p)

GEN Flx_sqr(GEN x, ulong p)

GEN Flx_divrem(GEN x, GEN y, ulong p, GEN *pr)

GEN Flx_div(GEN x, GEN y, ulong p)

GEN Flx_rem(GEN x, GEN y, ulong p)

GEN Flx_deriv(GEN z, ulong p)

GEN Flx_gcd(GEN a, GEN b, ulong p) returns a (not necessarily monic) greatest common divisor of x and y.

GEN Flx_gcd_i(GEN a, GEN b, ulong p), same as Flx_gcd without collecting garbage.

GEN Flx_extgcd(GEN a, GEN b, ulong p, GEN *ptu, GEN *ptv)

GEN Flx_pow(GEN x, long n, ulong p)

Miscellaneous operations
GEN Flx_normalize(GEN z, ulong p), as FpX_normalize.

GEN Flx_Fl_mul(GEN y, ulong x, ulong p)

GEN Flx_recip(GEN x), returns the reciprocal polynomial

ulong Flx_resultant(GEN a, GEN b, ulong p), returns the resultant of a and b

ulong Flx_extresultant(GEN a, GEN b, ulong p, GEN *ptU, GEN *ptV) returns the resultant and sets Bezout coefficients (if the resultant is 0, the latter are not set).

GEN Flx_invmontgomery(GEN T, ulong p), returns the Montgomery inverse of T, i.e. truncate(x / polrecip(T)+O(x^n). Assumes T(0) != 0.

GEN Flx_rem_montgomery(GEN x, GEN mg, GEN T, ulong p), returns x modulo T, where mg is the Montgomery inverse of T.

GEN Flx_renormalize(GEN x, long l), as FpX_renormalize, where l = lg(x), in place.

GEN Flx_shift(GEN T, long n), returns T multiplied by x^n.

long Flx_valuation(GEN x) returns the valuation of x, i.e. the multiplicity of the 0 root.

GEN FlxYqQ_pow(GEN x, GEN n, GEN S, GEN T, ulong p), as FpXQYQ_pow.

GEN Flx_div_by_X_x(GEN A, ulong a, ulong p, ulong *rem), returns the Euclidean quotient of the Flx A by X - a, and sets rem to the remainder A(a).

ulong Flx_eval(GEN x, ulong y, ulong p), as FpX_eval.

GEN FlxV_Flc_mul(GEN V, GEN W, ulong p), as FpXV_FpC_mul.

int Flx_is_squarefree(GEN z, ulong p)

long Flx_nbfact(GEN z, ulong p), as FpX_nbfact.

long Flx_nbroots(GEN f, ulong p), as FpX_nbroots.

GEN Flv_polint(GEN x, GEN y, ulong p, long sv) as FpV_polint, returning an Flx in variable v.

GEN Flv_roots_to_pol(GEN a, ulong p, long sv) as FpV_roots_to_pol returning an Flx in variable v.

Flxq.

See FpXQ operations.

GEN Flxq_mul(GEN y, GEN x, GEN pol, ulong p)

GEN Flxq_sqr(GEN y, GEN pol, ulong p)

GEN Flxq_inv(GEN x, GEN T, ulong p)

GEN Flxq_invsafe(GEN x, GEN T, ulong p)

GEN Flxq_pow(GEN x, GEN n, GEN pol, ulong p)

GEN Flxq_powers(GEN x, long l, GEN T, ulong p)

GEN FlxqV_roots_to_pol(GEN V, GEN T, ulong p, long v) as FqV_roots_to_pol returning an FlxqX in variable v.

FlxX.

See FpXX operations.

GEN FlxX_add(GEN P, GEN Q, ulong p)

GEN FlxX_renormalize(GEN x, long l), as normalizepol, where l = lg(x), in place.

GEN FlxX_shift(GEN a, long n)

FlxqX.

See FpXQX operations.

GEN FlxqX_mul(GEN x, GEN y, GEN T, ulong p)

GEN FlxqX_Flxq_mul(GEN P, GEN U, GEN T, ulong p)

GEN FlxqX_normalize(GEN z, GEN T, ulong p)

GEN FlxqX_sqr(GEN x, GEN T, ulong p)

GEN FlxqX_divrem(GEN x, GEN y, GEN T, ulong p, GEN *pr)

GEN FlxqX_red(GEN z, GEN T, ulong p)

GEN FlxqXV_prod(GEN V, GEN T, ulong p)

GEN FlxqXQ_pow(GEN x, GEN n, GEN S, GEN T, ulong p)

Flv, Flm.

See FpV, FpM operations.

GEN Flm_Flc_mul(GEN x, GEN y, ulong p)

GEN Flm_deplin(GEN x, ulong p)

GEN Flm_gauss(GEN a, GEN b, ulong p)

GEN Flm_indexrank(GEN x, ulong p)

GEN Flm_inv(GEN x, ulong p)

GEN Flm_ker(GEN x, ulong p)

GEN Flm_ker_sp(GEN x, ulong p, long deplin), as Flm_ker, in place (destroys x).

GEN Flm_mul(GEN x, GEN y, ulong p)

FlxqV, FlxqM.

See FqV, FqM operations.

GEN FlxqM_ker(GEN x, GEN T, ulong p)

QX.

GEN QXQ_inv(GEN A, GEN B) returns the inverse of A modulo B where A and B are QXs.

RgX.

GEN RgX_add(GEN x,GEN y) adds x and y.

GEN RgX_sub(GEN x,GEN y) subtracts x and y.

GEN RgX_neg(GEN x,GEN p) returns -x.

The functions above are currently implemented through the generic routines, but it might change in the future.

GEN RgX_mul(GEN x, GEN y) multiplies the two t_POL (in the same variable) x and y. Uses Karatsuba algorithm.

GEN RgX_mulspec(GEN a, GEN b, long na, long nb). Internal routine: a and b are arrays of coefficients representing polynomials sum_{i = 0}^{na} a[i] X^i and sum_{i = 0}^{nb} b[i] X^i. Returns their product (as a true GEN).

GEN RgX_sqr(GEN x) squares the t_POL x. Uses Karatsuba algorithm.

GEN RgX_sqrspec(GEN a, long na). Internal routine: a is an array of coefficients representing polynomial sum_{i = 0}^{na} a[i] X^i. Return its square (as a true GEN).

GEN RgX_divrem(GEN x, GEN y, GEN *r)

GEN RgX_div(GEN x, GEN y, GEN *r)

GEN RgX_div_by_X_x(GEN A, GEN a, GEN *r) returns the quotient of the RgX A by (X - a), and sets r to the remainder A(a).

GEN RgX_rem(GEN x, GEN y, GEN *r)

GEN RgX_mulXn(GEN x, long n) returns x * t^n. This may be a t_FRAC if n < 0 and the valuation of x is not large enough.

GEN RgX_shift(GEN x, long n) returns x * t^n if n >= 0, and x \t^{-n} otherwise.

GEN RgX_shift_shallow(GEN x, long n) as RgX_shift, but shallow (coefficients are not copied). This is not suitable for gerepile or gerepileupto.

GEN RgX_extgcd(GEN x, GEN y, GEN *u, GEN *v) returns d = {GCD}(x,y), and sets *u, *v to the Bezout coefficients such that *ux + *vy = d.

GEN RgXQ_mul(GEN y, GEN x, GEN T)

GEN RgXQ_norm(GEN x, GEN T) returns the norm of Mod(x, T).

GEN RgXQ_sqr(GEN x, GEN T)

GEN RgXQ_powers(GEN x, long n, GEN T, GEN p) returns [x^0, ..., x^n] as a t_VEC of RgXQs.

GEN RgXQC_red(GEN z, GEN T) z a vector whose coefficients are RgXs (arbitrary GENs in fact), reduce them to RgXQs (applying grem coefficientwise) in a t_COL.

GEN RgXQV_red(GEN z, GEN T) z a t_POL whose coefficients are RgXs (arbitrary GENs in fact), reduce them to RgXQs (applying grem coefficientwise) in a t_VEC.

GEN RgXQX_red(GEN z, GEN T) z a t_POL whose coefficients are RgXs (arbitrary GENs in fact), reduce them to RgXQs (applying grem coefficientwise).

GEN RgXQX_mul(GEN x, GEN y, GEN T)

GEN RgX_Rg_mul(GEN y, GEN x) multiplies the RgX y by the scalar x.

GEN RgX_Rg_div(GEN y, GEN x) divides the RgX y by the scalar x.

GEN RgXQX_RgXQ_mul(GEN x, GEN y, GEN T) multiplies the RgXQX y by the scalar (RgXQ) x.

GEN RgXQX_sqr(GEN x, GEN T)

GEN RgXQX_divrem(GEN x, GEN y, GEN T, GEN *pr)

GEN RgXQX_div(GEN x, GEN y, GEN T, GEN *r)

GEN RgXQX_rem(GEN x, GEN y, GEN T, GEN *r)

GEN RgX_rescale(GEN P, GEN h) returns h^{ deg (P)} P(x/h). P is an RgX and h is non-zero. (Leaves small objects on the stack. Suitable but inefficient for gerepileupto.)

GEN RgX_unscale(GEN P, GEN h) returns P(h x). (Leaves small objects on the stack. Suitable but inefficient for gerepileupto.)

GEN RgXV_unscale(GEN v, GEN h)

Conversions involving single precision objects

To single precision
GEN ZX_to_Flx(GEN x, ulong p) reduce ZX x modulo p (yielding an Flx).

GEN ZV_to_Flv(GEN x, ulong p) reduce ZV x modulo p (yielding an Flv).

GEN ZXV_to_FlxV(GEN v, ulong p), as ZX_to_Flx, repeatedly called on the vector's coefficients.

GEN ZXX_to_FlxX(GEN B, ulong p, long v), as ZX_to_Flx, repeatedly called on the polynomial's coefficients.

GEN ZXXV_to_FlxXV(GEN V, ulong p, long v), as ZXX_to_FlxX, repeatedly called on the vector's coefficients.

GEN ZM_to_Flm(GEN x, ulong p) reduce ZM x modulo p (yielding an Flm).

GEN ZV_to_zv(GEN z), converts coefficients using itos

GEN ZV_to_nv(GEN z), converts coefficients using itou

GEN ZM_to_zm(GEN z), converts coefficients using itos

GEN FqC_to_FlxC(GEN x, GEN T, GEN p), converts coefficients in Fq to coefficient in Flx, result being a column vector.

GEN FqV_to_FlxV(GEN x, GEN T, GEN p), converts coefficients in Fq to coefficient in Flx, result being a line vector.

GEN FqM_to_FlxM(GEN x, GEN T, GEN p), converts coefficients in Fq to coefficient in Flx.

From single precision
GEN Flx_to_ZX(GEN z), converts to ZX (t_POL of non-negative t_INTs in this case)

GEN Flx_to_ZX_inplace(GEN z), same as Flx_to_ZX, in place (z is destroyed).

GEN FlxX_to_ZXX(GEN B), converts an FlxX to a polynomial with ZX or t_INT coefficients (repeated calls to Flx_to_ZX).

GEN FlxC_to_ZXC(GEN x), converts a vector of Flx to a column vector of polynomials with t_INT coefficients (repeated calls to Flx_to_ZX).

GEN FlxM_to_ZXM(GEN z), converts a matrix of Flx to a matrix of polynomials with t_INT coefficients (repeated calls to Flx_to_ZX).

GEN zx_to_ZX(GEN z), as Flx_to_ZX, without assuming coefficients are non-negative.

GEN Flc_to_ZC(GEN z), converts to ZC (t_COL of non-negative t_INTs in this case)

GEN Flv_to_ZV(GEN z), converts to ZV (t_VEC of non-negative t_INTs in this case)

GEN Flm_to_ZM(GEN z), converts to ZM (t_MAT with non-negative t_INTs coefficients in this case)

GEN zc_to_ZC(GEN z) as Flc_to_ZC, without assuming coefficients are non-negative.

GEN zv_to_ZV(GEN z) as Flv_to_ZV, without assuming coefficients are non-negative.

GEN zm_to_ZM(GEN z) as Flm_to_ZM, without assuming coefficients are non-negative.

Mixed precision linear algebra Assumes dimensions are compatible.
Multiply a multiprecision object by a single-precision one.

GEN RgM_zc_mul(GEN x, GEN y)

GEN RgM_zm_mul(GEN x, GEN y)

GEN RgV_zc_mul(GEN x, GEN y)

GEN RgV_zm_mul(GEN x, GEN y)

GEN ZM_zc_mul(GEN x, GEN y)

GEN ZM_zm_mul(GEN x, GEN y)

Miscellaneous
GEN zero_Flx(long sv) returns a zero Flx in variable v.

GEN zero_zx(long sv) as zero_Flx

GEN polx_Flx(long sv) returns the variable v as degree 1 Flx.

GEN polx_zx(long sv) as polx_Flx

GEN Fl_to_Flx(ulong x, long sv) converts a unsigned long to a scalar Flx in shifted variable sv.

GEN Z_to_Flx(GEN x, ulong p, long v) converts a t_INT to a scalar polynomial in variable v.

GEN Flx_to_Flv(GEN x, long n) converts from Flx to Flv with n components (assumed larger than the number of coefficients of x).

GEN zx_to_zv(GEN x, long n) as Flx_to_Flv.

GEN Flv_to_Flx(GEN x, long sv) converts from vector (coefficient array) to (normalized) polynomial in variable v.

GEN zv_to_zx(GEN x, long n) as Flv_to_Flx.

GEN matid_Flm(long n) returns an Flm which is an n x n identity matrix.

GEN Flm_to_FlxV(GEN x, long sv) converts the colums of Flm x to an array of Flx (repeated calls to Flv_to_Flx).

GEN zm_to_zxV(GEN x, long n) as Flm_to_FlxV.

GEN Flm_to_FlxX(GEN x, long sv,long w) converts the colums of Flm x to the coefficient of an FlxX, and normalize the result.

GEN FlxV_to_Flm(GEN v, long n) reverse Flm_to_FlxV, to obtain an Flm with n rows (repeated calls to Flx_to_Flv).

GEN FlxX_to_Flm(GEN v, long n) reverse Flm_to_FlxX, to obtain an Flm with n rows (repeated calls to Flx_to_Flv).


Operations on general PARI objects

Assignment

void gaffsg(long s, GEN x) assigns the long s into the object x.

void gaffect(GEN x, GEN y) assigns the object x into the object y.

Conversions

Scalars
double rtodbl(GEN x) applied to a t_REAL x, converts x into a double if possible.

GEN dbltor(double x) converts the double x into a t_REAL.

double gtodouble(GEN x) if x is a real number (not necessarily a t_REAL), converts x into a double if possible.

long gtolong(GEN x) if x is an integer (not necessarily a t_INT), converts x into a long if possible.

GEN fractor(GEN x, long l) applied to a t_FRAC x, converts x into a t_REAL of length prec.

GEN quadtoc(GEN x, long l) applied to a t_QUAD x, converts x into a t_REAL or t_COMPLEX depending on the sign of the discriminant of x, to precision l BIL-bit words. line brk at hyphen here [GN]

GEN ctofp(GEN x, long prec) converts the t_COMPLEX x to a a complex whose real and imaginary parts are t_REAL of length prec, using gtofp;

GEN gtofp(GEN x, long prec) converts the complex number x (t_INT, t_REAL, t_FRAC, t_QUAD or t_COMPLEX) to either a t_REAL or t_COMPLEX whose components are t_REAL of length prec.

GEN gcvtop(GEN x, GEN p, long l) converts x into a t_PADIC p-adic number of precision l.

GEN gprec(GEN x, long l) returns a copy of x whose precision is changed to l digits. The precision change is done recursively on all components of x. Digits means decimal, p-adic and X-adic digits for t_REAL, t_SER, t_PADIC components, respectively.

GEN gprec_w(GEN x, long l) returns a shallow copy of x whose t_REAL components have their precision changed to l words. This is often more useful than gprec. Shallow copy means that unaffected components are not copied; in particular, this funciton is not suitable for gerepileupto.

GEN gprec_wtrunc(GEN x, long l) returns a shallow copy of x whose t_REAL components have their precision truncated to l words. Contrary to gprec_w, this function may never increase the precision of x. Shallow copy means that unaffected components are not copied; in particular, this funciton is not suitable for gerepileupto.

Modular objects
GEN gmodulo(GEN x, GEN y) creates the object Mod(x,y) on the PARI stack, where x and y are either both t_INTs, and the result is a t_INTMOD, or x is a scalar or a t_POL and y a t_POL, and the result is a t_POLMOD.

GEN gmodulgs(GEN x, long y) same as gmodulo except y is a long.

GEN gmodulss(long x, long y) same as gmodulo except both x and y are longs.

Between polynomials and coefficient arrays
GEN gtopoly(GEN x, long v) converts or truncates the object x into a t_POL with main variable number v. A common application would be the conversion of coefficient vectors (coefficients are given by decreasing degree). E.g. [2,3] goes to 2*v + 3

GEN gtopolyrev(GEN x, long v) converts or truncates the object x into a t_POL with main variable number v, but vectors are converted in reverse order compared to gtopoly (coefficients are given by increasing degree). E.g. [2,3] goes to 3*v + 2. In other words the vector represents a polynomial in the basis (1,v,v^2,v^3,...).

GEN normalizepol(GEN x) applied to an unnormalized t_POL x (with all coefficients correctly set except that leading_term(x) might be zero), normalizes x correctly in place and returns x. For internal use.

The following routines do not copy coefficients on the stack (they only move pointers around), hence are very fast but not suitable for gerepile calls. Recall that an RgV (resp. an RgX, resp. an RgM) is a t_VEC or t_COL (resp. a t_POL, resp. a t_MAT) with arbitrary components. Similarly, an RgXV is a t_VEC or t_COL with RgX components, etc.

GEN RgV_to_RgX(GEN x, long v) converts the RgV x to a (normalized) polynomial in variable v (as gtopolyrev, without copy).

GEN RgX_to_RgV(GEN x, long N) converts the t_POL x to a t_COL v with N components. Other types than t_POL are allowed for x, which is then considered as a constant polynomial. Coefficients of x are listed by increasing degree, so that y[i] is the coefficient of the term of degree i-1 in x.

GEN RgM_to_RgXV(GEN x, long v) converts the RgM x to a t_VEC of RgX, by repeated calls to RgV_to_RgX.

GEN RgXV_to_RgM(GEN v, long N) converts the vector of RgX v to a t_MAT with N rows, by repeated calls to RgX_to_RgV.

GEN RgM_to_RgXX(GEN x, long v,long w) converts the RgM x into a t_POL in variable v, whose coefficients are t_POLs in variable w. This is a shortcut for

    RgV_to_RgX( RgM_to_RgXV(x, w), v );

There are no consistency checks with respect to variable priorities: the above is an invalid object if varncmp(v, w) >= 0.

GEN RgXX_to_RgM(GEN x, long N) converts the t_POL x with RgX (or constant) coefficients to a matrix with N rows.

GEN RgXY_swap(GEN P, long n, long w) converts the bivariate polynomial P(u,v) (a t_POL with t_POL coefficients) to P(pol_x[w],u), assuming n is an upper bound for deg _v(P).

GEN greffe(GEN x, long l, int use_stack) applied to a t_POL x, creates a t_SER of length l starting with x, but without actually copying the coefficients, just the pointers. If use_stack is 0, this is created through malloc, and must be freed after use. Intended for internal use only.

GEN gtoser(GEN x, long v) converts the object x into a t_SER with main variable number v.

GEN gtocol(GEN x) converts the object x into a t_COL

GEN gtomat(GEN x) converts the object x into a t_MAT.

GEN gtovec(GEN x) converts the object x into a t_VEC.

GEN gtovecsmall(GEN x) converts the object x into a t_VECSMALL.

GEN normalize(GEN x) applied to an unnormalized t_SER x (i.e. type t_SER with all coefficients correctly set except that x[2] might be zero), normalizes x correctly in place. Returns x. For internal use.

Clean Constructors

GEN zeropadic(GEN p, long n) creates a 0 t_PADIC equal to O(p^n).

GEN zeroser(long v, long n) creates a 0 t_SER in variable v equal to O(X^n).

GEN scalarser(GEN x, long v, long prec) creates a constant t_SER in variable v and precision prec, whose constant coefficient is (a copy of) x, in other words x + O(v^prec). Assumes that x is non-zero.

GEN zeropol(long v) creates a 0 t_POL in variable v.

GEN scalarpol(GEN x, long v) creates a constant t_POL in variable v, whose constant coefficient is (a copy of) x.

GEN zerocol(long n) creates a t_COL with n components set to gen_0.

GEN zerovec(long n) creates a t_VEC with n components set to gen_0.

GEN col_ei(long n, long i) creates a t_COL with n components set to gen_0, but the i-th one which us set to gen_1 (i-th vector in the canonical basis).

GEN vec_ei(long n, long i) creates a t_VEC with n components set to gen_0, but the i-th one which us set to gen_1 (i-th vector in the canonical basis).

GEN zeromat(long m, long n) creates a t_MAT with m x n components set to gen_0. Note that the result allocates a single column, so modifying an entry in one column modifies it in all columns. To fully allocate a matrix initialized with zero entries, use zeromatcopy.

GEN zeromatcopy(long m, long n) creates a t_MAT with m x n components set to gen_0. Note that

See also next section for analogs of the following functions:

GEN mkcolcopy(GEN x) creates a 1-dimensional t_COL containing x.

GEN mkmatcopy(GEN x) creates a 1-by-1 t_MAT containing x.

GEN mkveccopy(GEN x) creates a 1-dimensional t_VEC containing x.

GEN mkvec2copy(GEN x, GEN y) creates a 2-dimensional t_VEC equal to [x,y].

GEN mkvecs(long x) creates a 1-dimensional t_VEC containing stoi(x).

GEN mkvec2s(long x, long y) creates a 2-dimensional t_VEC containing [stoi(x), stoi(y)].

GEN mkvec3s(long x, long y, long z) creates a 3-dimensional t_VEC containing [stoi(x), stoi(y), stoi(z)].

GEN mkvecsmall(long x) creates a 1-dimensional t_VECSMALL containing x.

GEN mkvecsmall2(long x, long y) creates a 2-dimensional t_VECSMALL containing [x, y].

GEN mkvecsmall3(long x, long y, long z) creates a 3-dimensional t_VECSMALL containing [x, y, z].

Unclean Constructors

Contrary to the policy of general PARI functions, the functions in this subsection do not copy their arguments, nor do they produce an object a priori suitable for gerepileupto. In particular, they are faster than their clean equivalent (which may not exist). If you restrict their arguments to universal objects (e.g gen_0), then the above warning does not apply.

GEN mkcomplex(GEN x, GEN y) creates the t_COMPLEX x + iy.

GEN mkfrac(GEN x, GEN y) creates the t_FRAC x/y. Assumes that y > 1 and (x,y) = 1.

GEN mkrfrac(GEN x, GEN y) creates the t_RFRAC x/y. Assumes that y is a t_POL, x a compatible type whose variable has lower or same priority, with (x,y) = 1.

GEN mkcol(GEN x) creates a 1-dimensional t_COL containing x.

GEN mkintmod(GEN x, GEN y) creates the t_INTMOD Mod(x, y). The input must be t_INTs satisfying 0 <= x < y.

GEN mkpolmod(GEN x, GEN y) creates the t_POLMOD Mod(x, y). The input must satisfy deg x < deg y with respect to the main variable of the t_POL y. x may be a scalar.

GEN mkmat(GEN x) creates a 1-by-1 t_MAT containing x.

GEN mkvec(GEN x) creates a 1-dimensional t_VEC containing x.

GEN mkvec2(GEN x, GEN y) creates a 2-dimensional t_VEC equal to [x,y].

GEN mkvec3(GEN x, GEN y, GEN z) creates a 3-dimensional t_VEC equal to [x,y,z].

GEN mkvec4(GEN x, GEN y, GEN z, GEN t) creates a 4-dimensional t_VEC equal to [x,y,z,t].

GEN mkintn(long n, ...) returns the non-negative t_INT whose development in base 2^{32} is given by the following n words (unsigned long). It is assumed that all such arguments are less than 2^{32} (the actual sizeof(long) is irrelevant, the behaviour is also as above on 64-bit machines).

    mkintn(3, a2, a1, a0);

returns a_2 2^{64} + a_1 2^{32} + a_0.

GEN mkpoln(long n, ...) Returns the t_POL whose n coefficients (GEN) follow, in order of decreasing degree.

    mkpoln(3, gen_1, gen_2, gen_0);

returns the polynomial X^2 + 2X (in variable 0, use setvarn if you want other variable numbers). Beware that n is the number of coefficients, hence one more than the degree.

GEN mkvecn(long n, ...) returns the t_VEC whose n coefficients (GEN) follow.

GEN mkcoln(long n, ...) returns the t_COL whose n coefficients (GEN) follow.

Integer parts

GEN gfloor(GEN x) creates the floor of x, i.e. the (true) integral part.

GEN gfrac(GEN x) creates the fractional part of x, i.e. x minus the floor of x.

GEN gceil(GEN x) creates the ceiling of x.

GEN ground(GEN x) rounds towards + oo the components of x to the nearest integers.

GEN grndtoi(GEN x, long *e) same as ground, but in addition sets *e to the binary exponent of x - ground(x). If this is positive, all significant bits are lost. This kind of situation raises an error message in ground but not in grndtoi.

GEN gtrunc(GEN x) truncates x. This is the false integer part if x is a real number (i.e. the unique integer closest to x among those between 0 and x). If x is a t_SER, it is truncated to a t_POL; if x is a t_RFRAC, this takes the polynomial part.

GEN gcvtoi(GEN x, long *e) same as grndtoi except that rounding is replaced by truncation.

Valuation and shift

GEN gshift[z](GEN x, long n[, GEN z]) yields the result of shifting (the components of) x left by n (if n is non-negative) or right by -n (if n is negative). Applies only to t_INT and vectors/matrices of such. For other types, it is simply multiplication by 2^{n}.

GEN gmul2n[z](GEN x, long n[, GEN z]) yields the product of x and 2^{n}. This is different from gshift when n is negative and x is a t_INT: gshift truncates, while gmul2n creates a fraction if necessary.

long ggval(GEN x, GEN p) returns the greatest exponent e such that p^e divides x, when this makes sense.

long gval(GEN x, long v) returns the highest power of the variable number v dividing the t_POL x.

long polvaluation(GEN P, GEN *z) returns the valuation v of the t_POL P with respect to its main variable X. Check whether coefficients are 0 using gcmp0. If z is non-NULL, set it to P / X^v.

long polvaluation_inexact(GEN P, GEN *z) as polvaluation but use isexactzero instead of gcmp0.

long ZX_valuation(GEN P, GEN *z) as polvaluation, but assumes P has t_INT coefficients.

Comparison operators

int isexactzero(GEN x) returns 1 (true) if x is exactly equal to 0, 0 (false) otherwise. Note that many PARI functions return a pointer to gen_0 when they are aware that the result they return is an exact zero, so it is almost always faster to test for pointer equality first, and call isexactzero (or gcmp0) only when the first test fails.

int isinexact(GEN x) returns 0 (false) if x has an inexact component, and 1 (true) otherwise.

int isint(GEN x, GEN *n) returns 0 (false) if x does not round to an integer. Otherwise, returns 1 (true) and set n to the rounded value.

int issmall(GEN x, long *n) returns 0 (false) if x does not round to a small integer (suitable for itos). Otherwise, returns 1 (true) and set n to the rounded value.

int gcmp0(GEN x) returns 1 (true) if x is equal to 0, 0 (false) otherwise.

int gcmp1(GEN x) returns 1 (true) if x is equal to 1, 0 (false) otherwise.

int gcmp_1(GEN x) returns 1 (true) if x is equal to -1, 0 (false) otherwise.

long gcmp(GEN x, GEN y) comparison of x with y (returns the sign of x-y).

long gcmpsg(long s, GEN x) comparison of the long s with x.

long gcmpgs(GEN x, long s) comparison of x with the long s.

long lexcmp(GEN x, GEN y) comparison of x with y for the lexicographic ordering.

long gequal(GEN x, GEN y) returns 1 (true) if x is equal to y, 0 otherwise. A priori, this makes sense only if x and y have the same type. When the types are different, a true result means that x - y was successfully computed and found equal to 0 (by gcmp0). In particular

    gequal(cgetg(1, t_VEC), gen_0)

is true, and the relation is not transitive. E.g. an empty t_COL and an empty t_VEC are not equal but are both equal to gen_0.

long gequalsg(long s, GEN x) returns 1 (true) if the long s is equal to x, 0 otherwise.

long gequalgs(GEN x, long s) returns 1 (true) if x is equal to the long s, 0 otherwise.

long iscomplex(GEN x) returns 1 (true) if x is a complex number (of component types embeddable into the reals) but is not itself real, 0 if x is a real (not necessarily of type t_REAL), or raises an error if x is not embeddable into the complex numbers.

long ismonome(GEN x) returns 1 (true) if x is a non-zero monomial in its main variable, 0 otherwise.

Generic unary operators

GEN gneg[\key{z}](GEN x[, GEN z]) yields -x.

GEN gabs[\key{z}](GEN x[, GEN z]) yields |x|.

GEN gsqr(GEN x) creates the square of x.

GEN ginv(GEN x) creates the inverse of x.

Divisibility, Euclidean division

GEN gdivexact(GEN x, GEN y) returns the quotient x / y, assuming y divides x.

int gdvd(GEN x, GEN y) returns 1 (true) if y divides x, 0 otherwise.

GEN gdiventres(GEN x, GEN y) creates a 2-component vertical vector whose components are the true Euclidean quotient and remainder of x and y.

GEN gdivent[z](GEN x, GEN y[, GEN z]) yields the true Euclidean quotient of x and the t_INT or t_POL y.

GEN gdiventsg[z](long s, GEN y[, GEN z]), as gdivent except that x is a long.

GEN gdiventgs[z](GEN x, long s[, GEN z]), as gdivent except that y is a long.

GEN gmod[z](GEN x, GEN y[, GEN z]) yields the true remainder of x modulo the t_INT or t_POL y. A t_REAL or t_FRAC y is also allowed, in which case the remainder is the unique real r such that 0 <= r < |y| and y = qx + r for some (in fact unique) integer q.

GEN gmodsg[z](long s, GEN y[, GEN z]) as gmod, except x is a long.

GEN gmodgs[z](GEN x, long s[, GEN z]) as gmod, except y is a long.

GEN gdivmod(GEN x, GEN y, GEN *r) If r is not equal to NULL or ONLY_REM, creates the (false) Euclidean quotient of x and y, and puts (the address of) the remainder into *r. If r is equal to NULL, do not create the remainder, and if r is equal to ONLY_REM, create and output only the remainder. The remainder is created after the quotient and can be disposed of individually with a cgiv(r).

GEN poldivrem(GEN x, GEN y, GEN *r) same as gdivmod but specifically for t_POLx and y, not necessarily in the same variable. Either of x and y may also be scalars (treated as polynomials of degree 0)

GEN gdeuc(GEN x, GEN y) creates the Euclidean quotient of the t_POLx and y. Either of x and y may also be scalars (treated as polynomials of degree 0)

GEN grem(GEN x, GEN y) creates the Euclidean remainder of the t_POL x divided by the t_POL y.

GEN gdivround(GEN x, GEN y) if x and y are t_INT, as diviiround. Operate componentwise if x is a t_COL, t_VEC or t_MAT. Otherwise as gdivent.

GEN centermod_i(GEN x, GEN y, GEN y2), as centermodii, componentwise.

GEN centermod(GEN x, GEN y), as centermod_i, except that y2 is computed (and left on the stack for efficiency).

GEN ginvmod(GEN x, GEN y) creates the inverse of x modulo y when it exists. y must be of type t_INT (in which case x is of type t_INT) or t_POL (in which case x is either a scalar type or a t_POL).

GCD, content and primitive part

GEN subres(GEN x, GEN y) creates the resultant of the t_POLs x and y computed using the subresultant algorithm. Either of x and y may also be scalars (treated as polynomials of degree 0)

GEN ggcd(GEN x, GEN y) creates the GCD of x and y.

GEN glcm(GEN x, GEN y) creates the LCM of x and y.

GEN gbezout(GEN x,GEN y, GEN *u,GEN *v) creates the GCD of x and y, and puts (the addresses of) objects u and v such that ux+vy = gcd(x,y) into *u and *v.

GEN bezoutpol(GEN a,GEN b, GEN *u,GEN *v), returns the GCD d of t_INTs a and b and sets u, v to the Bezout coefficients such that au + bv = d.

GEN content(GEN x) creates the GCD of all the components of x.

GEN primitive_part(GEN x, GEN *c), sets c to content(x) and returns the primitive part x / c.

GEN primpart(GEN x) as primitive_part but the content is lost. (For efficiency, the content remains on the stack.)

Generic binary operators.

Let ``op'' be a binary operation among

op = add: addition (x + y).

op = sub: subtraction (x - y).

op = mul: multiplication (x * y).

op = div: division (x / y).

op = max: maximum (max(x, y))

op = min: minimum (min(x, y))

The names and prototypes of the functions corresponding to op are as follows:

GEN gop[z](GEN x, GEN y[, GEN z])

GEN gopgs[z](GEN x, long s[, GEN z])

GEN gopsg[z](long s, GEN y[, GEN z])

GEN gpow(GEN x, GEN y, long l) creates x^{y}. If y is a t_INT, return powgi(x,y) (the precision l is not taken into account). Otherwise, the result is exp (y* log (x)) computed to precision l.

GEN gpowgs(GEN x, long n) creates x^{n} using binary powering.

GEN powgi(GEN x, GEN y) creates x^{y}, where y is a t_INT, using left-shift binary powering.

GEN gsubst(GEN x, long v, GEN y) substitutes the object y into x for the variable number v.

Miscellaneous functions

const char* type_name(long t) given a type number t this routine returns a string containing its symbolic name. E.g type_name(t_INT) returns "t_INT". The return value is read-only.


Further type specific functions

Vectors and Matrices

See Label se:clean and Label se:unclean for various useful constructors. Coefficients are accessed and set using gel, gcoeff, see Label se:accessors. There are many internal functions to extract or manipulate subvectors or submatrices but, like the accessors above, none of them are suitable for gerepileupto. Worse, there are no type verification, nor bound checking, so use at your own risk.

Note. In the function names below, i stands for interval and p for permutation.

GEN shallowcopy(GEN x) returns a t_GEN whose components are the components of x (no copy is made). The result may now be used to compute in place without destroying x. This is essentially equivalent to

    GEN y = cgetg(lg(x), typ(x));
    for (i = 1; i < lg(x); i++) y[i] = x[i];
    return y;

except that t_POLMOD (resp. t_MAT) are treated specially since a dummy copy of the representative (resp. all columns) is also made.

GEN shallowtrans(GEN x) returns the transpose of x, without copying its components, i. e., it returns a GEN whose components are (physically) the components of x. This is the internal function underlying gtrans.

GEN shallowconcat(GEN x, GEN y) concatenate x and y, without copying compoents, i. e., it returns a GEN whose components are (physically) the comonents of x and y.

GEN vconcat(GEN A, GEN B) concatenate vertically the two t_MAT A and B of compatible dimensions. A NULL pointer is accepted for an empty matrix. See shallowconcat.

GEN row(GEN A, long i) return A[i,], the i-th row of the t_MAT A.

GEN row_i(GEN A, long i, long j1, long j2) return part of the i-th row of t_MAT A: A[i,j_1], A[i,j_1+1]...,A[i,j_2]. Assume j_1 E<lt>= j_2.

GEN rowslice(GEN A, long i1, long i2) return the t_MAT formed by the i_1-th through i_2-th rows of t_MAT A. Assume i_1 E<lt>= i_2.

GEN rowpermute(GEN A, GEN p), p being a t_VECSMALL representing a list [p_1,...,p_n] of rows of t_MAT A, returns the matrix whose rows are A[p_1,],..., A[p_n,].

GEN rowslicepermute(GEN A, GEN p, long x1, long x2), short for

    rowslice(rowpermute(A,p), x1, x2)

(more efficient).

GEN vecslice(GEN A, long j1, long j2), return A[j_1],..., A[j_2]. If A is a t_MAT, these correspond to columns of A. The object returned has the same type as A (t_VEC, t_COL or t_MAT). Assume j_1 <= j_2.

GEN vecpermute(GEN A, GEN p) p is a t_VECSMALL representing a list [p_1,...,p_n] of indices. Returns a GEN which has the same type as A (t_VEC, t_COL or t_MAT), and whose components are A[p_1],...,A[p_n]. If A is a t_MAT, these are the columns of A.

GEN vecslicepermute(GEN A, GEN p, long y1, long y2) short for

    vecslice(vecpermute(A,p), y1, y2)

(more efficient).

Low-level vectors and columns functions

Theses functions handle t_VEC as an abstract container type of GENs. No specific meaning is attached to the content.

They accept both t_VEC and t_COL as input, but col functions always return t_COL and vec functions alway return t_VEC.

Note. All the functions below are shallow.

GEN const_col(long n, long c) returns a t_COL of n components equal to c.

GEN const_vec(long n, long c) returns a t_VEC of n components equal to c.

int vec_isconst(GEN v) Returns 1 if all the components of v are equal, else returns 0.

int vec_is1to1(GEN v) Returns 1 if the components of v are pair-wise distinct, i.e. if i|--->v[i] is a 1-to-1 mapping, else returns 0.

GEN vec_shorten(GEN v, long n) shortens the vector v to n components.

GEN vec_lengthen(GEN v, long n) lengthens the vector v to n components. The extra components are not initialised.

Function to handle t_VECSMALL

Theses functions handle t_VECSMALL as an abstract container type of small signed integers. No specific meaning is attached to the content.

GEN const_vecsmall(long n, long c) returns a t_VECSMALL of n components equal to c.

GEN vec_to_vecsmall(GEN z) identical to ZV_to_zv(z).

GEN vecsmall_to_vec(GEN z) identical to zv_to_ZV(z).

GEN vecsmall_to_col(GEN z) identical to zv_to_ZC(z).

GEN vecsmall_copy(GEN x) makes a copy of x on the stack.

GEN vecsmall_shorten(GEN v, long n) shortens the t_VECSMALL v to n components.

GEN vecsmall_lengthen(GEN v, long n) lengthens the t_VECSMALL v to n components. The extra components are not initialised.

GEN vecsmall_indexsort(GEN x) performs an indirect sort of the components of the t_VECSMALL x and return a permutation stored in a t_VECSMALL.

void vecsmall_sort(GEN v) sorts the t_VECSMALL v in place.

GEN vecsmall_uniq(GEN v) given a sorted t_VECSMALL v, return the vector of unique occurences.

int vecsmall_lexcmp(GEN x, GEN y) compares two t_VECSMALL lexically

int vecsmall_prefixcmp(GEN x, GEN y) truncate the longest t_VECSMALL to the length of the shortest and compares them lexicographically.

GEN vecsmall_prepend(GEN V, long s) prepend s to the t_VECSMALL V.

GEN vecsmall_append(GEN V, long s) append s to the t_VECSMALL V.

GEN vecsmall_concat(GEN u, GEN v) concat the t_VECSMALL u and v.

long vecsmall_coincidence(GEN u, GEN v) returns the numbers of indices where u and v agree.

long vecsmall_pack(GEN v, long base, long mod) handles the t_VECSMALL v as the digit of a number in base base and return this number modulo mod. This can be used as an hash function.

Functions to handle bits-vectors

Theses functions manipulate vectors of bits (stored in t_VECSMALL). Bits are numbered from 0.

GEN bitvec_alloc(long n) allocates a bits-vector suitable for n bits.

GEN bitvec_shorten(GEN bitvec, long n) shortens a bits-vector bitvec to n bits.

long bitvec_test(GEN bitvec, long b) returns the bit of index b of bitvec.

void bitvec_set(GEN bitvec, long b) (in place) sets the bit of index b of bitvec.

void bitvec_clear(GEN bitvec, long b) (in place) clears the bit of index b of bitvec.

Functions to handle vectors of t_VECSMALL

Theses functions manipulate vectors of t_VECSMALL (vecvecsmall).

GEN vecvecsmall_sort(GEN x) sorts lexicographically the components of the vector x.

GEN vecvecsmall_indexsort(GEN x) performs an indirect lexicographic sorting of the components of the vector x and return a permutation stored in a t_VECSMALL.

long vecvecsmall_search(GEN x, GEN y, long flag) x being a sorted vecvecsmall and y a t_VECSMALL, search y inside x. fla has the same meaning as for setsearch.