mersenneforum.org

mersenneforum.org (https://www.mersenneforum.org/index.php)
-   Programming (https://www.mersenneforum.org/forumdisplay.php?f=29)
-   -   Porting pari/gp routines into a C program?? (https://www.mersenneforum.org/showthread.php?t=17343)

EdH 2012-10-27 22:08

Porting pari/gp routines into a C program??
 
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()
)
);
)
}
[/code]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
?
[/code]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;
}
[/code]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;
}
[/code]and, then I compiled and ran the program:
[code]
Process returned -1 (0xFFFFFFFF) execution time: 0.004 s
Press ENTER to continue.
[/code]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
[/code]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=CRGreathouse;316232]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.[/QUOTE]
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)
)
}
[/code]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
?
[/code]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;
}
[/code]Debugger message:
[code]
Program received signal SIGSEGV, Segmentation fault.
In gcopy () (usr/lib64/libpari-gmp.so.2)
[/code]
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 );[/CODE]before you execute any pari code.

It's also good form to call:

[CODE] pari_close();
[/CODE]at the end.

Another problem is the casts:

[CODE] j=(GEN)100000;
k=(GEN)200000;
[/CODE]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=jcrombie;316254]Another problem is the casts:

[CODE] j=(GEN)100000;
k=(GEN)200000;
[/CODE]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[/QUOTE]


Ah, good catch -- I was just looking at the code that gp2c generated. Those should be
[code]j = utoipos(100000);
k = utoipos(200000);[/code]

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)
[/code]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)
)
}
[/code]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
?
[/code]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;
}
[/code]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;
}
[/code]c program session output:
[code]
Process returned -1 (0xFFFFFFFF) execution time : 0.004 s
Press ENTER to continue.
[/code]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)
[/code]My program console gives me the following message:
[code]
warning: GDB: Failed to set controlling terminal: Operation not permitted
[/code]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 ??() (??:??)
[/code]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;
}
[/CODE]

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
[/CODE]

This runs as expected.

frmky 2012-10-28 22:15

[QUOTE=EdH;316242]If/when this does work, will the c version have a noticeable difference in speed over running it via gp2c-run?[/QUOTE]
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.


All times are UTC. The time now is 13:50.

Powered by vBulletin® Version 3.8.11
Copyright ©2000 - 2021, Jelsoft Enterprises Ltd.