Incorporating Pari into a C program

Hi Pari Users

I do a lot of things with graph polynomials and their roots and for that I use GP a lot, but only in a naïve way – essentially just calling polrootsreal on billions of polynomials (using a Perl script to actually format the relevant GP commands).

Now I need to do something different and although I have tried to read the manual, I am lost.

What I need is a C function that will operate as a “plugin” to an existing C program (not written by me).

The C program traverses a huge search tree, and at every node it calls the plugin function (whose name is given to the outer program at compile time) and whose arguments are fixed. So I write a routine, something like

int decision(int *a, int n)

which is called by main program at every node (and there will be billions of nodes).

The response given by my plugin determines whether that branch of the search tree should be extended or pruned.

What does my function do?

It uses the input array to construct an n x n real symmetric matrix, and tests whether its largest eigenvalue is greater than 2, and then returns 1 or 0 accordingly.

The matrices are of moderate-ish size (20x20 – 30x30) and the entries are small rational numbers.

At the moment, (just in the testing phase) I create the matrix using floating point arithmetic and use a simple eigenvalue program that had made available on some academic’s home page.

Now I want to do it rigorously – that is, I need to make sure that when my routine says “the largest eigenvalue is greater than 2” it is definitely correct. I am willing to have a routine that makes an occasional mistake in the other direction – it can incorrectly report that the eigenvalue is less than 2, but only very occasionally.

(I also want to do it faster (of course) and the obvious place to gain speed is to only calculate the largest eigenvalue.)

There are various ways that I could try using Pari, but I cannot understand the documentation about using Pari from a C program – in particular, the instructions about initializing the Pari stack and allocating and freeing memory are beyond me.

Can anyone provide a simple example of how I could use a Pari routine from within my plugin function? And how I should manage the initialization step… the plugin mechanism does allow me define a second function that is called just once at the beginning of the search, and presumably I would put the pari-initialize stuff in there?

Unfortunately, my C is sufficiently rusty that I am not sure how the plugin can access things created by the “initialize plugin” function.

I have also considered that perhaps Pari is overkill for what I need, and that it would be better to use, say GMP, to create my rational matrix, then simply run a few iterations of the power method starting with some vector v0, so v1=Bv0, v2=Bv1, v3=Bv2, … and then check if the Rayleigh quotient v3^TBv^3/v3^Tv is greater than 2.

Apologies for the long post, but any thoughts and comments welcome…

Gordon