Porting pari/gp routines into a C program??

EdH · 2012-10-27, 22:08

This will be detailed; possibly much more information than is needed, but I'm getting frustrated (again).

Background: I'm trying to learn how to use some pari/gp routines within some c programs, via an IDE called Code::Blocks. Unfortunately, I don't understand some of the declarations for the inner workings of pari/gp. I didn't think I needed to, to construct a simple executable. Apparently, I was wrong (yet, again).

I have the following function, called cycle02.gp, that looks for aliquot cycles:

Code:

findcycle(a,b)=
{
  for(e=a,b,
    i=e;
    for(c=1,10,
      i=sigma(i)-i;
      if(i<e,
        break()
      );
      if(i==e,
        print(i,"  ",c);
        break()
      )
    );
  )
}

I can load this into gp and run it:

Code:

[user@comp pari]$ gp2c-run -g cycle02.gp
                  GP/PARI CALCULATOR Version 2.3.5 (released)
          amd64 running linux (x86-64/GMP-4.3.2 kernel) 64-bit version
      compiled: Feb  8 2011, gcc-4.6.0 20110205 (Red Hat 4.6.0-0.6) (GCC) 
    (readline v6.2 enabled [was v6.1 in Configure], extended help available)

                     Copyright (C) 2000-2006 The PARI Group

PARI/GP is free software, covered by the GNU General Public License, and 
comes WITHOUT ANY WARRANTY WHATSOEVER.

Type ? for help, \q to quit.
Type ?12 for how to get moral (and possibly technical) support.

parisize = 8000000, primelimit = 500000
? findcycle(10000,20000)
10744  2
12285  2
12496  5
17296  2
?

gp2c creates the following file cycle02.gp.c:

Code:

/*-*- compile-command: "/usr/bin/gcc -c -o cycle02.gp.o -O3 -Wall -fno-strict-aliasing -fomit-frame-pointer -fPIC -I/usr/include cycle02.gp.c && /usr/bin/gcc -o cycle02.gp.so -shared -O3 -Wall -fno-strict-aliasing -fomit-frame-pointer -fPIC -Wl,-shared cycle02.gp.o -lc -ldl -lm -L/usr/lib64 -lgmp -L/usr/lib64 -lpari"; -*-*/
#include <pari/pari.h>
/*
GP;install("init_cycle02","v","init_cycle02","./cycle02.gp.so");
GP;install("findcycle","vD0,G,D0,G,","findcycle","./cycle02.gp.so");
*/
void init_cycle02(void);
void findcycle(GEN a, GEN b);
/*End of prototype*/

void
init_cycle02(void)      /* void */
{
  pari_sp ltop = avma;
  avma = ltop;
  return;
}

void
findcycle(GEN a, GEN b)      /* void */
{
  pari_sp ltop = avma;
  GEN i = pol_x[fetch_user_var("i")];
  {
    pari_sp btop = avma, st_lim = stack_lim(btop, 1);
    GEN e = gen_0;
    for (e = gcopy(a); gcmp(e, b) <= 0; e = gaddgs(e, 1))
    {
      i = e;
      {
        pari_sp btop = avma, st_lim = stack_lim(btop, 1);
        long c;
        for (c = 1; c <= 10; ++c)
        {
          i = gsub(gsumdiv(i), i);
          if (gcmp(i, e) < 0)
            break;
          if (gequal(i, e))
          {
            pariprintf("%Z  %ld\n", i, c);
            break;
          }
          if (low_stack(st_lim, stack_lim(btop, 1)))
            i = gerepilecopy(btop, i);
        }
      }
      if (low_stack(st_lim, stack_lim(btop, 1)))
        gerepileall(btop, 2, &e, &i);
    }
  }
  avma = ltop;
  return;
}

I took the above routines and tried to build a basic program in c using Code::Blocks:

Code:

#include <stdio.h>
#include <stdlib.h>
#include "pari.h"

void init_cycle02(void);
void findcycle(GEN a, GEN b);


int main(int argc, char *argv[])
{
    GEN j, k;

    j=(GEN)100000;
    k=(GEN)200000;

    findcycle(j, k);
    return 0;
}

