NAME

libPARI - The gp Calculator


Introduction

Originally, gp was designed as a debugging device for the PARI system library. Over the years, it has become a powerful user-friendly stand-alone calculator. The mathematical functions available in PARI and gp are described in the next chapter. In the present one, we describe the specific use of the gp programmable calculator.

 If you have GNU Emacs and use the PariEmacs package, you can work in a
special Emacs shell, described in L<Label se:emacs>. Specific features of
this Emacs shell are indicated by an EMACS sign in the left margin.

Startup

To start the calculator, the general command line syntax is:

gp [-D key = val] [files]

where items within brackets are optional. The [files] argument is a list of files written in the GP scripting language, which will be loaded on startup. There can be any number of arguments of the form -D key = val, setting some internal parameters of gp, or defaults: each sets the default key to the value val. See Label se:defaults below for a list and explanation of all defaults. These defaults can be changed by adding parameters to the input line as above, or interactively during a gp session, or in a preferences file also known as gprc.

If a preferences file (to be discussed in Label se:gprc) is found, gp then reads it and executes the commands it contains. This provides an easy way to customize gp. The files argument is processed right after the gprc.

A copyright banner then appears which includes the version number, and a lot of useful technical information. After the copyright, the computer writes the top-level help information, some initial defaults, and then waits after printing its prompt, which is '? ' by default . Whether extended on-line help and line editing are available or not is indicated in this gp banner, between the version number and the copyright message. Consider investigating the matter with the person who installed gp if they are not. Do this as well if there is no mention of the GMP kernel.

Getting help

To get help, type a ? and hit return. A menu appears, describing the main categories of available functions and how to get more detailed help. If you now type ?n with 1 <= n <= 11, you get the list of commands corresponding to category n and simultaneously to Section 3.n of this manual. If you type ?functionname where functionname is the name of a PARI function, you will get a short explanation of this function.

If extended help (see Label se:exthelp) is available on your system, you can double or triple the ? sign to get much more: respectively the complete description of the function (e.g. ??sqrt), or a list of gp functions relevant to your query (e.g.  ???"elliptic curve" or ???"quadratic field").

If gp was properly installed (see Appendix A), a line editor is available to correct the command line, get automatic completions, and so on. See Label se:readline or ??readline for a short summary of the line editor's commands.

If you type ?\ you will get a short description of the metacommands (keyboard shortcuts).