void
init_cycle02(void)      /* void */
{
  pari_sp ltop = avma;
  avma = ltop;
  return;
}

void
findcycle(GEN a, GEN b)      /* void */
{
  pari_sp ltop = avma;
  GEN i = pol_x[fetch_user_var("i")];
  {
    pari_sp btop = avma, st_lim = stack_lim(btop, 1);
    GEN e = gen_0;
    for (e = gcopy(a); gcmp(e, b) <= 0; e = gaddgs(e, 1))
    {
      i = e;
      {
        pari_sp btop = avma, st_lim = stack_lim(btop, 1);
        long c;
        for (c = 1; c <= 10; ++c)
        {
          i = gsub(gsumdiv(i), i);
          if (gcmp(i, e) < 0)
            break;
          if (gequal(i, e))
          {
            pariprintf("%Z  %ld\n", i, c);
            break;
          }
          if (low_stack(st_lim, stack_lim(btop, 1)))
            i = gerepilecopy(btop, i);
        }
      }
      if (low_stack(st_lim, stack_lim(btop, 1)))
        gerepileall(btop, 2, &e, &i);
    }
  }
  avma = ltop;
  return;
}

and, then I compiled and ran the program:

Code:

Process returned -1 (0xFFFFFFFF)   execution time: 0.004 s
Press ENTER to continue.

So, then I ran the debugger:

Code:

Program received signal SIGSEGV, Segmentation fault.
In fetch_named_var () (usr/lib64/libpari-gmp.so.2)
Debugger finished with status 0

I am obviously missing something...

Would anyone like to point me in a direction where I can find the missing details on how to use pari/gp functions within c programs?

Thanks!

CRGreathouse · 2012-10-28, 01:27

Presumably it's trying to access the global variable i, but there are no global variables since you're not in gp. Try changing

i=e;

into

my(i=e);

in the GP script to make the variable local, then see if it works. Also don't forget to use the -g option when running gp2c.

EdH · 2012-10-28, 03:00

Quote:

Originally Posted by CRGreathouse

Presumably it's trying to access the global variable i, but there are no global variables since you're not in gp. Try changing

i=e;

into

my(i=e);

in the GP script to make the variable local, then see if it works. Also don't forget to use the -g option when running gp2c.

Thanks for the help, but I got no detectable change. I do normally run the -g option, but simply skipped it for the small sample. I did use it this time.

I should be able to do what I'm try to, shouldn't I?

Maybe I'll try something more simple...

If/when this does work, will the c version have a noticeable difference in speed over running it via gp2c-run?

EdH · 2012-10-28, 03:30

I tried again with a version that's about as simple as I can get. I used the following gp routine:

Code:

pTest(a,b)=
{
  for(e=a,b,
    print("e is ",e)
  )
}

Terminal result:

Code:

[user@comp pari]$ gp2c-run -g pariTest.gp
                  GP/PARI CALCULATOR Version 2.3.5 (released)

...

Type ? for help, \q to quit.
Type ?12 for how to get moral (and possibly technical) support.

parisize = 8000000, primelimit = 500000
? pTest(15,32)
e is 15
e is 16
...
e is 31
e is 32
?

The PariTest.gp.c file:

Code:

/*-*- compile-command: "/usr/bin/gcc -c -o pariTest.gp.o -O3 -Wall -fno-strict-aliasing -fomit-frame-pointer -fPIC -I/usr/include pariTest.gp.c && /usr/bin/gcc -o pariTest.gp.so -shared -O3 -Wall -fno-strict-aliasing -fomit-frame-pointer -fPIC -Wl,-shared pariTest.gp.o -lc -ldl -lm -L/usr/lib64 -lgmp -L/usr/lib64 -lpari"; -*-*/
#include <pari/pari.h>
/*
GP;install("init_pariTest","v","init_pariTest","./pariTest.gp.so");
GP;install("pTest","vD0,G,D0,G,","pTest","./pariTest.gp.so");
*/
void init_pariTest(void);
void pTest(GEN a, GEN b);
/*End of prototype*/

void
init_pariTest(void)      /* void */
{
  pari_sp ltop = avma;
  avma = ltop;
  return;
}