Finally, typing ?. will return the list of available (pre-defined) member functions. These are functions attached to specific kind of objects, used to retrieve easily some information from complicated structures (you can define your own but they won't be shown here). We will soon describe these commands in more detail.

More generally, commands starting with the symbols \b or ?, are not computing commands, but are metacommands which allow you to exchange information with gp. The available metacommands can be divided into default setting commands (explained below) and simple commands (or keyboard shortcuts, to be dealt with in Label se:meta).

Input

Just type in an instruction, e.g. 1 + 1, or Pi. No action is undertaken until you hit the < Return > key. Then computation starts, and a result is eventually printed. To suppress printing of the result, end the expression with a ; sign. Note that many systems use ; to indicate end of input. Not so in gp: a final semicolon means the result should not be printed. (Which is certainly useful if it occupies several screens.)

Interrupt, Quit

Typing quit at the prompt ends the session and exits gp. At any point you can type Ctrl-C (that is press simultaneously the Control and C keys): the current computation is interrupted and control given back to you at the gp prompt, together with a message like

    ***   at top-level: gcd(a,b)
    ***                 ^--------
    *** gcd: user interrupt after 236 ms.

telling you how much time elapsed since the last command was typed in and in which GP function the computation was aborted. It does not mean that that much time was spent in the function, only that the evaluator was busy processing that specific function when you stopped it.


The general gp input line

The gp calculator uses a purely interpreted language GP. The structure of this language is reminiscent of LISP with a functional notation, f(x,y) rather than (f x y): all programming constructs, such as if, while, etc...are functionsFOOTNOTE<<< Not exactly, since not all their arguments need be evaluated. For instance it would be stupid to evaluate both branches of an if statement: since only one will apply, only this one is evaluated. >>>, and the main loop does not really execute, but rather evaluates (sequences of) expressions. Of course, it is by no means a true LISP, and has been strongly influenced by C and Perl since then.

Introduction

User interaction with a gp session proceeds as follows. First, one types a sequence of characters at the gp prompt; see Label se:readline for a description of the line editor. When you hit the < Return > key, gp gets your input, evaluates it, then prints the result and assigns it to an ``history'' array.

More precisely, the input is case-sensitive and, outside of character strings, blanks are completely ignored. Inputs are either metacommands or sequences of expressions. Metacommands are shortcuts designed to alter gp's internal state, such as the working precision or general verbosity level; we shall describe them in Label se:meta, and ignore them for the time being.

The evaluation of a sequence of instructions proceeds in two phases: your input is first digested (byte-compiled) to a bytecode suitable for fast evaluation, in particular loop bodies are compiled only once but a priori evaluated many times; then the bytecode is evaluated.

An expression is formed by combining constants, variables, operator symbols, functions and control statements. It is evaluated using the conventions about operator priorities and left to right associativity. An expression always has a value, which can be any PARI object:

  ? 1 + 1
  %1 = 2          \\ an ordinary integer
  ? x
  %2 = x          \\ a polynomial of degree 1 in the unknown C<x>
  ? print("Hello")
  Hello           \\ C<void> return value
  ? f(x) = x^2
  %3 = (x)->x^2   \\ a user function

@3In the third example, Hello is printed as a side effect, but is not the return value. The print command is a procedure, which conceptually returns nothing. But in fact procedures return a special void object, meant to be ignored (but which evaluates to 0 in a numeric context, and stored as 0 in the history or results). The final example assigns to the variable f the function x|--->x^2, the alternative form f = x- > x^2 achieving the same effect; the return value of a function definition is, unsurprisingly, a function object (of type t_CLOSURE).

Several expressions are combined on a single line by separating them with semicolons (';'). Such an expression sequence will be called a seq. A seq also has a value, which is the value of the last expression in the sequence. Under gp, the value of the seq, and only this last value, becomes an history entry. The values of the other expressions in the seq are discarded after the execution of the seq is complete, except of course if they were assigned into variables. In addition, the value of the seq is printed if the line does not end with a semicolon ;.

The gp history of results

This is not to be confused with the history of your commands, maintained by readline. The gp history contains the results they produced, in sequence.

The successive elements of the history array are called %1, %2, ...As a shortcut, the latest computed expression can also be called %, the previous one %`, the one before that %`` and so on.

When you suppress the printing of the result with a semicolon, it is still stored in the history, but its history number will not appear either. It is a better idea to assign it to a variable for later use than to mentally recompute what its number is. Of course, on the next line, you may just use %.

The time used to compute that history entry is also stored as part of the entry and can be recovered using the %# operator: %#1, %#2, %#`; %# by itself returns the time needed to compute the last result (the one returned by %).

@3Remark. The history ``array'' is in fact better thought of as a queue: its size is limited to 5000 entries by default, after which gp starts forgetting the initial entries. So %1 becomes unavailable as gp prints %5001. You can modify the history size using histsize.

Special editing characters

A GP program can of course have more than one line. Since your commands are executed as soon as you have finished typing them, there must be a way to tell gp to wait for the next line or lines of input before doing anything. There are three ways of doing this.

The first one is to use the backslash character \ at the end of the line that you are typing, just before hitting < Return > . This tells gp that what you will write on the next line is the physical continuation of what you have just written. In other words, it makes gp forget your newline character. You can type a \ anywhere. It is interpreted as above only if (apart from ignored whitespace characters) it is immediately followed by a newline. For example, you can type

  ? 3 + \
  4

@3instead of typing 3 + 4.

The second one is a variation on the first, and is mostly useful when defining a user function (see Label se:user_defined): since an equal sign can never end a valid expression, gp disregards a newline immediately following an = .

  ? a =
  123
  %1 = 123

The third one is in general much more useful, and uses braces { and }. An opening brace { signals that you are typing a multi-line command, and newlines are ignored until you type a closing brace }. There are two important, but easily obeyed, restrictions: first, braces do not nest; second, inside an open brace-close brace pair, all input lines are concatenated, suppressing any newlines. Thus, all newlines should occur after a semicolon (;), a comma (,) or an operator (for clarity's sake, never split an identifier over two lines in this way). For instance, the following program

  {
    a = b
    b = c
  }

@3would silently produce garbage, since this is interpreted as a = bb = c which assigns the value of c to both bb and a. It should have been written

  {
    a = b;
    b = c;
  }


The PARI types

We see here how to input values of the different data types known to PARI. Recall that blanks are ignored in any expression which is not a string (see below).

@3A note on efficiency. The following types are provided for convenience, not for speed: t_INTMOD, t_FRAC, t_PADIC, t_QUAD, t_POLMOD, t_RFRAC. Indeed, they always perform a reduction of some kind after each basic operation, even though it is usually more efficient to perform a single reduction at the end of some complex computation. For instance, in a convolution product sum_{i+j = n} x_i y_j in Z/NZ --- common when multiplying polynomials! ---, it is quite wasteful to perform n reductions modulo N. In short, basic individual operations on these types are fast, but recursive objects with such components could be handled more efficiently: programming with libpari will save large constant factors here, compared to GP.

Integers (t_INT)

After an (optional) leading + or -, type in the decimal digits of your integer. No decimal point!

  ? 1234567
  %1 = 1234567
  ? -3
  %2 = -3
  ? 1.         \\ oops, not an integer
  %3 = 1.000000000000000000000000000

Real numbers (t_REAL)

Real numbers are represented (approximately) in a floating point system, internally in base 2, but converted to base 10 for input / output purposes. A t_REAL object has a given accuracy (or precision) \ell >= 0; it comprises

@3* a sign s: +1, -1 or 0;

@3* a mantissa m: a multiprecision integer, 0 <= m < 10^\ell;

@3* an exponent e: a small integer in [-E,E], where E ~ 2^B log _{10} 2, and B = 32 on a 32-bit machine and 64 otherwise.

This data may represent any real number x such that

  |x - s m 10^e| < 10^{e-\ell}.

We consider that a t_REAL with sign s = 0 has accuracy \ell = 0, so that its mantissa is useless, but it still has an exponent e and acts like a machine epsilon for all accuracies < e.

After an (optional) leading + or -, type a number with a decimal point. Leading zeroes may be omitted, up to the decimal point, but trailing zeroes are important: your t_REAL is assigned an internal precision, which is the supremum of the input precision, one more than the number of decimal digits input, and the default precision. For example, if the default precision is 28 digits, typing 2. yields a precision of 28 digits, but 2.0...0 with 45 zeros gives a number with internal precision at least 45, although less may be printed.

You can also use scientific notation with the letter E or e. As usual, en is interpreted as x 10^n for all integers n. Since the result is converted to a t_REAL, you may often omit the decimal point in this case: 6.02 E 23 or 1e-5 are fine, but e10 is not.

By definition, 0.E n returns a real 0 of exponent n, whereas 0. returns a real 0 ``of default precision'' (of exponent -realprecision), see Label se:whatzero, behaving like the machine epsilon for the current default accuracy: any float of smaller absolute value is indistinguishable from 0.

@3Note on output formats. A zero real number is printed in e format as 0.Exx where xx is the (usually negative) decimal exponent of the number (cf. Label se:whatzero). This allows the user to check the accuracy of that particular zero.

When the integer part of a real number x is not known exactly because the exponent of x is greater than the internal precision, the real number is printed in e format.

Intmods (t_INTMOD)

To create the image of the integer a in Z/bZ (for some non-zero integer b), type Mod(a,b); not a%b. Internally, all operations are done on integer representatives belonging to [0,b-1].

Note that this type is available for convenience, not for speed: each elementary operation involves a reduction modulo b.

If x is a t_INTMOD Mod(a,b), the following member function is defined:

x.mod: return the modulus b.

Rational numbers (t_FRAC)

All fractions are automatically reduced to lowest terms, so it is impossible to work with reducible fractions. To enter n/m just type it as written. As explained in Label se:gdiv, floating point division is not performed, only reduction to lowest terms.

Note that rational computation are almost never the fastest method to proceed: in the PARI implementation, each elementary operation involves computing a gcd. It is generally a little more efficient to cancel denominators and work with integers only:

  ? P = Pol( vector(10^3,i, 1/i) ); \\ big polynomial with small rational coeffs
  ? P^2
  time = 1,392 ms.
  ? c = content(P); c^2 * (P/c)^2;  \\ same computation in integers
  time = 1,116 ms.

And much more efficient (but harder to setup) to use homomorphic imaging schemes and modular computations. As the simple example below indicates, if you only need modular information, it is very worthwhile to work with t_INTMODs directly, rather than deal with t_FRACs all the way through:

  ? p = nextprime(10^7);
  ? sum(i=1, 10^5, 1/i) % p
  time = 13,288 ms.
  %1 = 2759492
  ? sum(i=1, 10^5, Mod(1/i, p))
  time = 60 ms.
  %2 = Mod(2759492, 10000019)

Finite field elements (t_FFELT)

Let T belongs to F_p[X] be a monic irreducible polynomial defining your finite field over F_p, for instance obtained using ffinit. Then the ffgen function creates a generator of the finite field as an F_p-algebra, namely the class of X in F_p[X]/(T), from which you can build all other elements. For instance, to create the field F_{2^8}, we write

  ? T = ffinit(2, 8);
  ? y = ffgen(T, 'y);
  ? y^0    \\ the unit element in the field
  %3 = 1
  ? y^8
  %4 = y^6 + y^5 + y^4 + y^3 + y + 1

The second (optional) parameter to ffgen is only used to display the result; it is customary to use the name of the variable we assign the generator to. If g is a t_FFELT, the following member functions are defined:

g.pol: the polynomial (with reduced integer coefficients) expressing g in term of the field generator.

g.p: the characteristic of the finite field.

g.f: the dimension of the definition field over its prime field; the cardinality of the definition field is thus p^f.

g.mod: the minimal polynomial (with reduced integer coefficients) of the field generator.

Complex numbers (t_COMPLEX)

To enter x+iy, type x + I*y. (That's I, not i!) The letter I stands for sqrt {-1}. The ``real'' and ``imaginary'' parts x and y can be of type t_INT, t_REAL, t_INTMOD, t_FRAC, or t_PADIC.

p-adic numbers (t_PADIC):

Typing O(p^k), where p and k are integers, yields a p-adic 0 of accuracy k, representing any p-adic number whose valuation is >= k. To input a general non-0 p-adic number, write a suitably precise rational or integer approximation and add O(p^ k) to it.

Note that it is not checked whether p is indeed prime but results are undefined if this is not the case: you can work on 10-adics if you want, but disasters will happen as soon as you do something non-trivial like taking a square root. Note that O(25) is not the same as O(5^2); you want the latter!

For example, you can type in the 7-adic number

2*7^(-1) + 3 + 4*7 + 2*7^2 + O(7^3)

@3exactly as shown, or equivalently as 905/7 + O(7^3).

If a is a t_PADIC, the following member functions are defined:

a.mod: returns the modulus p^k.

a.p: returns p.

Note that this type is available for convenience, not for speed: internally, t_PADICs are stored as p-adic units modulo some p^k. Each elementary operation involves updating p^k (multiplying or dividing by powers of p) and a reduction mod p^k. In particular, additions are slow.

      ? n = 1+O(2^20);   for (i=1,10^6, n++)
      time = 841 ms.
      ? n = Mod(1,2^20); for (i=1,10^6, n++)
      time = 441 ms.
      ? n = 1;           for (i=1,10^6, n++)
      time = 328 ms.

@3The penalty associated with maintaining p^k decreases steeply as p increases (and updates become very rare). But t_INTMODs remain at least 25% more efficient. (But they do not have denominators!)

Quadratic numbers (t_QUAD)

This type is used to work in the quadratic order of discriminant d, where d is a non-square integer congruent to 0 or 1 (modulo 4). The command

      w = quadgen(d)

assigns to w the ``canonical'' generator for the integer basis of the order of discriminant d, i.e. w = sqrt {d}/2 if d = 0 mod 4, and w = (1+ sqrt {d})/2 if d = 1 mod 4. The name w is of course just a suggestion, but corresponds to traditional usage. You can use any variable name that you like, but quadgen(d) is always printed as w, regardless of the discriminant. So beware, two t_QUADs can be printed in the same way and not be equal; however, gp will refuse to add or multiply them for example.

Since the order is Z + wZ, any other element can be input as x+y*w for some integers x and y. In fact, you may work in its fraction field Q( sqrt {d}) and use t_FRAC values for x and y.

Polmods (t_POLMOD)

Exactly as for intmods, to enter x mod y (where x and y are polynomials), type Mod(x,y), not x%y. Note that when y is an irreducible polynomial in one variable, polmods whose modulus is y are simply algebraic numbers in the finite extension defined by the polynomial y. This allows us to work easily in number fields, finite extensions of the p-adic field Q_p, or finite fields.

Note that this type is available for convenience, not for speed: each elementary operation involves a reduction modulo y. If p is a t_POLMOD, the following member functions are defined:

p.pol: return a representative of the polynomial class of minimal degree.

p.mod: return the modulus.

\misctitle{Important remark} Mathematically, the variables occurring in a polmod are not free variables. But internally, a congruence class in R[t]/(y) is represented by its representative of lowest degree, which is a t_POL in R[t], and computations occur with polynomials in the variable t. PARI will not recognize that Mod(y, y^2 + 1) is ``the same'' as Mod(x, x^2 + 1), since x and y are different variables.

To avoid inconsistencies, polmods must use the same variable in internal operations (i.e. between polmods) and variables of lower priority for external operations, typically between a polynomial and a polmod. See Label se:priority for a definition of ``priority'' and a discussion of (PARI's idea of) multivariate polynomial arithmetic. For instance:

      ? Mod(x, x^2+ 1) + Mod(x, x^2 + 1)
      %1 = Mod(2*x, x^2 + 1)    \\ 2i (or -2i), with i^2 = -1
      ? x + Mod(y, y^2 + 1)
      %2 = x + Mod(y, y^2 + 1)  \\ in B<Q>(i)[x]
      ? y + Mod(x, x^2 + 1)
      %3 = Mod(x + y, x^2 + 1)  \\ in B<Q>(y)[i]

The first two are straightforward, but the last one may not be what you want: y is treated here as a numerical parameter, not as a polynomial variable.

If the main variables are the same, it is allowed to mix t_POL and t_POLMODs. The result is the expected t_POLMOD. For instance

      ? x + Mod(x, x^2 + 1)
      %1 = Mod(2*x, x^2 + 1)

Polynomials (t_POL)

Type the polynomial in a natural way, not forgetting to put a ``*'' between a coefficient and a formal variable;

  ? 1 + 2*x + 3*x^2
  %1 = 3*x^2 + 2*x + 1

This assumes that x is still a ''free variable''.

  ? x = 1; 1 + 2*x + 3*x^2
  %2 = 6

generates an integer, not a polynomial! It is good practice to never assign values to polynomial variables to avoid the above problem, but a foolproof construction is available using 'x instead of x: 'x is a constant evaluating to the free variable with name x, independently of the current value of x.

  ? x = 1; 1 + 2*'x + 3*'x^2
  %3 = 1 + 2*x + 3*x^2
  ? x = 'x; 1 + 2*x + 3*x^2
  %4 = 1 + 2*x + 3*x^2

You may also use the functions Pol or Polrev:

  ? Pol([1,2,3])       \\ C<Pol> creates a polynomial in C<x> by default
  %1 = x^2 + 2*x + 3
  ? Polrev([1,2,3])
  %2 = 3*x^2 + 2*x + 1
  ? Pol([1,2,3], 'y)   \\ we use C<'y>, safer than C<y>
  %3 = y^2 + 2*y + 3

@3The latter two are much more efficient constructors than an explicit summation (the latter is quadratic in the degree, the former linear):

  ? for (i=1, 10^4, Polrev( vector(100, i,i) ) )
  time = 124ms
  ? for (i=1, 10^4, sum(i = 1, 100, (i+1) * 'x^i) )
  time = 3,985ms

Polynomials are always printed as univariate polynomials, with monomials sorted by decreasing degree:

  ? (x+y+1)^2
  %1 = x^2 + (2*y + 2)*x + (y^2 + 2*y + 1)

(Univariate polynomial in x whose coefficients are polynomials in y.) See Label se:varsymb for valid variable names, and a discussion of multivariate polynomial rings.

Power series (t_SER)

Typing O(X^k), where k is an integer, yields an X-adic 0 of accuracy k, representing any power series in X whose valuation is >= k. Of course, X can be replaced by any other variable name! To input a general non-0 power series, type in a polynomial or rational function (in X, say), and add O(X^k) to it. The discussion in the t_POL section about variables remains valid; a constructor Ser replaces Pol and Polrev.

@3Caveat. Power series with inexact coefficients sometimes have a non-intuitive behavior: if k significant terms are requested, an inexact zero is counted as significant, even if it is the coefficient of lowest degree. This means that useful higher order terms may be disregarded.

If a series with a zero leading coefficient must be inverted, then as a desperation measure that coefficient is discarded, and a warning is issued:

  ? C = 0. + y + O(y^2);
  ? 1/C
    *** _/_: Warning: normalizing a series with 0 leading term.
  %2 = y^-1 + O(1)

The last output could be construed as a bug since it is a priori impossible to deduce such a result from the input (0. represents any sufficiently small real number). But it was thought more useful to try and go on with an approximate computation than to raise an early exception.

If the series precision is insufficient, errors may occur (mostly division by 0), which could have been avoided by a better global understanding of the computation:

  ? A = 1/(y + 0.); B = 1. + O(y);
  ? B * denominator(A)
  %2 = 0.E-28 + O(y)
  ? A/B
    *** _/_: Warning: normalizing a series with 0 leading term.
  %3 = 1.000000000000000000000000000*y^-1 + O(1)
  ? A*B
    *** _*_: Warning: normalizing a series with 0 leading term.
  %4 = 1.000000000000000000000000000*y^-1 + O(1)

Rational functions (t_RFRAC)

As for fractions, all rational functions are automatically reduced to lowest terms. All that was said about fractions in Label se:FRAC remains valid here.

Binary quadratic forms of positive or negative discriminant (t_QFR and t_QFI)

These are input using the function Qfb. For example Qfb(1,2,3) creates the binary form x^2+2xy+3y^2. It is imaginary (of internal type t_QFI) since its discriminant 2^2 - 4 x 3 = -8 is negative. Although imaginary forms could be positive or negative definite, only positive definite forms are implemented.

In the case of forms with positive discriminant (t_QFR), you may add an optional fourth component (related to the regulator, more precisely to Shanks and Lenstra's ``distance''), which must be a real number. See also the function qfbprimeform which directly creates a prime form of given discriminant.

Row and column vectors (t_VEC and t_COL)

) To enter a row vector, type the components separated by commas ``,'', and enclosed between brackets ``[

  '' and ``

]'', e.g. [1,2,3]. To enter a column vector, type the vector horizontally, and add a tilde ``~'' to transpose. [ ] yields the empty (row) vector. The function Vec can be used to transform any object into a vector (see Chapter 3). The construction [i..j], where i <= j are two integers returns the vector [i, i+1,..., j-1, j]

  ? [1,2,3]
  %1 = [1, 2, 3]
  ? [-2..3]
  %2 = [-2, -1, 0, 1, 2, 3]

Let the variable v contain a (row or column) vector:

@3* v[m] refers to its m-th entry; you can assign any value to v[m], i.e. write something like v[m] = expr.

@3* v[i..j], where i <= j, returns the vector slice containing elements v[i],..., v[j]; you can not assign a result to v[i..j].

@3* v[^i] returns the vector whose i-th entry has been removed; you can not assign a result to v[^i].

@3In the last two constructions v[i..j] and v[^i], i and j are allowed to be negative integers, in which case, we start counting from the end of the vector: e.g., -1 is the index of the last element.

  ? v = [1,2,3,4];
  ? v[2..4]
  %2 = [2, 3, 4]
  ? v[^3]
  %3 = [1, 2, 4]
  ? v[^-1]
  %3 = [1, 2, 3]
  ? v[-3..-1]
  %4 = [2, 3, 4]

@3Remark. vector is the standard constructor for row vectors whose i-th entry is given by a simple function of i; vectorv is similar for column vectors:

   ? vector(10, i, i^2+1)
   %1 = [2, 5, 10, 17, 26, 37, 50, 65, 82, 101]

The functions Vec and Col convert objects to row and column vectors respectively (as well as Vecrev and Colrev, which revert the indexing):

   ? T = poltchebi(5)   \\ 5-th Chebyshev polynomial
   %1 = 16*x^5 - 20*x^3 + 5*x
   ? Vec(T)
   %2 = [16, 0, -20, 0, 5, 0]  \\ coefficients of T
   ? Vecrev(T)
   %3 = [0, 5, 0, -20, 0, 16]  \\ ... in reverse order

@3Remark. For v a t_VEC, t_COL, t_LIST or t_MAT, the alternative set-notations

   [g(x) | x <- v, f(x)]
   [x | x <- v, f(x)]
   [g(x) | x <- v]

are available as shortcuts for

   apply(g, select(f, Vec(v)))
   select(f, Vec(v))
   apply(g, Vec(v))

@3respectively, and may serve as t_VEC constructors:

   ? [ p | p <- primes(10), isprime(p+2) ]
   %2 = [3, 5, 11, 17, 29]

@3returns the primes p (among the first 10 primes) such that (p, p+2) is a twin pair;

   ? [ p^2 | p <- primes(10), p % 4 == 1 ]
   %1 = [25, 169, 289, 841]

@3returns the squares of the primes congruent to 1 modulo 4, where p runs among the first 10 primes.

Matrices (t_MAT)

To enter a matrix, type the components row by row, the components being separated by commas ``,'', the rows by semicolons ``;'', and everything enclosed in brackets ``[

  '' and ``

]'', e.g. [x,y; z,t; u,v]. [;] yields an empty (0 x 0) matrix. The function Mat transforms any object into a matrix, and matrix creates matrices whose (i,j)-th entry is described by a function f(i,j):

  ? Mat(1)
  %1 =
  [1]
  ? matrix(2,2, i,j, 2*i+j)
  %2 =
  [3 4]
  [5 6]

@3Let the variable M contain a matrix, and let i,j,k,l denote four integers:

@3* M[i,j] refers to its (i,j)-th entry; you can assign any result to M[i,j].

@3* M[i,] refers to its i-th row; you can assign a t_VEC of the right dimension to M[i,].

@3* M[,j] refers to its j-th column; you can assign a t_COL of the right dimension to M[,j].

@3But M[i] is meaningless and triggers an error. The ``range'' i..j and ``caret'' ^c notations are available as for vectors; you can not assign to any of these:

@3* M[i..j, k..l], i <= j, k <= l, returns the submatrix built from the rows i to j and columns k to l of M. assign to M[i..j, j..l].

@3* M[i..j,] returns the submatrix built from the rows i to j of M.

@3* M[,i..j] returns the submatrix built from the columns i to j of M.

@3* M[i..j,^k], i <= j, returns the submatrix built from the rows i to j and column k removed.

@3* M[^k,] returns the submatrix with row k removed.

@3* M[,^k] returns the submatrix with column k removed.

@3Finally,

@3* M[i..j, k] returns the t_COL built from the k-th column (entries i to j).

@3* M[^i, k] returns the t_COL built from the k-th column (entry i removed).

@3* M[k, i..j] returns the t_VEC built from the k-th row (entries i to j).

@3* M[k,^i] returns the t_VEC built from the k-th row (entry i removed).

  ? M = [1,2,3;4,5,6;7,8,9];
  ? M[1..2, 2..3]
  %2 =
  [2 3]
  [5 6]
  ? M[1..2,]
  %3 =
  [1 2 3]
  [4 5 6]
  ? M[,2..3]
  %4 =
  [2 3]
  [5 6]
  [8 9]

All this is recursive, so if M is a matrix of matrices of..., an expression such as M[1,1][,3][4] = 1 is perfectly valid (and actually identical to M[1,1][4,3] = 1), assuming that all matrices along the way have compatible dimensions.

@3Technical note (design flaw). Matrices are internally represented as a vector of columns. All matrices with 0 columns are thus represented by the same object (internally, an empty vector), and there is no way to distinguish between them. Thus it is not possible to create or represent matrices with zero columns and an actual nonzero number of rows. The empty matrix [;] is handled as though it had an arbitrary number of rows, exactly as many as needed for the current computation to make sense:

  ? [1,2,3; 4,5,6] * [;]
  %1 = [;]

The empty matrix on the first line is understood as a 3 x 0 matrix, and the result as a 2 x 0 matrix. On the other hand, it is possible to create matrices with a given positive number of columns, each of which has zero rows, e.g. using Mat as above or using the matrix function.

Note that although the internal representation is essentially the same, a row vector of column vectors is not a matrix; for example, multiplication will not work in the same way. It is easy to go from one representation to the other using Vec / Mat, though:

  ? [1,2,3;4,5,6]
  %1 =
  [1 2 3]
  [4 5 6]
  ? Vec(%)
  %2 = [[1, 4]~, [2, 5]~, [3, 6]~]
  ? Mat(%)
  %3 =
  [1 2 3]
  [4 5 6]

Lists (t_LIST)

Lists can be input directly, as in List([1,2,3,4]); but in most cases, one creates an empty list, then appends elements using listput:

    ? a = List(); listput(a,1); listput(a,2);
    ? a
    %2 = List([1, 2])

Elements can be accessed directly as with the vector types described above.

Strings (t_STR)

To enter a string, enclose it between double quotes ", like this: "this is a string". The function Str can be used to transform any object into a string.

Small vectors (t_VECSMALL)

This is an internal type, used to code in an efficient way vectors containing only small integers, such as permutations. Most gp functions will refuse to operate on these objects.

Functions (t_CLOSURE)

We will explain this at length in Label se:user_defined. For the time being, suffice it to say that functions can be assigned to variables, as any other object, and the following equivalent basic forms are available to create new ones

    f = (x,y) -> x^2 + y^2
    f(x,y) = x^2 + y^2

Error contexts (t_ERROR)

An object of this type is created whenever an error occurs: it contains some information about the error and the error context. Usually, an appropriate error is printed immediately, the computation is aborted, and GP enters the ``break loop'':

    ? 1/0; 1 + 1
      ***   at top-level: 1/0;1+1
      ***                  ^------
      *** _/_: division by a non-invertible object
      ***   Break loop: type 'break' to go back to the GP prompt

@3Here the computation is aborted as soon as we try to evaluate 1/0, and 1 + 1 is never executed. Exceptions can be trapped using iferr, however: we can evaluate some expression and either recover an ordinary result (no error occurred), or an exception (an error did occur).

    ? i = Mod(6,12); iferr(1/i, E, print(E)); 1 + 1
    error("impossible inverse modulo: Mod(6, 12).")
    %1 = 2

@3One can ignore the exception, print it as above, or extract non trivial information from the error context:

    ? i = Mod(6,12); iferr(1/i, E, print(component(E,1)));
    Mod(6, 12)

@3We can also rethrow the exception: error(E).


GP operators

@3Loosely speaking, an operator is a function, usually associated to basic arithmetic operations, whose name contains only non-alphanumeric characters. For instance + or -, but also = or += , or even [ ] (the selection operator). As all functions, operators take arguments, and return a value; assignment operators also have side effects: besides returning a value, they change the value of some variable.

Each operator has a fixed and unchangeable priority, which means that, in a given expression, the operations with the highest priority is performed first. Unless mentioned otherwise, operators at the same priority level are left-associative (performed from left to right), unless they are assignments, in which case they are right-associative. Anything enclosed between parenthesis is considered a complete subexpression, and is resolved recursively, independently of the surrounding context. For instance,

    a + b + c    -->   (a + b) + c     \\ left-associative
    a = b = c    -->   a = (b = c)     \\ right-associative

Assuming that op_1, op_2, op_3 are binary operators with increasing priorities (think of +, *, ^),

   x op_1 y op_2 z op_2 x op_3 y

is equivalent to

   x op_1 ((y op_2 z) op_2 (x op_3 y)).

GP contains many different operators, either unary (having only one argument) or binary, plus a few special selection operators. Unary operators are defined as either prefix or postfix, meaning that they respectively precede (op x) and follow (x op) their single argument. Some symbols are syntactically correct in both positions, like !, but then represent different operators: the ! symbol represents the negation and factorial operators when in prefix and postfix position respectively. Binary operators all use the (infix) syntax x op y.

Most operators are standard (+, %, = ), some are borrowed from the C language (++, << ), and a few are specific to GP (\, #). Beware that some GP operators differ slightly from their C counterparts. For instance, GP's postfix ++ returns the new value, like the prefix ++ of C, and the binary shifts << , >> have a priority which is different from (higher than) that of their C counterparts. When in doubt, just surround everything by parentheses; besides, your code will be more legible.

@3Here is the list of available operators, ordered by decreasing priority, binary and left-associative unless mentioned otherwise. An expression is an lvalue if something can be assigned to it. (The name comes from left-value, to the left of a = operator; e.g. x, or v[1] are lvalues, but x + 1 is not.)

\point{Priority 14}

: as in x:small, is used to indicate to the GP2C compiler that the variable on the left-hand side always contains objects of the type specified on the right hand-side (here, a small integer) in order to produce more efficient or more readable C code. This is ignored by GP.

\point{Priority 13} ( ) is the function call operator. If f is a closure and args is a comma-separated list of arguments (possibly empty), f(args) evaluates f on those arguments.

@3* Priority 12

++ and -- (unary, postfix): if x is an lvalue, x++ assigns the value x+1 to x, then returns the new value of x. This corresponds to the C statement ++x: there is no prefix ++ operator in GP. x-- does the same with x-1. These operators are not associative, i.e. x++++ is invalid, since x++ is not an lvalue.

@3* Priority 11

.member (unary, postfix): x.member extracts member from structure x (see Label se:member).

[ ] is the selection operator. x[i] returns the i-th component of vector x; x[i,j], x[,j] and x[i,] respectively return the entry of coordinates (i,j), the j-th column, and the i-th row of matrix x. If the assignment operator ( = ) immediately follows a sequence of selections, it assigns its right hand side to the selected component. E.g x[1][1] = 0 is valid; but beware that (x[1])[1] = 0 is not (because the parentheses force the complete evaluation of x[1], and the result is not modifiable).

@3* Priority 10

' (unary, postfix): derivative with respect to the main variable. If f is a function (t_CLOSURE), f' is allowed and defines a new function, which will perform numerical derivation when evaluated at a scalar x; this is defined as (f(x+varepsilon) - f(x-varepsilon)) / 2varepsilon for a suitably small epsilon depending on current precision.

  ? (x^2 + y*x + y^2)'  \\ derive with respect to main variable C<x>
  %1 = 2*x + y
  ? SIN = cos'
  %2 = cos'
  ? SIN(Pi/6)         \\ numerical derivation
  %3 = -0.5000000000000000000000000000
  ? cos'(Pi/6)        \\ works directly: no need for intermediate C<SIN>
  %4 = -0.5000000000000000000000000000

~ (unary, postfix): vector/matrix transpose.

! (unary, postfix): factorial. x! = x(x-1)...1.

! (unary, prefix): logical not. !x returns 1 if x is equal to 0 (specifically, if gequal0(x) == 1), and 0 otherwise.

@3* Priority 9

# (unary, prefix): cardinality; #x returns length(x).

@3* Priority 8

^: powering. This operator is right associative: 2^3^4 is understood as 2^(3^4).

@3* Priority 7

+, - (unary, prefix): - toggles the sign of its argument, + has no effect whatsoever.

@3* Priority 6

*: multiplication.

/: exact division (3/2 yields 3/2, not 1.5).

\, %: Euclidean quotient and remainder, i.e. if x = qy + r, then x \ y = q, x%y = r. If x and y are scalars, then q is an integer and r satisfies 0 <= r < y; if x and y are polynomials, then q and r are polynomials such that deg r < deg y and the leading terms of r and x have the same sign.

\/: rounded Euclidean quotient for integers (rounded towards + oo when the exact quotient would be a half-integer).

<< , >> : left and right binary shift. By definition, x << n = x * 2^n if n > 0, and truncate(x 2^{-n}) otherwise. Right shift is defined by x >> n = x << (-n).

@3* Priority 5

+, -: addition/subtraction.

@3* Priority 4

< , > , <= , >= : the usual comparison operators, returning 1 for true and 0 for false. For instance, x <= 1 returns 1 if x <= 1 and 0 otherwise.

< > , != : test for (exact) inequality.

== : test for (exact) equality. t_QFR having the same coefficients but a different distance component are tested as equal.

== = : test whether two objects are identical component-wise. This is stricter than == : for instance, the integer 0, a 0 polynomial or a vector with 0 entries, are all tested equal by == , but they are not identical.

@3* Priority 3

&&: logical and.

||: logical (inclusive) or. Any sequence of logical or and and operations is evaluated from left to right, and aborted as soon as the final truth value is known. Thus, for instance,

    x == 0 || test(1/x)

will never produce an error since test(1/x) is not even evaluated when the first test is true (hence the final truth value is true). Similarly

    type(p) == "t_INT" && isprime(p)

does not evaluate isprime(p) if p is not an integer.

@3* Priority 2

= (assignment, lvalue = expr). The result of x = y is the value of the expression y, which is also assigned to the variable x. This assignment operator is right-associative. This is not the equality test operator; a statement like x = 1 is always true (i.e. non-zero), and sets x to 1; the equality test would be x == 1. The right hand side of the assignment operator is evaluated before the left hand side.

It is crucial that the left hand-side be an lvalue there, it avoids ambiguities in expressions like 1 + x = 1. The latter evaluates as 1 + (x = 1), not as (1 + x) = 1, even though the priority of = is lower than the priority of +: 1 + x is not an lvalue.

If the expression cannot be parsed in a way where the left hand side is an lvalue, raise an error.

  ? x + 1 = 1
    ***   unused characters: x+1=1
    ***                         ^--

op = , where op is any binary operator among +, -, *, %, /, \, \/, << , or >> (composed assignment lvalue op = expr). The expression x op = y assigns (x op y) to x, and returns the new value of x. The result is not an lvalue; thus

    (x += 2) = 3

is invalid. These assignment operators are right-associative:

    ? x = 'x; x += x *= 2
    %1 = 3*x

@3* Priority 1 - > (function definition): (vars)- > expr returns a function object, of type t_CLOSURE.

@3Remark. Use the op = operators as often as possible since they make complex assignments more legible: one needs not parse complicated expressions twice to make sure they are indeed identical. Compare

  v[i+j-1] = v[i+j-1] + 1    -->    v[i+j-1]++
  M[i,i+j] = M[i,i+j] * 2    -->    M[i,i+j] *= 2

@3Remark. Less important but still interesting. The ++, -- and op = operators are slightly more efficient:

  ? a = 10^6;
  ? i = 0; while(i<a, i=i+1)
  time = 365 ms.
  ? i = 0; while(i<a, i++)
  time = 352ms.

@3For the same reason, the shift operators should be preferred to multiplication:

  ? a = 1<<(10^5);
  ? i = 1; while(i<a, i=i*2);
  time = 1,052 ms.
  ? i = 1; while(i<a, i<<=1);
  time = 617 ms.


Variables and symbolic expressions

In this section we use variable in the standard mathematical sense, symbols representing algebraically independent elements used to build rings of polynomials and power series, and explain the all-important concept of variable priority. In the next Label se:scope, we shall no longer consider only free variables, but adopt the viewpoint of computer programming and assign values to these symbols: (bound) variables are names associated to values in a given scope.

Variable names

A valid name starts with a letter, followed by any number of keyword characters: _ or alphanumeric characters ([A-Za-z0-9]). The built-in function names are reserved and cannot be used; see the list with \c, including the constants Pi, Euler, Catalan and I = sqrt {-1}.

GP names are case sensitive. For instance, the symbol i is perfectly safe to use, and will not be mistaken for I = sqrt {-1}; analogously, o is not synonymous to O.

In GP you can use up to 16383 variable names (up to 65535 on 64-bit machines). If you ever need thousands of variables and this becomes a serious limitation, you should probably be using vectors instead: e.g. instead of variables X1, X2, X3,..., you might equally well store their values in X[1], X[2], X[3],...

Variables and polynomials

What happens when you use a valid variable name,t say, for the first time before assigning a value into it? This registers a new free variable with the interpreter (which will be written as t), and evaluates to a monomial of degree 1 in the said variable t. It is important to understand that PARI/GP is not a symbolic manipulation package: even free variables already have default values FOOTNOTE<<< More generally, any expression has a value, and is replaced by its value as soon as it is read; it never stays in an abstract form. >>> , there is no such thing as an ``unbound'' variable in GP. You have access to this default value using the quote operator: 't always evaluates to the above monomial of degree 1, independently of assignments made since then (e.g. t = 1).

  ? t^2 + 1
  %1 = t^2 + 1
  ? t = 2; t^2 + 1
  %2 = 5
  ? %1
  %3 = t^2 + 1
  ? eval(%1)
  %4 = 5

In the above, t is initially a free variable, later bound to 2. We see that assigning a value to a variable does not affect previous expressions involving it; to take into account the new variable's value, one must force a new evaluation, using the function eval (see Label se:eval). It is preferable to leave alone your ``polynomial variables'', never assigning values to them, and to use subst and its more powerful variants rather than eval. You will avoid the following kind of problems:

  ? p = t^2 + 1; subst(p, t, 2)
  %1 = 5
  ? t = 2;
  ? subst(p, t, 3)    \\ C<t> is no longer free: it evaluates to 2
    ***   at top-level: subst(p,t,3)
    ***                         ^----
    ***   variable name expected.
  ? subst(p, 't, 3)   \\ OK
  %3 = 10

A statement like x = 'x in effect restores x as a free variable.

Variable priorities, multivariate objects

A multivariate polynomial in PARI is just a polynomial (in one variable), whose coefficients are themselves polynomials, arbitrary but for the fact that they do not involve the main variable. (PARI currently has no sparse representation for polynomials, listing only non-zero monomials.) All computations are then done formally on the coefficients as if the polynomial was univariate.

This is not symmetrical. So if I enter x + y in a clean session, what happens? This is understood as

   x^1 + (y^1 + 0*y^0)*x^0 belongs to (Z[y])[x]

but how do we know that x is ``more important'' than y ? Why not y^1 + x*y^0, which is the same mathematical entity after all?

The answer is that variables are ordered implicitly by the interpreter: when a new identifier (e.g x, or y as above) is input, the corresponding variable is registered as having a strictly lower priority than any variable in use at this pointFOOTNOTE<<< This is not strictly true: the variable x is predefined and always has the highest possible priority. >>> . To see the ordering used by gp at any given time, type variable().

Given such an ordering, multivariate polynomials are stored so that the variable with the highest priority is the main variable. And so on, recursively, until all variables are exhausted. A different storage pattern (which could only be obtained via libpari programming and low-level constructors) would produce an invalid object, and eventually a disaster.

In any case, if you are working with expressions involving several variables and want to have them ordered in a specific manner in the internal representation just described, the simplest is just to write down the variables one after the other under gp before starting any real computations. You could also define variables from your gprc to have a consistent ordering of common variable names in all your gp sessions, e.g read in a file variables.gp containing

  x;y;z;t;a;b;c;d;

@3Important note. PARI allows Euclidean division of multivariate polynomials, but assumes that the computation takes place in the fraction field of the coefficient ring (if it is not an integral domain, the result will a priori not make sense). This can become tricky; for instance assume x has highest priority (which is always the case), then y:

  ? x % y
  %1 = 0
  ? y % x
  %2 = y             \\ these two take place in B<Q>(y)[x]
  ? x * Mod(1,y)
  %3 = Mod(1, y)*x   \\ in (B<Q>(y)/yB<Q>(y))[x]  ~  B<Q>[x]
  ? Mod(x,y)
  %4 = 0

@3In the last example, the division by y takes place in Q(y)[x], hence the Mod object is a coset in (Q(y)[x]) / (yQ(y)[x]), which is the null ring since y is invertible! So be very wary of variable ordering when your computations involve implicit divisions and many variables. This also affects functions like numerator/denominator or content:

  ? denominator(x / y)
  %1 = 1
  ? denominator(y / x)
  %2 = x
  ? content(x / y)
  %3 = 1/y
  ? content(y / x)
  %4 = y
  ? content(2 / x)
  %5 = 2

@3Can you see why? Hint: x/y = (1/y) * x is in Q(y)[x] and denominator is taken with respect to Q(y)(x); y/x = (y*x^0) / x is in Q(y)(x) so y is invertible in the coefficient ring. On the other hand, 2/x involves a single variable and the coefficient ring is simply Z.

These problems arise because the variable ordering defines an implicit variable with respect to which division takes place. This is the price to pay to allow % and / operators on polynomials instead of requiring a more cumbersome divrem(x, y, var) (which also exists). Unfortunately, in some functions like content and denominator, there is no way to set explicitly a main variable like in divrem and remove the dependence on implicit orderings. This will hopefully be corrected in future versions.

Multivariate power series

Just like multivariate polynomials, power series are fundamentally single-variable objects. It is awkward to handle many variables at once, since PARI's implementation cannot handle multivariate error terms like O(x^i y^j). (It can handle the polynomial O(y^j) x x^i which is a very different thing, see below.)

The basic assumption in our model is that if variable x has higher priority than y, then y does not depend on x: setting y to a function of x after some computations with bivariate power series does not make sense a priori. This is because implicit constants in expressions like O(x^i) depend on y (whereas in O(y^j) they can not depend on x). For instance

    ? O(x) * y
    %1 = O(x)
    ? O(y) * x
    %2 = O(y)*x

Here is a more involved example:

    ? A = 1/x^2 + 1 + O(x); B = 1/x + 1 + O(x^3);
    ? subst(z*A, z, B)
    %2 = x^-3 + x^-2 + x^-1 + 1 + O(x)
    ? B * A
    %3 = x^-3 + x^-2 + x^-1 + O(1)
    ? z * A
    %4 = z*x^-2 + z + O(x)

The discrepancy between %2 and %3 is surprising. Why does %2 contain a spurious constant term, which cannot be deduced from the input? Well, we ignored the rule that forbids to substitute an expression involving high-priority variables to a low-priority variable. The result %4 is correct according to our rules since the implicit constant in O(x) may depend on z. It is obviously wrong if z is allowed to have negative valuation in x. Of course, the correct error term should be O(xz), but this is not possible in PARI.


Variables and Scope

This section is rather technical, and strives to explain potentially confusing concepts. Skip to the last subsection for practical advice, if the next discussion does not make sense to you. After learning about user functions, study the example in Label se:bewarescope then come back.

@3Definitions.

A scope is an enclosing context where names and values are associated. A user's function body, the body of a loop, an individual command line, all define scopes; the whole program defines the global scope. The argument of eval is evaluated in the enclosing scope.

Variables are bound to values within a given scope. This is traditionally implemented in two different ways:

@3* lexical (or static) scoping: the binding makes sense within a given block of program text. The value is private to the block and may not be accessed from outside. Where to find the value is determined at compile time.

@3* dynamic scoping: introducing a local variable, say x, pushes a new value on a stack associated to the name x (possibly empty at this point), which is popped out when the control flow leaves the scope. Evaluating x in any context, possibly outside of the given block, always yields the top value on this dynamic stack.

GP implements both lexical and dynamic scoping, using the keywords FOOTNOTE<<< The names are borrowed from the Perl scripting language. >>> my (lexical) and local (dynamic):

    x = 0;
    f() = x
    g() =    my(x = 1); f()
    h() = local(x = 1); f()

The function g returns 0 since the global x binding is unaffected by the introduction of a private variable of the same name in g. On the other hand, h returns 1; when it calls f(), the binding stack for the x identifier contains two items: the global binding to 0, and the binding to 1 introduced in h, which is still present on the stack since the control flow has not left h yet.

Scoping rules

Named parameters in a function definition, as well as all loop indicesFOOTNOTE<<< More generally, in all iterative constructs which use a variable name (for, prod, sum, vector, matrix, plot, etc.) the given variable is lexically scoped to the construct's body. >>>, have lexical scope within the function body and the loop body respectively.

  p = 0;
  forprime (p = 2, 11, print(p)); p   \\ prints 0 at the end
  x = 0;
  f(x) = x++;
  f(1)  \\ returns 2, and leave global x unaffected (= 0)

@3 If you exit the loop prematurely, e.g. using the break statement, you must save the loop index in another variable since its value prior the loop will be restored upon exit. For instance

    for(i = 1, n,
      if (ok(i), break);
    );
    if (i > n, return(failure));

@3 is incorrect, since the value of i tested by the (i > n) is quite unrelated to the loop index. One ugly workaround is

    for(i = 1, n,
      if (ok(i), isave = i; break);
    );
    if (isave > n, return(failure));

@3 But it is usually more natural to wrap the loop in a user function and use return instead of break:

  try() =
  {
    for(i = 1, n,
      if (ok(i), return (i));
    );
    0 \\ failure
  }

A list of variables can be lexically or dynamically scoped (to the block between the declaration and the end of the innermost enclosing scope) using a my or local declaration:

  for (i = 1, 10,
    my(x, y, z, i2 = i^2); \\ temps needed within the loop body
    ...
  )

Note how the declaration can include (optional) initial values, i2 = i^2 in the above. Variables for which no explicit default value is given in the declaration are initialized to 0. It would be more natural to initialize them to free variables, but this would break backward compatibility. To obtain this behavior, you may explicitly use the quoting operator:

  my(x = 'x, y = 'y, z = 'z);

A more complicated example:

  for (i = 1, 3,
    print("main loop");
    my(x = i);          \\ local to the outermost loop
    for (j = 1, 3,
      my (y = x^2);     \\ local to the innermost loop
      print (y + y^2);
      x++;
    )
  )

When we leave the loops, the values of x, y, i, j are the same as before they were started.

Note that eval is evaluated in the given scope, and can access values of lexical variables:

  ? x = 1;
  ? my(x = 0); eval("x")
  %2 = 0    \\ we see the local C<x> scoped to this command line, not the global one

Variables dynamically scoped using local should more appropriately be called temporary values since they are in fact local to the function declaring them and any subroutine called from within. In practice, you almost certainly want true private variables, hence should use almost exclusively my.

We strongly recommended to explicitly scope (lexically) all variables to the smallest possible block. Should you forget this, in expressions involving such ``rogue'' variables, the value used will be the one which happens to be on top of the value stack at the time of the call; which depends on the whole calling context in a non-trivial way. This is in general not what you want.


User defined functions

The most important thing to understand about user-defined functions is that they are ordinary GP objects, bound to variables just like any other object. Those variables are subject to scoping rules as any other: while you can define all your functions in global scope, it is usually possible and cleaner to lexically scope your private helper functions to the block of text where they will be needed.

Whenever gp meets a construction of the form expr(argument list) and the expression expr evaluates to a function (an object of type t_CLOSURE), the function is called with the proper arguments. For instance, constructions like funcs[i](x) are perfectly valid, assuming funcs is an array of functions.

Defining a function

A user function is defined as follows:

(list of formal variables) - > seq.

@3The list of formal variables is a comma-separated list of distinct variable names and allowed to be empty. It there is a single formal variable, the parentheses are optional. This list corresponds to the list of parameters you will supply to your function when calling it.

In most cases you want to assign a function to a variable immediately, as in

  R = (x,y) -> sqrt( x^2+y^2 );
  sq = x -> x^2;  \\ or equivalently C<(x) - > x^2>

but it is quite possible to define (a priori short-lived) anonymous functions. The trailing semicolon is not part of the definition, but as usual prevents gp from printing the result of the evaluation, i.e. the function object. The construction

f(list of formal variables) = seq

@3is available as an alias for

f = (list of formal variables) - > seq

@3Using that syntax, it is not possible to define anonymous functions (obviously), and the above two examples become:

  R(x,y) = sqrt( x^2+y^2 );
  sq(x) = x^2;

The semicolon serves the same purpose as above: preventing the printing of the resulting function object; compare

  ? sq(x) = x^2;  \\ no output
  ? sq(x) = x^2   \\ print the result: a function object
  %2 = (x)->x^2

@3Of course, the sequence seq can be arbitrarily complicated, in which case it will look better written on consecutive lines, with properly scoped variables:part { f(x_0, x_1,...) = my(t_0, t_1,...); \@com variables lexically scoped to the function body ... } @eprog@3Note that the following variant would also work:part f(x_0, x_1,...) = { my(t_0, t_1,...); \@com variables lexically scoped to the function body ... } @eprog

(the first newline is disregarded due to the preceding = sign, and the others because of the enclosing braces). The my statements can actually occur anywhere within the function body, scoping the variables to more restricted blocks than the whole function body.

Arguments are passed by value, not as variables: modifying a function's argument in the function body is allowed, but does not modify its value in the calling scope. In fact, a copy of the actual parameter is assigned to the formal parameter when the function is called. Formal parameters are lexically scoped to the function body. It is not allowed to use the same variable name for different parameters of your function:

  ? f(x,x) = 1
    ***   variable declared twice: f(x,x)=1
    ***                                ^----

@3Finishing touch. You can add a specific help message for your function using addhelp, but the online help system already handles it. By default ?name will print the definition of the function name: the list of arguments, as well as their default values, the text of seq as you input it. Just as \c prints the list of all built-in commands, \u outputs the list of all user-defined functions.

@3Backward compatibility (lexical scope). Lexically scoped variables were introduced in version 2.4.2. Before that, the formal parameters were dynamically scoped. If your script depends on this behavior, you may use the following trick: replace the initial f(x) = by

  f(x_orig) = local(x = x_orig)

@3Backward compatibility (disjoint namespaces). Before version 2.4.2, variables and functions lived in disjoint namespaces and it was not possible to have a variable and a function share the same name. Hence the need for a kill function allowing to reuse symbols. This is no longer the case.

There is now no distinction between variable and function names: we have PARI objects (functions of type t_CLOSURE, or more mundane mathematical entities, like t_INT, etc.) and variables bound to them. There is nothing wrong with the following sequence of assignments:

  ? f = 1       \\ assigns the integer 1 to C<f>
  %1 = 1;
  ? f() = 1     \\ a function with a constant value
  %2 = ()->1
  ? f = x^2     \\ C<f> now holds a polynomial
  %3 = x^2
  ? f(x) = x^2  \\...and now a polynomial function
  %4 = (x)->x^2
  ? g(fun) = fun(Pi);\\ a function taking a function as argument
  ? g(cos)
  %6 = -1.000000000000000000000000000

@3 Previously used names can be recycled as above: you are just redefining the variable. The previous definition is lost of course.

@3Important technical note. Built-in functions are a special case since they are read-only (you cannot overwrite their default meaning), and they use features not available to user functions, in particular pointer arguments. In the present version 2.7.1, it is possible to assign a built-in function to a variable, or to use a built-in function name to create an anonymous function, but some special argument combinations may not be available:

  ? issquare(9, &e)
  %1 = 1
  ? e
  %2 = 3
  ? g = issquare;
  ? g(9)
  %4 = 1
  ? g(9, &e)  \\ pointers are not implemented for user functions
    ***   unexpected &: g(9,&e)
    ***                     ^---

Function call, Default arguments

You may now call your function, as in f(1,2), supplying values for the formal variables. The number of parameters actually supplied may be less than the number of formal variables in the function definition. An uninitialized formal variable is given an implicit default value of (the integer) 0, i.e. after the definition

  f(x, y) = ...

you may call f(1, 2), supplying values for the two formal parameters, or for example

  f(2) equivalent to f(2,0),

  f() f(0,0),

  f(,3) f(0,3). (``Empty argument'' trick) @3This implicit default value of 0, is actually deprecated and setting

    default(strictargs, 1)

@3allows to disable it (see \secref{se:def,strictargs}).

The recommended practice is to explicitly set a default value: in the function definition, you can append = expr to a formal parameter, to give that variable a default value. The expression gets evaluated the moment the function is called, and may involve the preceding function parameters: a default value for x_i may involve x_j for j < i. For instance, after

  f(x = 1, y = 2, z = y+1) = ....

typing in f(3,4) would give you f(3,4,5). In the rare case when you want to set some far away argument, and leave the defaults in between as they stand, use the ``empty argument'' trick: f(6,,1) would yield f(6,2,1). Of course, f() by itself yields f(1,2,3) as was to be expected.

In short, the argument list is filled with user supplied values, in order. A comma or closing parenthesis, where a value should have been, signals we must use a default value. When no input arguments are left, the defaults are used instead to fill in remaining formal parameters. A final example:

  f(x, y=2, z=3) = print(x, ":", y, ":", z);

@3defines a function which prints its arguments (at most three of them), separated by colons.

  ? f(6,7)
  6:7:3
  ? f(,5)
  0:5:3
  ? f()
  0:2:3

@3If strictargs is set (recommended), x is now a mandatory argument, and the above becomes:

  ? default(strictargs,1)
  ? f(6,7)
  6:7:3
  ? f(,5)
    ***   at top-level: f(,5)
    ***                 ^-----
    ***   in function f: x,y=2,z=3
    ***                  ^---------
    ***   missing mandatory argument 'x' in user function.

@3Example. We conclude with an amusing example, intended to illustrate both user-defined functions and the power of the sumalt function. Although the Riemann zeta-function is included (as zeta) among the standard functions, let us assume that we want to check other implementations. Since we are highly interested in the critical strip, we use the classical formula

   (2^{1-s} - 1)zeta(s) = sum_{n >= 1} (-1)^n n^{-s}, Re s > 0.

The implementation is obvious:

  ZETA(s) = sumalt(n=1, (-1)^n*n^(-s)) / (2^(1-s) - 1)

Note that n is automatically lexically scoped to the sumalt ``loop'', so that it is unnecessary to add a my(n) declaration to the function body. Surprisingly, this gives very good accuracy in a larger region than expected:

  ? check = z -> ZETA(z) / zeta(z);
  ? check(2)
  %1 = 1.000000000000000000000000000
  ? check(200)
  %2 = 1.000000000000000000000000000
  ? check(0)
  %3 = 0.9999999999999999999999999994
  ? check(-5)
  %4 = 1.00000000000000007549266557
  ? check(-11)
  %5 = 0.9999752641047824902660847745
  ? check(1/2+14.134*I)  \\ very close to a non-trivial zero
  %6 = 1.000000000000000000003747432 + 7.62329066 E-21*I
  ? check(-1+10*I)
  %7 = 1.000000000000000000000002511 + 2.989950968 E-24*I

@3Now wait a minute; not only are we summing a series which is certainly no longer alternating (it has complex coefficients), but we are also way outside of the region of convergence, and still get decent results! No programming mistake this time: sumalt is a ``magic'' functionFOOTNOTE<<< sumalt is heuristic, but its use can be rigorously justified for a given function, in particular our zeta(s) formula. Indeed, Peter Borwein (An efficient algorithm for the Riemann zeta function, CMS Conf. Proc. 27 (2000), pp. 29--34) proved that the formula used in sumalt with n terms computes (1-2^{1-s})zeta(s) with a relative error of the order of (3+ sqrt {8})^{-n}|Gamma(s)|^{-1}. >>>, providing very good convergence acceleration; in effect, we are computing the analytic continuation of our original function. To convince ourselves that sumalt is a non-trivial implementation, let us try a simpler example:

  ? sum(n=1, 10^7, (-1)^n/n, 0.) / (-log(2)) \\ approximates the well-known formula
  time = 7,417 ms.
  %1 = 0.9999999278652515622893405457
  ? sumalt(n=1, (-1)^n/n) / (-log(2))        \\ accurate and fast
  time = 0 ms.
  %2 = 1.000000000000000000000000000

@3No, we are not using a powerful simplification tool here, only numerical computations. Remember, PARI is not a computer algebra system!

Beware scopes

Be extra careful with the scopes of variables. What is wrong with the following definition?

  FirstPrimeDiv(x) =
  { my(p);
    forprime(p=2, x, if (x%p == 0, break));
    p
  }
  ? FirstPrimeDiv(10)
  %1 = 0

@3\misctitle{Hint} The function body is equivalent to

  { my(newp = 0);
    forprime(p=2, x, if (x%p == 0, break));
    newp
  }

@3Detailed explanation. The index p in the forprime loop is lexically scoped to the loop and is not visible to the outside world. Hence, it will not survive the break statement. More precisely, at this point the loop index is restored to its preceding value. The initial my(p), although well-meant, adds to the confusion: it indeed scopes p to the function body, with initial value 0, but the forprime loop introduces another variable, unfortunately also called p, scoped to the loop body, which shadows the one we wanted. So we always return 0, since the value of the p scoped to the function body never changes and is initially 0.

To sum up, the routine returns the p declared local to it, not the one which was local to forprime and ran through consecutive prime numbers. Here is a corrected version:

  ? FirstPrimeDiv(x) = forprime(p=2, x, if (x%p == 0, return(p)))

Recursive functions

Recursive functions can easily be written as long as one pays proper attention to variable scope. Here is an example, used to retrieve the coefficient array of a multivariate polynomial (a non-trivial task due to PARI's unsophisticated representation for those objects):

  coeffs(P, nbvar) =
  {
    if (type(P) != "t_POL",
      for (i=1, nbvar, P = [P]);
      return (P)
    );
    vector(poldegree(P)+1, i, coeffs(polcoeff(P, i-1), nbvar-1))
  }

@3If P is a polynomial in k variables, show that after the assignment v = coeffs(P,k), the coefficient of x_1^{n_1}... x_k^{n_k} in P is given by v[n_1+1][...][n_k+1].

The operating system automatically limits the recursion depth:

  ? dive(n) = dive(n+1)
  ? dive(0);
    ***   [...] at: dive(n+1)
    ***             ^---------
    ***   in function dive: dive(n+1)
    ***                     ^---------
    \\ (last 2 lines repeated 19 times)
    ***   deep recursion.

There is no way to increase the recursion limit (which may be different on your machine) from within gp. To increase it before launching gp, you can use ulimit or limit, depending on your shell, and raise the process available stack space (increase stacksize).

Function which take functions as parameters This is done as follows

  ? calc(f, x) = f(x)
  ? calc(sin, Pi)
  %2 = -5.04870979 E-29
  ? g(x) = x^2;
  ? calc(g, 3)
  %4 = 9

@3If we do not need g elsewhere, we should use an anonymous function here, calc(x- > x^2, 3). Here is a variation:

  ? funs = [cos, sin, tan, x->x^3+1]; \\ an array of functions
  ? call(i, x) = funs[i](x)

evaluates the appropriate function on argument x, provided 1 <= i <= 4. Finally, a more useful example:

  APPLY(f, v) = vector(#v, i, f(v[i]))

applies the function f to every element in the vector v. (The built-in function apply is more powerful since it also applies to lists and matrices.)

Defining functions within a function

Defining a single function is easy:

  init(x) = (add = y -> x+y);

@3Basically, we are defining a global variable add whose value is the function y- > x+y. The parentheses were added for clarity and are not mandatory.

  ? init(5);
  ? add(2)
  %2 = 7

A more refined approach is to avoid global variables and return the function:

  init(x) = y -> x+y
  add = init(5)

@3Then add(2) still returns 7, as expected! Of course, if add is in global scope, there is no gain, but we can lexically scope it to the place where it is useful:

    my ( add = init(5) );

How about multiple functions then? We can use the last idea and return a vector of functions, but if we insist on global variables? The first idea

  init(x) = add(y) = x+y; mul(y) = x*y;

@3does not work since in the construction f() = seq, the function body contains everything until the end of the expression. Hence executing init defines the wrong function add (itself defining a function mul). The way out is to use parentheses for grouping, so that enclosed subexpressions will be evaluated independently:

  ? init(x) = ( add(y) = x+y ); ( mul(y) = x*y );
  ? init(5);
  ? add(2)
  %3 = 7
  ? mul(3)
  %4 = 15

@3This defines two global functions which have access to the lexical variables private to init! The following would work in exactly the same way:

  ? init5() = my(x = 5); ( add(y) = x+y ); ( mul(y) = x*y );

Closures as Objects

Contrary to what you might think after the preceding examples, GP's closures may not be used to simulate true ``objects'', with private and public parts and methods to access and manipulate them. In fact, closures indeed incorporate an existing context (they may access lexical variables that existed at the time of their definition), but then may not change it. More precisely, they access a copy, which they are welcome to change, but a further function call still accesses the original context, as it existed at the time the function was defined:

  init() =
  { my(count = 0);
    inc()=count++;
    dec()=count--;
  }
  ? inc()
  %1 = 1
  ? inc()
  %2 = 1
  ? inc()
  %3 = 1


Member functions

Member functions use the `dot' notation to retrieve information from complicated structures. The built-in structures are bid, ell, galois, ff, nf, bnf, bnr and prid, which will be described at length in Chapter 3. The syntax structure.member is taken to mean: retrieve member from structure, e.g. E.j returns the j-invariant of the elliptic curve E, or outputs an error message if E is not a proper ell structure. To define your own member functions, use the syntax

 I<var>C<.>I<member>C< = >I<seq>,

@3where the formal variable var is scoped to the function body seq. This is of course reminiscent of a user function with a single formal variable var. For instance, the current implementation of the ell type is a vector, the j-invariant being the thirteenth component. It could be implemented as

  x.j =
  {
    if (type(x) != "t_VEC" || #x < 14, error("not an elliptic curve: " x));
    x[13]
  }

@3As for user functions, you can redefine your member functions simply by typing new definitions. On the other hand, as a safety measure, you cannot redefine the built-in member functions, so attempting to redefine x.j as above would in fact produce an error; you would have to call it e.g. x.myj in order for gp to accept it.

@3Rationale. In most cases, member functions are simple accessors of the form

    x.a = x[1];
    x.b = x[2];
    x.c = x[3];

where x is a vector containing relevant data. There are at least three alternative approaches to the above member functions: 1) hardcode x[1], etc. in the program text, 2) define constant global variables AINDEX = 1, BINDEX = 2 and hardcode x[AINDEX], 3) user functions a(x) = x[1] and so on.

Even if 2) improves on 1), these solutions are neither elegant nor flexible, and they scale badly. 3) is a genuine possibility, but the main advantage of member functions is that their namespace is independent from the variables (and functions) namespace, hence we can use very short identifiers without risk. The j-invariant is a good example: it would clearly not be a good idea to define j(E) = E[13], because clashes with loop indices are likely.

@3Note. Typing \um will output all user-defined member functions.

@3Member function names. A valid name starts with a letter followed by any number of keyword characters: _ or alphanumeric characters ([A-Za-z0-9]). The built-in member function names are reserved and cannot be used (see the list with ?.). Finally, names starting with e or E followed by a digit are forbidden, due to a clash with the floating point exponent notation: we understand 1.e2 as 100.000..., not as extracting member e2 of object 1.


Strings and Keywords

Strings

GP variables can hold values of type character string (internal type t_STR). This section describes how they are actually used, as well as some convenient tricks (automatic concatenation and expansion, keywords) valid in string context.

As explained above, the general way to input a string is to enclose characters between quotes ". This is the only input construct where whitespace characters are significant: the string will contain the exact number of spaces you typed in. Besides, you can ``escape'' characters by putting a \ just before them; the translation is as follows

     \e: <Escape>
     \n: <Newline>
     \t: <Tab>

For any other character x, \x is expanded to x. In particular, the only way to put a " into a string is to escape it. Thus, for instance, "\"a\"" would produce the string whose content is ``a''. This is definitely not the same thing as typing "a", whose content is merely the one-letter string a.

You can concatenate two strings using the concat function. If either argument is a string, the other is automatically converted to a string if necessary (it will be evaluated first).

  ? concat("ex", 1+1)
  %1 = "ex2"
  ? a = 2; b = "ex"; concat(b, a)
  %2 = "ex2"
  ? concat(a, b)
  %3 = "2ex"

Some functions expect strings for some of their arguments: print would be an obvious example, Str is a less obvious but useful one (see the end of this section for a complete list). While typing in such an argument, you will be said to be in string context. The rest of this section is devoted to special syntactical tricks which can be used with such arguments (and only here; you will get an error message if you try these outside of string context):

@3* Writing two strings alongside one another will just concatenate them, producing a longer string. Thus it is equivalent to type in "a " "b" or "a b". A little tricky point in the first expression: the first whitespace is enclosed between quotes, and so is part of a string; while the second (before the "b") is completely optional and gp actually suppresses it, as it would with any number of whitespace characters at this point (i.e. outside of any string).

@3* If you insert any expression when a string is expected, it gets ``expanded'': it is evaluated as a standard GP expression, and the final result (as would have been printed if you had typed it by itself) is then converted to a string, as if you had typed it directly. For instance "a" 1+1 "b" is equivalent to "a2b": three strings get created, the middle one being the expansion of 1+1, and these are then concatenated according to the rule described above. Another tricky point here: assume you did not assign a value to aaa in a GP expression before. Then typing aaa by itself in a string context will actually produce the correct output (i.e. the string whose content is aaa), but in a fortuitous way. This aaa gets expanded to the monomial of degree one in the variable aaa, which is of course printed as aaa, and thus will expand to the three letters you were expecting.

@3Warning. Expression involving strings are not handled in a special way; even in string context, the largest possible expression is evaluated, hence print("a"[1]) is incorrect since "a" is not an object whose first component can be extracted. On the other hand print("a", [1]) is correct (two distinct argument, each converted to a string), and so is print("a" 1) (since "a"1 is not a valid expression, only "a" gets expanded, then 1, and the result is concatenated as explained above).

Keywords

Since there are cases where expansion is not desirable, we now distinguish between ``Keywords'' and ``Strings''. String is what has been described so far. Keywords are special relatives of Strings which are automatically assumed to be quoted, whether you actually type in the quotes or not. Thus expansion is never performed on them. They get concatenated, though. The analyzer supplies automatically the quotes you have ``forgotten'' and treats Keywords just as normal strings otherwise. For instance, if you type "a"b+b in Keyword context, you will get the string whose contents are ab+b. In String context, on the other hand, you would get a2*b.

All GP functions have prototypes (described in Chapter 3 below) which specify the types of arguments they expect: either generic PARI objects (GEN), or strings, or keywords, or unevaluated expression sequences. In the keyword case, only a very small set of words will actually be meaningful (the default function is a prominent example).

@3Reference. The arguments of the following functions are processed in string context:

  Str

  addhelp (second argument)

  default (second argument)

  error

  extern

  plotstring (second argument)

  plotterm (first argument)

  read and readvec

  system

  all the printxxx functions

  all the writexxx functions

@3The arguments of the following functions are processed as keywords:

  alias

  default (first argument)

  install (all arguments but the last)

  trap (first argument)

  whatnow

Useful examples

The function Str converts its arguments into strings and concatenate them. Coupled with eval, it is very powerful. The following example creates generic matrices:

  ? genmat(u,v,s="x") = matrix(u,v,i,j, eval( Str(s,i,j) ))
  ? genmat(2,3) + genmat(2,3,"m")
  %1 =
  [x11 + m11 x12 + m12 x13 + m13]
  [x21 + m21 x22 + m22 x23 + m23]

Two last examples: hist(10,20) returns all history entries from %10 to %20 neatly packed into a single vector; histlast(10) returns the last 10 history entries:

    hist(a,b) = vector(b-a+1, i, eval(Str("%", a-1+i)))
    histlast(n) = vector(n, i, eval(Str("%", %#-i+1)))


Errors and error recovery

Errors

Your input program is first compiled to a more efficient bytecode; then the latter is evaluated, calling appropriate functions from the PARI library. Accordingly, there are two kind of errors: syntax errors produced by the compiler, and runtime errors produced by the PARI library either by the evaluator itself, or in a mathematical function. Both kinds are fatal to your computation: gp will report the error, perform some cleanup (restore variables modified while evaluating the erroneous command, close open files, reclaim unused memory, etc.), and output its usual prompt.

When reporting a syntax error, gp gives meaningful context by copying (part of) the expression it was trying to compile, indicating where the error occurred with a caret ^-, as in

  ? factor()
    ***   too few arguments: factor()
    ***                             ^-
  ? 1+
    ***   syntax error, unexpected $end: 1+
    ***                                   ^-

possibly enlarged to a full arrow given enough trailing context

  ? if (isprime(1+, do_something())
    ***   syntax error, unexpected {\tt ','}: if(isprime(1+,do_something()))
    ***                                              ^----------------

These error messages may be mysterious, because gp cannot guess what you were trying to do, and the error may occur once gp has been sidetracked. The first error is straightforward: factor has one mandatory argument, which is missing.

The other two are simple typos involving an ill-formed addition 1 + missing its second operand. The error messages differ because the parsing context is slightly different: in the first case we reach the end of input ($end) while still expecting a token, and in the second one, we received an unexpected token (the comma).

Here is a more complicated one:

  ? factor(x
    ***   syntax error, unexpected C<end, expecting )-> or {\tt ','} or {\tt ')'}: factor(x
    ***                                                                      ^-

The error is a missing parenthesis, but from gp's point of view, you might as well have intended to give further arguments to factor (this is possible and useful, see the description of the function). In fact gp expected either a closing parenthesis, or a second argument separated from the first by a comma. And this is essentially what the error message says: we reached the end of the input (end>) while expecting a ')' or a ','.

Actually, a third possibility is mentioned in the error message )- > , which could never be valid in the above context, but a subexpression like (x)- > sin(x), defining an inline closure would be valid, and the parser is not clever enough to rule that out, so we get the same message as in

  ? (x
    ***   syntax error, unexpected $end, expecting )-> or {\tt ','} or {\tt ')'}: (x
    ***                                                                ^-

where all three proposed continuations would be valid.

Runtime errors from the evaluator are nicer because they answer a correctly worded query, otherwise the bytecode compiler would have protested first; here is a slightly pathological case:

  ? if (siN(x) < eps, do_something())
    ***   at top-level: if(siN(x)<eps,do_someth
    ***                    ^--------------------
    ***   not a function in function call

@3(no arrow!) The code is syntactically correct and compiled correctly, even though the siN function, a typo for sin, was not defined at this point. When trying to evaluate the bytecode, however, it turned out that siN is still undefined so we cannot evaluate the function call siN(x).

Library runtime errors are even nicer because they have more mathematical content, which is easier to grasp than a parser's logic:

  ? 1/Mod(2,4)
    ***   at top-level: 1/Mod(2,4)
    ***                  ^---------
    *** _/_: impossible inverse in Fp_inv: Mod(2, 4).

@3telling us that a runtime error occurred while evaluating the binary / operator (the _ surrounding the operator are placeholders), more precisely the Fp_inv library function was fed the argument Mod(2,4) and could not invert it. More context is provided if the error occurs deep in the call chain:

  ? f(x) = 1/x;
  ? g(N) = for(i = -N, N, f(i + O(5)));
  ? g(10)
    ***   at top-level: g(10)
    ***                 ^-----
    ***   in function g: for(i=-N,N,f(i))
    ***                             ^-----
    ***   in function f: 1/x
    ***                   ^--
    *** _/_: impossible inverse in ginv: O(5).

@3In this example, the debugger reports (at least) 3 enclosed frames: last (innermost) is the body of user function f, the body of g, and the top-level (global scope). In fact, the for loop in g's body defines an extra frame, since there exist variables scoped to the loop body.

Error recovery

It is annoying to wait for a program to finish and find out the hard way that there was a mistake in it (like the division by 0 above), sending you back to the prompt. First you may lose some valuable intermediate data. Also, correcting the error may not be obvious; you might have to change your program, adding a number of extra statements and tests to narrow down the problem.

A different situation, still related to error recovery, is when you actually foresee that some error may occur, are unable to prevent it, but quite capable of recovering from it, given the chance. Examples include lazy factorization, where you knowingly use a pseudo prime N as if it were prime; you may then encounter an ``impossible'' situation, but this would usually exhibit a factor of N, enabling you to refine the factorization and go on. Or you might run an expensive computation at low precision to guess the size of the output, hence the right precision to use. You can then encounter errors like ``precision loss in truncation'', e.g when trying to convert 1E1000, known to 28 digits of accuracy, to an integer; or ``division by 0'', e.g inverting 0E1000 when all accuracy has been lost, and no significant digit remains. It would be enough to restart part of the computation at a slightly higher precision.

We now describe error trapping, a useful mechanism which alleviates much of the pain in the first situation (the break loop debugger), and provides satisfactory ways out of the second one (the iferr exception handler).

Break loop

A break loop is a special debugging mode that you enter whenever a user interrupt (Control-C) or runtime error occurs, freezing the gp state, and preventing cleanup until you get out of the loop. By runtime error, we mean an error from the evaluator, the library or a user error (from error), not syntax errors. When a break loop starts, a prompt is issued (break > ). You can type in a gp command, which is evaluated when you hit the < Return > key, and the result is printed as during the main gp loop, except that no history of results is kept. Then the break loop prompt reappears and you can type further commands as long as you do not exit the loop. If you are using readline, the history of commands is kept, and line editing is available as usual. If you type in a command that results in an error, you are sent back to the break loop prompt: errors do not terminate the loop.

To get out of a break loop, you can use next, break, return, or type C-d (EOF), any of which will let gp perform its usual cleanup, and send you back to the gp prompt. Note that C-d is slightly dangerous, since typing it twice will not only send you back to the gp prompt, but to your shell prompt! (Since C-d at the gp prompt exits the gp session.)

If the break loop was started by a user interrupt Control-C, and not by an error, inputting an empty line, i.e hitting the < Return > key at the break > prompt, resumes the temporarily interrupted computation. A single empty line has no effect in case of a fatal error, to avoid getting get out of the loop prematurely, thereby losing valuable debugging data. Any of next, break, return, or C-d will abort the computation and send you back to the gp prompt as above.

Break loops are useful as a debugging tool. You may inspect the values of gp variables to understand why an error occurred, or change gp's state in the middle of a computation (increase debugging level, start storing results in a log file, set variables to different values...): hit C-c, type in your modifications, then let the computation go on as explained above. A break loop looks like this:

  ? v = 0; 1/v
    ***   at top-level: v=0;1/v
    ***                      ^--
    *** _/_: impossible inverse in gdiv: 0.
    ***   Break loop (type 'break' to go back to the GP prompt)
  break>

@3So the standard error message is printed first. The break > at the bottom is a prompt, and hitting v then < Return > , we see:

  break> v
  0

@3explaining the problem. We could have typed any gp command, not only the name of a variable, of course. Lexically-scoped variables are accessible to the evaluator during the break loop:

  ? for(v = -2, 2, print(1/v))
  -1/2
  -1
    ***   at top-level: for(v=-2,2,print(1/v))
    ***                                   ^----
    *** _/_: impossible inverse in gdiv: 0.
    ***   Break loop (type 'break' to go back to the GP prompt)
  break> v
  0

Even though loop indices are automatically lexically scoped and no longer exist when the break loop is run, enough debugging information is retained in the bytecode to reconstruct the evaluation context. Of course, when the error occurs in a nested chain of user function calls, lexically scoped variables are available only in the corresponding frame:

  ? f(x) = 1/x;
  ? g(x) = for(i = 1, 10, f(x+i));
  ? for(j = -5,5, g(j))
    ***   at top-level: for(j=-5,5,g(j))
    ***                            ^-----
    ***   in function g: for(i=1,10,f(x+i))
    ***                             ^-------
    ***   in function f: 1/x
    ***                   ^--
    *** _/_: impossible inverse in gdiv: 0.
    ***   Break loop: type 'break' to go back to GP prompt
  break> [i,j,x]     \\  the x in f's body.
  [i, j, 0]
  break> dbg_up      \\  go up one frame
    ***   at top-level: for(j=-5,5,g(j))
    ***                            ^-----
    ***   in function g: for(i=1,10,f(x+i))
    ***                             ^-------
  break> [i,j,x]      \\  the x in g's body, i in the for loop.
  [5, j, -5]

The following GP commands are available during a break loop to help debugging:

dbg_up(n): go up n frames, as seen above.

dbg_down(n): go down n frames, cancelling previous dbg_up's.

dbg_x(t): examine t, as \x but more flexible.

dbg_err(): returns the current error context t_ERROR. The error components often provide useful additional information:

    ? O(2) + O(3)
      ***   at top-level: O(2)+O(3)
      ***                     ^-----
      *** _+_: inconsistent addition t_PADIC + t_PADIC.
      ***   Break loop: type 'break' to go back to GP prompt
    break> E = dbg_err()
    error("inconsistent addition t_PADIC + t_PADIC.")
    break> Vec(E)
    ["e_OP", "+", O(2), O(3)]

@3Note. The debugger is enabled by default, and fires up as soon as a runtime error occurs. If you do not like this behavior, you may disable it by setting the default breakloop to 0 in for gprc. A runtime error will send you back to the prompt. Note that the break loop is automatically disabled when running gp in non interactive mode, i.e. when the program's standard input is not attached to a terminal.

@3Technical Note. When you enter a break loop due to a PARI stack overflow, the PARI stack is reset so that you can run commands. Otherwise the stack would immediately overflow again! Still, as explained above, you do not lose the value of any gp variable in the process.

Protecting code

The expression

iferr(statements, ERR, recovery)

@3evaluates and returns the value of statements, unless an error occurs during the evaluation in which case the value of recovery is returned. As in an if/else clause, with the difference that statements has been partially evaluated, with possible side effects. We shall give a lot more details about the ERR argument shortly; it is the name of a variable, lexically scoped to the recovery expression sequence, whose value is set by the exception handler to help the recovery code decide what to do about the error.

For instance one can define a fault tolerant inversion function as follows:

  ? inv(x) = iferr(1/x, ERR, "oo")    \\ ERR is unused...
  ? for (i=-1,1, print(inv(i)))
  -1
  oo
  1

@3Protected codes can be nested without adverse effect. Let's now see how ERR can be used; as written, inv is too tolerant:

  ? inv("blah")
  %2 = "oo"

Let's improve it by checking that we caught a ``division by 0'' exception, and not an unrelated one like the type error 1 / "blah".

  ? inv2(x) = {
    iferr(1/x,
          ERR, if (errname(ERR) != "e_INV", error(ERR)); "oo")
  }
  ? inv2(0)
  %3 = "oo"  \\ as before
  ? inv2("blah")
    ***   at top-level: inv2("blah")
    ***                 ^------------
    ***   in function inv2: ...f(errname(ERR)!="e_INV",error(ERR));"oo")
    ***                                                 ^-----------------
    *** error: forbidden division t_INT / t_STR.

@3In the inv2("blah") example, the error type was not expected, so we rethrow the exception: error(ERR) triggers the original error that we mistakenly trapped. Since the recovery code should always check whether the error is the one expected, this construction is very common and can be simplified to

  ? inv3(x) = iferr(1/x,
                    ERR, "oo",
                    errname(ERR) == "e_INV")

@3More generally

iferr(statements, ERR, recovery, predicate)

@3only catches the exception if predicate (allowed to check various things about ERR, not only its name) is non-zero.

Rather than trapping everything, then rethrowing whatever we do not like, we advise to only trap errors of a specific kind, as above. Of course, sometimes, one just want to trap everything because we do not know what to expect. The following function check whether install works correctly in your gp:

  broken_install() =
  { \\ can we install?
    iferr(install(addii,GG),
          ERR, return ("OS"));
    \\ can we use the installed function?
    iferr(if (addii(1,1) != 2, return("BROKEN")),
          ERR, return("USE"));
    return (0);
  }

@3The function returns OS if the operating system does not support install, USE if using an installed function triggers an error, BROKEN if the installed function did not behave as expected, and 0 if everything works.

The ERR formal parameter contains more useful data than just the error name, which we recovered using errname(ERR). In fact, a t_ERROR object usually has extra components, which can be accessed as component(ERR,1), component(ERR,2), and so on. Or globally by casting the error to a t_VEC: Vec(ERR) returns the vector of all components at once. See Label se:iferr for the list of all exception types, and the corresponding contents of ERR.


Interfacing GP with other languages

The PARI library was meant to be interfaced with C programs. This specific use is dealt with extensively in the User's guide to the PARI library. Of course, gp itself provides a convenient interpreter to execute rather intricate scripts (see Label se:programming).

Scripts, when properly written, tend to be shorter and clearer than C programs, and are certainly easier to write, maintain or debug. You don't need to deal with memory management, garbage collection, pointers, declarations, and so on. Because of their intrinsic simplicity, they are more robust as well. They are unfortunately somewhat slower. Thus their use will remain complementary: it is suggested that you test and debug your algorithms using scripts, before actually coding them in C if speed is paramount. The GP2C compiler often eases this part.

The install command (see Label se:install) efficiently imports foreign functions for use under gp, which can of course be written using other libraries than PARI. Thus you may code only critical parts of your program in C, and still maintain most of the program as a GP script.

We are aware of three PARI-related Free Software packages to embed PARI in other languages. We neither endorse nor support any of them, but you may want to give them a try if you are familiar with the languages they are based on. The first is William Stein's Python-based SAGEFOOTNOTE<<< see http://sagemath.org/ >>> system. The second is the Math::Pari Perl module (see any CPAN mirror), written by Ilya Zakharevich. Finally, Michael Stoll has integrated PARI into CLISPFOOTNOTE<<< see http://clisp.cons.org/ >>>, which is a Common Lisp implementation by Bruno Haible, Marcus Daniels and others; this interface has been updated for pari-2 by Sam Steingold.

These provide interfaces to gp functions for use in python, perl, or Lisp programs, respectively.


Defaults

@3There are many internal variables in gp, defining how the system will behave in certain situations, unless a specific override has been given. Most of them are a matter of basic customization (colors, prompt) and will be set once and for all in your preferences file (see Label se:gprc), but some of them are useful interactively (set timer on, increase precision, etc.).

The function used to manipulate these values is called default, which is described in Label se:default. The basic syntax is

default(def, value),

@3which sets the default def to value. In interactive use, most of these can be abbreviated using gp metacommands (mostly, starting with \b), which we shall describe in the next section.

Here we will only describe the available defaults and how they are used. Just be aware that typing default by itself will list all of them, as well as their current values (see \d). Just after the default name, we give between parentheses the initial value when gp starts, assuming you did not tamper with factory settings using command-line switches or a gprc.

@3Note. The suffixes k, M or G can be appended to a value which is a numeric argument, with the effect of multiplying it by 10^3, 10^6 and 10^9 respectively. Case is not taken into account there, so for instance 30k and 30K both stand for 30000. This is mostly useful to modify or set the default parisize which typically involve a lot of trailing zeroes.

@3(somewhat technical) Note. As we saw in Label se:strings, the second argument to default is subject to string context expansion, which means you can use run-time values. In other words, something like

    a = 3;
    default(logfile, "file" a ".log")

logs the output in file3.log.

Some special defaults, corresponding to file names and prompts, expand further the resulting value at the time they are set. Two kinds of expansions may be performed:

@3* time expansion: the string is sent through the library function strftime. This means that %char combinations have a special meaning, usually related to the time and date. For instance, %H = hour (24-hour clock) and %M = minute [00,59] (on a Unix system, you can try man strftime at your shell prompt to get a complete list). This is applied to prompt, psfile, and logfile. For instance,

default(prompt,"(%H:%M) ? ")

will prepend the time of day, in the form (hh:mm) to gp's usual prompt.

@3* environment expansion: When the string contains a sequence of the form SOMEVAR, e.g. HOME, the environment is searched and if SOMEVAR is defined, the sequence is replaced by the corresponding value. Also the ~ symbol has the same meaning as in many shells --- ~ by itself stands for your home directory, and ~user is expanded to user's home directory. This is applied to all file names.

Available defaults are described in the reference guide, Label se:gp_defaults.


Simple metacommands

Simple metacommands are meant as shortcuts and should not be used in GP scripts (see Label se:programming). Beware that these, as all of gp input, are case sensitive. For example, \Q is not identical to \q. In the following list, braces are used to denote optional arguments, with their default values when applicable, e.g. {n = 0} means that if n is not there, it is assumed to be 0. Whitespace (or spaces) between the metacommand and its arguments and within arguments is optional. (This can cause problems only with \w, when you insist on having a file name whose first character is a digit, and with \r or \w, if the file name itself contains a space. In such cases, just use the underlying read or write function; see Label se:write).

\subseckbd{?{command}} The gp on-line help interface. If you type ?n where n is a number from 1 to 11, you will get the list of functions in Section 3.n of the manual (the list of sections being obtained by simply typing ?).

These names are in general not informative enough. More details can be obtained by typing ?function, which gives a short explanation of the function's calling convention and effects. Of course, to have complete information, read Chapter 3 of this manual (the source code is at your disposal as well, though a trifle less readable).

If the line before the copyright message indicates that extended help is available (this means perl is present on your system and the PARI distribution was correctly installed), you can add more ? signs for extended functionality:

?? keyword yields the function description as it stands in this manual, usually in Chapter 2 or 3. If you're not satisfied with the default chapter chosen, you can impose a given chapter by ending the keyword with @ followed by the chapter number, e.g. ?? Hello@2 will look in Chapter 2 for section heading Hello (which doesn't exist, by the way).

All operators (e.g. +, &&, etc.) are accepted by this extended help, as well as a few other keywords describing key gp concepts, e.g. readline (the line editor), integer, nf (``number field'' as used in most algebraic number theory computations), ell (elliptic curves), etc.

In case of conflicts between function and default names (e.g log, simplify), the function has higher priority. To get the default help, use

    ?? default(log)
    ?? default(simplify)

??? pattern produces a list of sections in Chapter 3 of the manual related to your query. As before, if pattern ends by @ followed by a chapter number, that chapter is searched instead; you also have the option to append a simple @ (without a chapter number) to browse through the whole manual.

If your query contains dangerous characters (e.g ? or blanks) it is advisable to enclose it within double quotes, as for GP strings (e.g ??? "elliptic curve").

Note that extended help is much more powerful than the short help, since it knows about operators as well: you can type ?? * or ?? &&, whereas a single ? would just yield a not too helpful

  &&: unknown identifier.}

@3message. Also, you can ask for extended help on section number n in Chapter 3, just by typing ?? n (where ?n would yield merely a list of functions). Finally, a few key concepts in gp are documented in this way: metacommands (e.g ?? "??"), defaults (e.g ?? psfile) and type names (e.g t_INT or integer), as well as various miscellaneous keywords such as edit (short summary of line editor commands), operator, member, "user defined", nf, ell,...

Last but not least: ?? without argument will open a dvi previewer (xdvi by default, $GPXDVI if it is defined in your environment) containing the full user's manual. ??tutorial and ??refcard do the same with the tutorial and reference card respectively.

@3Technical note. This functionality is provided by an external perl script that you are free to use outside any gp session (and modify to your liking, if you are perl-knowledgeable). It is called gphelp, lies in the doc subdirectory of your distribution (just make sure you run Configure first, see Appendix A) and is really two programs in one. The one which is used from within gp is gphelp which runs TeX on a selected part of this manual, then opens a previewer. gphelp -detex is a text mode equivalent, which looks often nicer especially on a colour-capable terminal (see misc/gprc.dft for examples). The default help selects which help program will be used from within gp. You are welcome to improve this help script, or write new ones (and we would like to know about it so that we may include them in future distributions). By the way, outside of gp you can give more than one keyword as argument to gphelp.

\subseckbd{/*...*/} A comment. Everything between the stars is ignored by gp. These comments can span any number of lines.

\subseckbd{\} A one-line comment. The rest of the line is ignored by gp.

\a {n}

Prints the object number n (%n) in raw format. If the number n is omitted, print the latest computed object (%).

\c

Prints the list of all available hardcoded functions under gp, not including operators written as special symbols (see Label se:operators). More information can be obtained using the ? metacommand (see above). For user-defined functions / member functions, see \u and \um.

\d

Prints the defaults as described in the previous section (shortcut for default(), see Label se:default).

\e {n}

Switches the echo mode on (1) or off (0). If n is explicitly given, set echo to n.

\g {n}

Sets the debugging level debug to the non-negative integer n.

\gf {n}

Sets the file usage debugging level debugfiles to the non-negative integer n.

\gm {n}

Sets the memory debugging level debugmem to the non-negative integer n.

\h {m-n}

Outputs some debugging info about the hashtable. If the argument is a number n, outputs the contents of cell n. Ranges can be given in the form m-n (from cell m to cell n, $ = last cell). If a function name is given instead of a number or range, outputs info on the internal structure of the hash cell this function occupies (a struct entree in C). If the range is reduced to a dash ('-'), outputs statistics about hash cell usage.

\l {logfile}

Switches log mode on and off. If a logfile argument is given, change the default logfile name to logfile and switch log mode on.

\m

As \a, but using prettymatrix format.

\o {n} Sets output mode to n (0: raw, 1

prettymatrix, 3: external prettyprint).

\p {n}

Sets realprecision to n decimal digits. Prints its current value if n is omitted.

\ps {n}

Sets seriesprecision to n significant terms. Prints its current value if n is omitted.

\q

Quits the gp session and returns to the system. Shortcut for quit() (see Label se:quit).

\r {filename}

Reads into gp all the commands contained in the named file as if they had been typed from the keyboard, one line after the other. Can be used in combination with the \w command (see below). Related but not equivalent to the function read (see Label se:read); in particular, if the file contains more than one line of input, there will be one history entry for each of them, whereas read would only record the last one. If filename is omitted, re-read the previously used input file (fails if no file has ever been successfully read in the current session). If a gp binary file (see Label se:writebin) is read using this command, it is silently loaded, without cluttering the history.

Assuming gp figures how to decompress files on your machine, this command accepts compressed files in compressed (.Z) or gzipped (.gz or .z) format. They will be uncompressed on the fly as gp reads them, without changing the files themselves.

\s

Prints the state of the PARI stack and heap. This is used primarily as a debugging device for PARI.

\t

Prints the internal longword format of all the PARI types. The detailed bit or byte format of the initial codeword(s) is explained in Chapter 4, but its knowledge is not necessary for a gp user.

\u

Prints the definitions of all user-defined functions.

\um

Prints the definitions of all user-defined member functions.

\v

Prints the version number and implementation architecture (680x0, Sparc, Alpha, other) of the gp executable you are using.

\w {n} {filename}

Writes the object number n ( %n ) into the named file, in raw format. If the number n is omitted, writes the latest computed object ( % ). If filename is omitted, appends to logfile (the GP function write is a trifle more powerful, as you can have arbitrary file names).

\x {n}

Prints the complete tree with addresses and contents (in hexadecimal) of the internal representation of the object number n ( %n ). If the number n is omitted, uses the latest computed object in gp. As for \s, this is used primarily as a debugging device for PARI, and the format should be self-explanatory. The underlying GP function dbg_x is more versatile, since it can be applied to other objects than history entries.

\y {n}

Switches simplify on (1) or off (0). If n is explicitly given, set simplify to n.

\subseckbd{#} Switches the timer on or off.

\subseckbd{##} Prints the time taken by the latest computation. Useful when you forgot to turn on the timer.


The preferences file

This file, called gprc in the sequel, is used to modify or extend gp default behavior, in all gp sessions: e.g customize default values or load common user functions and aliases. gp opens the gprc file and processes the commands in there, before doing anything else, e.g. creating the PARI stack. If the file does not exist or cannot be read, gp will proceed to the initialization phase at once, eventually emitting a prompt. If any explicit command line switches are given, they override the values read from the preferences file.

Syntax

The syntax in the gprc file (and valid in this file only) is simple-minded, but should be sufficient for most purposes. The file is read line by line; as usual, white space is ignored unless surrounded by quotes and the standard multiline constructions using braces, \, or = are available (multiline comments between /* ... */ are also recognized).

Preprocessor: Two types of lines are first dealt with by a preprocessor

@3* comments are removed. This applies to all text surrounded by /* ... */ as well as to everything following \\ on a given line.

@3* lines starting with #if boolean are treated as comments if boolean evaluates to false, and read normally otherwise. The condition can be negated using either #if not (or #if !). If the rest of the current line is empty, the test applies to the next line (same behavior as = under gp). Only three tests can be performed:

EMACS: true if gp is running in an Emacs or TeXmacs shell (see Label se:emacs).

READL: true if gp is compiled with readline support (see Label se:readline).

VERSION op number: where op is in the set { > , < , <= , >= }, and number is a PARI version number of the form Major.Minor.patch, where the last two components can be omitted (i.e. 1 is understood as version 1.0.0). This is true if gp's version number satisfies the required inequality.

Commands: After preprocessing, the remaining lines are executed as sequence of expressions (as usual, separated by ; if necessary). Only two kinds of expressions are recognized

@3* default = value, where default is one of the available defaults (see Label se:defaults), which will be set to value on actual startup. Don't forget the quotes around strings (e.g. for prompt or help).

@3* read "some_GP_file" where some_GP_file is a regular GP script this time, which will be read just before gp prompts you for commands, but after initializing the defaults. In particular, file input is delayed until the gprc has been fully loaded. This is the right place to input files containing alias commands, or your favorite macros.

@3For instance you could set your prompt in the following portable way:

  \\ self modifying prompt looking like \rm(18:03) B<gp>C<E<gt> >
  prompt   = "(%H:%M) \e[1mgp\e[m > "
  \\ readline wants non-printing characters to be braced between ^A/^B pairs
  #if READL prompt = "(%H:%M) ^A\e[1m^Bgp^A\e[m^B > "
  \\ escape sequences not supported under emacs
  #if EMACS prompt = "(%H:%M) gp > "

@3Note that any of the last two lines could be broken in the following way

  #if EMACS
    prompt = "(%H:%M) gp > "

@3since the preprocessor directive applies to the next line if the current one is empty.

A sample gprc file called misc/gprc.dft is provided in the standard distribution. It is a good idea to have a look at it and customize it to your needs. Since this file does not use multiline constructs, here is one (note the terminating ; to separate the expressions):

  #if VERSION > 2.2.3
  {
    read "my_scripts";     \\ syntax errors in older versions
    new_galois_format = 1; \\ default introduced in 2.2.4
  }
  #if ! EMACS
  {
    colors = "9, 5, no, no, 4, 1, 2";
    help   = "gphelp -detex -ch 4 -cb 0 -cu 2";
  }

The gprc location

When gp is started, it looks for a customization file, or gprc in the following places (in this order, only the first one found will be loaded):

@3* gp checks whether the environment variable GPRC is set. On Unix, this can be done with something like:

  GPRC = /my/dir/anyname; export GPRC in sh syntax (for instance in your .profile),

  setenv GPRC /my/dir/anyname in csh syntax (in your .login or .cshrc file).

  env GPRC = /my/dir/anyname gp on the command line launching gp.

@3If so, the file named by $GPRC is the gprc.

@3* If GPRC is not set, and if the environment variable HOME is defined, gp then tries

$HOME/.gprc on a Unix system

$HOME\gprc.txt on a DOS, OS/2, or Windows system.

@3* If no gprc was found among the user files mentioned above we look for /etc/gprc for a system-wide gprc file (you will need root privileges to set up such a file yourself).

@3* Finally, we look in pari's datadir for a file named

.gprc on a Unix system

gprc.txt on a DOS, OS/2, or Windows system. If you are using our Windows installer, this is where the default preferences file is written.

@3Note that on Unix systems, the gprc's default name starts with a '.' and thus is hidden to regular ls commands; you need to type ls -a to list it.


Using readline

This very useful library provides line editing and contextual completion to gp. You are encouraged to read the readline user manual, but we describe basic usage here.

@3A (too) short introduction to readline. In the following, C- stands for ``the Control key combined with another'' and the same for M- with the Meta key; generally C- combinations act on characters, while the M- ones operate on words. The Meta key might be called Alt on some keyboards, will display a black diamond on most others, and can safely be replaced by Esc in any case.

Typing any ordinary key inserts text where the cursor stands, the arrow keys enabling you to move in the line. There are many more movement commands, which will be familiar to the Emacs user, for instance C-a/C-e will take you to the start/end of the line, M-b/M-f move the cursor backward/forward by a word, etc. Just press the < Return > key at any point to send your command to gp.

All the commands you type at the gp prompt are stored in a history, a multiline command being saved as a single concatenated line. The Up and Down arrows (or C-p/C-n) will move you through the history, M- < /M- > sending you to the start/end of the history. C-r/C-s will start an incremental backward/forward search. You can kill text (C-k kills till the end of line, M-d to the end of current word) which you can then yank back using the C-y key (M-y will rotate the kill-ring). C-_ will undo your last changes incrementally (M-r undoes all changes made to the current line). C-t and M-t will transpose the character (word) preceding the cursor and the one under the cursor.

Keeping the M- key down while you enter an integer (a minus sign meaning reverse behavior) gives an argument to your next readline command (for instance M-- C-k will kill text back to the start of line). If you prefer Vi--style editing, M-C-j will toggle you to Vi mode.

Of course you can change all these default bindings. For that you need to create a file named .inputrc in your home directory. For instance (notice the embedding conditional in case you would want specific bindings for gp):

  C<if Pari-GP
    set show-all-if-ambiguous
    "\C-h": backward-delete-char
    "\e\C-h": backward-kill-word
    "\C-xd": dump-functions
    (: "\C-v()\C-b"       # can be annoying when copy-pasting!
    [: "\C-v[]\C-b"
  >endif

@3C-x C-r will re-read this init file, incorporating any changes made to it during the current session.

@3Note. By default, ( and [ are bound to the function pari-matched-insert which, if ``electric parentheses'' are enabled (default: off) will automatically insert the matching closure (respectively ) and ]). This behavior can be toggled on and off by giving the numeric argument -2 to ( (M--2(), which is useful if you want, e.g to copy-paste some text into the calculator. If you do not want a toggle, you can use M--0 / M--1 to specifically switch it on or off).

@3Note. In some versions of readline (2.1 for instance), the Alt or Meta key can give funny results (output 8-bit accented characters for instance). If you do not want to fall back to the Esc combination, put the following two lines in your .inputrc:

    set convert-meta on
    set output-meta off

@3Command completion and online help. Hitting < TAB > will complete words for you. This mechanism is context-dependent: gp will strive to only give you meaningful completions in a given context (it will fail sometimes, but only under rare and restricted conditions).

For instance, shortly after a ~, we expect a user name, then a path to some file. Directly after default( has been typed, we would expect one of the default keywords. After whatnow( , we expect the name of an old function, which may well have disappeared from this version. After a '.', we expect a member keyword. And generally of course, we expect any GP symbol which may be found in the hashing lists: functions (both yours and GP's), and variables.

If, at any time, only one completion is meaningful, gp will provide it together with

@3* an ending comma if we are completing a default,

@3* a pair of parentheses if we are completing a function name. In that case hitting < TAB > again will provide the argument list as given by the online helpFOOTNOTE<<< recall that you can always undo the effect of the preceding keys by hitting C-_ >>>.

Otherwise, hitting < TAB > once more will give you the list of possible completions. Just experiment with this mechanism as often as possible, you will probably find it very convenient. For instance, you can obtain default(seriesprecision,10), just by hitting def < TAB > se < TAB > 10, which saves 18 keystrokes (out of 27).

Hitting M-h will give you the usual short online help concerning the word directly beneath the cursor, M-H will yield the extended help corresponding to the help default program (usually opens a dvi previewer, or runs a primitive tex-to-ASCII program). None of these disturb the line you were editing.


GNU Emacs and PariEmacs

If you install the PariEmacs package (see Appendix A), you may use gp as a subprocess in Emacs. You then need to include in your .emacs file the following lines:

    (autoload 'gp-mode "pari" nil t)
    (autoload 'gp-script-mode "pari" nil t)
    (autoload 'gp "pari" nil t)
    (autoload 'gpman "pari" nil t)
    (setq auto-mode-alist
      (cons '("\\.gp$" . gp-script-mode) auto-mode-alist))

@3which autoloads functions from the PariEmacs package and ensures that file with the .gp suffix are edited in gp-script mode.

Once this is done, under GNU Emacs if you type M-x gp (where as usual M is the Meta key), a special shell will be started launching gp with the default stack size and prime limit. You can then work as usual under gp, but with all the facilities of an advanced text editor. See the PariEmacs documentation for customizations, menus, etc.

\newpage