void
pTest(GEN a, GEN b)      /* void */
{
  pari_sp ltop = avma;
  {
    pari_sp btop = avma, st_lim = stack_lim(btop, 1);
    GEN e = gen_0;
    for (e = gcopy(a); gcmp(e, b) <= 0; e = gaddgs(e, 1))
    {
      pariprintf("e is %Z\n", e);
      if (low_stack(st_lim, stack_lim(btop, 1)))
        e = gerepilecopy(btop, e);
    }
  }
  avma = ltop;
  return;
}

Debugger message:

Code:

Program received signal SIGSEGV, Segmentation fault.
In gcopy () (usr/lib64/libpari-gmp.so.2)

Thanks for any assistance...

jcrombie · 2012-10-28, 05:06

EdH, One problem I see with your program is that there is no pari_init() at the beginning. You need something like:

Code:

  pari_init( 1e9, 1e8 );

before you execute any pari code.

It's also good form to call:

Code:

  pari_close();

at the end.

Another problem is the casts:

Code:

 j=(GEN)100000;
 k=(GEN)200000;

You are assigning integers into pointers here. I really don't know much about pari, but the proper coding techniques are to be found in the "libpari.dvi" documentation.

Cheers

CRGreathouse · 2012-10-28, 05:13

Quote:

Originally Posted by jcrombie

Another problem is the casts:

Code:

 j=(GEN)100000;
 k=(GEN)200000;

You are assigning integers into pointers here. I really don't know much about pari, but the proper coding techniques are to be found in the "libpari.dvi" documentation.

Cheers

Ah, good catch -- I was just looking at the code that gp2c generated. Those should be

Code:

j = utoipos(100000);
k = utoipos(200000);

EdH · 2012-10-28, 12:54

Thank you both for the help, but it still resists my efforts:

Code:

Program received signal SIGSEGV, Segmentation fault.
In pariputc () (/usr/lib64/libpari-gmp.so.2)

I had wondered what I was supposed to do with the created init_pariTest function, also. I tried simply inserting it at the beginning of main, but it had no effect, either. Does it belong somewhere?

Off to research libpari.dvi docs...

CRGreathouse · 2012-10-28, 17:22

Would you post the code you're using now, with the various modifications?

EdH · 2012-10-28, 18:31

This is my latest very basic testing setup from start to finish.

pariTest.gp:

Code:

pTest(a,b)=
{
  for(e=a,b,
    print("e is ",e)
  )
}

terminal session:

Code:

[user@comp pari]$ gp2c-run -g pariTest.gp
                  GP/PARI CALCULATOR Version 2.3.5 (released)
          amd64 running linux (x86-64/GMP-4.3.2 kernel) 64-bit version
      compiled: Feb  8 2011, gcc-4.6.0 20110205 (Red Hat 4.6.0-0.6) (GCC) 
    (readline v6.2 enabled [was v6.1 in Configure], extended help available)

                     Copyright (C) 2000-2006 The PARI Group

PARI/GP is free software, covered by the GNU General Public License, and 
comes WITHOUT ANY WARRANTY WHATSOEVER.

Type ? for help, \q to quit.
Type ?12 for how to get moral (and possibly technical) support.

parisize = 8000000, primelimit = 500000
? pTest(11,17)
e is 11
e is 12
e is 13
e is 14
e is 15
e is 16
e is 17
?

generated pariTest.gp.c

Code:

/*-*- compile-command: "/usr/bin/gcc -c -o pariTest.gp.o -O3 -Wall -fno-strict-aliasing -fomit-frame-pointer -fPIC -I/usr/include pariTest.gp.c && /usr/bin/gcc -o pariTest.gp.so -shared -O3 -Wall -fno-strict-aliasing -fomit-frame-pointer -fPIC -Wl,-shared pariTest.gp.o -lc -ldl -lm -L/usr/lib64 -lgmp -L/usr/lib64 -lpari"; -*-*/
#include <pari/pari.h>
/*
GP;install("init_pariTest","v","init_pariTest","./pariTest.gp.so");
GP;install("pTest","vD0,G,D0,G,","pTest","./pariTest.gp.so");
*/
void init_pariTest(void);
void pTest(GEN a, GEN b);
/*End of prototype*/

void
init_pariTest(void)      /* void */
{
  pari_sp ltop = avma;
  avma = ltop;
  return;
}

void
pTest(GEN a, GEN b)      /* void */
{
  pari_sp ltop = avma;
  {
    pari_sp btop = avma, st_lim = stack_lim(btop, 1);
    GEN e = gen_0;
    for (e = gcopy(a); gcmp(e, b) <= 0; e = gaddgs(e, 1))
    {
      pariprintf("e is %Z\n", e);
      if (low_stack(st_lim, stack_lim(btop, 1)))
        e = gerepilecopy(btop, e);
    }
  }
  avma = ltop;
  return;
}

c program:

Code:

#include <stdio.h>
#include <stdlib.h>
#include "pari.h"

void init_pariTest(void);
void pTest(GEN a, GEN b);


int main()
{
    GEN c, d;

    c=utoipos(11);
    d=utoipos(17);

    pari_init(1e9,1e8);
    init_pariTest();

    pTest(c,d);

    pari_close();

    printf("\nFinished!\n");
    return 0;
}

void
init_pariTest(void)      /* void */
{
  pari_sp ltop = avma;
  avma = ltop;
  return;
}

void
pTest(GEN a, GEN b)      /* void */
{
  pari_sp ltop = avma;
  {
    pari_sp btop = avma, st_lim = stack_lim(btop, 1);
    GEN e = gen_0;
    for (e = gcopy(a); gcmp(e, b) <= 0; e = gaddgs(e, 1))
    {
      pariprintf("e is %Z\n", e);
      if (low_stack(st_lim, stack_lim(btop, 1)))
        e = gerepilecopy(btop, e);
    }
  }
  avma = ltop;
  return;
}

c program session output:

Code:

Process returned -1 (0xFFFFFFFF)   execution time : 0.004 s
Press ENTER to continue.

debugger output:

Code:

Building to ensure sources are up-to-date
Build succeeded
Selecting target: 
Release
Adding source dir: /home/user/Programming/pariTest/
Adding source dir: /home/user/Programming/pariTest/
Adding file: bin/Release/pariTest
Starting debugger: 
done
Registered new type: wxString
Registered new type: STL String
Registered new type: STL Vector
Setting breakpoints
(no debugging symbols found)...done.
Debugger name and version: GNU gdb (GDB) Fedora (7.3.50.20110722-16.fc16)
Program received signal SIGSEGV, Segmentation fault.
In pariputc () (/usr/lib64/libpari-gmp.so.2)

My program console gives me the following message:

Code:

warning: GDB: Failed to set controlling terminal: Operation not permitted

Call stack:

Code:

#0 0x36e27f53fd    pariputc() (/usr/lib64/libpari-gmp.so.2:??)
#1 0x36e27fe74a    pari_err() (/usr/lib64/libpari-gmp.so.2:??)
#2 (    0x0000000000400976 in ??() (??:??)
#3 (    0x000000000040099f in ??() (??:??)
#4 (    0x00000000004009db in ??() (??:??)
#5 (    0x0000000000400a15 in ??() (??:??)
#6 0x36e122169d    __libc_start_main() (/lib64/libc.so.6:??)
#7 (    0x0000000000400839 in ??() (??:??)
#8 (    0x00007fffffffe708 in ??() (??:??)
#9 (    0x000000000000001c in ??() (??:??)
#10 (    0x0000000000000001 in ??() (??:??)
#11 (    0x00007fffffffea03 in ??() (??:??)
#12 (    0x0000000000000000 in ??() (??:??)

Thanks for any help. I have the aforementioned documents and am looking into them. I haven't tried constructing this in the referenced Makefile manner, since I've never programmed with makefiles. I have started looking into them, though...

frmky · 2012-10-28, 22:08

I took the output of gp2c -g and added a short main function to produce

Code:

/*-*- compile-command: "cc -c -o ptest.gp.o -g -O3 -Wall -fomit-frame-pointer -fno-strict-aliasing -fPIC -I"/usr/include" ptest.gp.c && cc -o ptest.gp.so -shared -g -O3 -Wall -fomit-frame-pointer -fno-strict-aliasing -fPIC -Wl,-shared ptest.gp.o -lc -lm -L/usr/lib -lpari"; -*-*/
#include <pari/pari.h>
/*
GP;install("init_ptest","v","init_ptest","./ptest.gp.so");
GP;install("pTest","vD0,G,D0,G,","pTest","./ptest.gp.so");
*/
void init_ptest(void);
void pTest(GEN a, GEN b);
/*End of prototype*/

void
init_ptest(void)          /* void */
{
  pari_sp ltop = avma;
  avma = ltop;
  return;
}

void
pTest(GEN a, GEN b)       /* void */
{
  pari_sp ltop = avma;
  {
    pari_sp btop = avma, st_lim = stack_lim(btop, 1);
    GEN e = gen_0;
    for (e = gcopy(a); gcmp(e, b) <= 0; e = gaddgs(e, 1))
    {
      pari_printf("e is %Ps\n", e);
      if (low_stack(st_lim, stack_lim(btop, 1)))
        e = gerepilecopy(btop, e);
    }
  }
  avma = ltop;
  return;
}

int main() {
  long a, b;
  GEN x, y;
  pari_init(1000000,2);
  init_ptest();
  printf("a = ");
  scanf("%ld",&a);
  printf("b = ");
  scanf("%ld",&b);
  x = stoi(a);
  y = stoi(b);
  pTest(x, y);
  pari_close();
  return 0;
}

This I compiled with

Code:

gcc -o ptest -g -O3 -Wall -fomit-frame-pointer -fno-strict-aliasing -fPIC -I"/usr/include" ptest.gp.c -lc -lm -L/usr/lib -lpari

This runs as expected.

frmky · 2012-10-28, 22:15

Quote:

Originally Posted by EdH

If/when this does work, will the c version have a noticeable difference in speed over running it via gp2c-run?

I doubt it. The idea of gp2c-run is to combine the speed of the c compiled function with the ease of use of the gp calculator.

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
Question regarding basic routines being used in Yafu.	storflyt32	YAFU	2	2015-06-29 23:25
Porting my factorization applet to Android	alpertron	Programming	2	2013-03-19 11:28
PARI vs GAP	skan	Miscellaneous Math	0	2012-12-16 00:13
porting Prime95 to GTK+?	ixfd64	Software	13	2012-02-21 07:36
Possibly stupid question about porting games to Linux.	jasong	Linux	4	2006-12-23 21:24

2012-10-28, 01:27	#2
CRGreathouse Aug 2006 1011101100011₂ Posts	Presumably it's trying to access the global variable i, but there are no global variables since you're not in gp. Try changing i=e; into my(i=e); in the GP script to make the variable local, then see if it works. Also don't forget to use the -g option when running gp2c.

2012-10-28, 05:06	#5
jcrombie "Jonathan" Jul 2010 In a tangled web... 5·43 Posts	EdH, One problem I see with your program is that there is no pari_init() at the beginning. You need something like: Code: pari_init( 1e9, 1e8 ); before you execute any pari code. It's also good form to call: Code: pari_close(); at the end. Another problem is the casts: Code: j=(GEN)100000; k=(GEN)200000; You are assigning integers into pointers here. I really don't know much about pari, but the proper coding techniques are to be found in the "libpari.dvi" documentation. Cheers Last fiddled with by jcrombie on 2012-10-28 at 05:11

2012-10-28, 12:54	#7
EdH "Ed Hall" Dec 2009 Adirondack Mtns 5224₁₀ Posts	Thank you both for the help, but it still resists my efforts: Code: Program received signal SIGSEGV, Segmentation fault. In pariputc () (/usr/lib64/libpari-gmp.so.2) I had wondered what I was supposed to do with the created init_pariTest function, also. I tried simply inserting it at the beginning of main, but it had no effect, either. Does it belong somewhere? Off to research libpari.dvi docs...

2012-10-28, 17:22	#8
CRGreathouse Aug 2006 5987₁₀ Posts	Would you post the code you're using now, with the various modifications?