0>[H[JPARI/GP Workshop      IHP    2004/09/14      Slide 1/12
How can I compile and run my scripts?
$ gp2c-run pari/examples/squfof.gp
         GP/PARI CALCULATOR Version 2.2.10 (development CHANGES-1.1117)
                i686 running linux (ix86 kernel) 32-bit version
         compiled: Mar  8 2005, gcc-2.95.4 20011002 (Debian prerelease)
               (readline v4.2a enabled, extended help available)

                       Copyright (C) 2003 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 = 4000000, primelimit = 500000
? squfof(3097180303181)
[419]
i = 596
Qfb(133225, 1719841, -261451, 0.E-28)
[0;35m%1 = [0;36m1691693[0m
? squfof(122294051504814979)
[20137]
? allocatemem()
? squfof(122294051504814979)
[20137]
[20137, 3445]
i = 46474
Qfb(321233929, 131349818, -367273962, 0.E-28)
[0;35m%2 = [0;36m73823023[0m
[0mGoodbye!



1>[H[JPARI/GP Workshop      IHP    2004/09/14      Slide 2/12
What happen behind the scene ?
$ ls pari/examples/squfof.gp*
pari/examples/squfof.gp
pari/examples/squfof.gp.c
pari/examples/squfof.gp.o
pari/examples/squfof.gp.so
/*-*- compile-command: "/usr/bin/gcc -c -o pari/examples/squfof.gp.o -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -I/home/bill/pari/i686/include pari/examples/squfof.gp.c && /usr/bin/gcc -o pari/examples/squfof.gp.so -shared -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -Wl,-shared pari/examples/squfof.gp.o -lc -ldl -lm -L/home/bill/pari/i686/lib -lpari"; -*-*/
#include <pari/pari.h>
/*
GP;install("squfof","D0,G,p","squfof","./pari/examples/squfof.gp.so");
GP;install("init_squfof","v","init_squfof","./pari/examples/squfof.gp.so");
*/
GEN squfof(GEN n, long prec);
void init_squfof(void);
/*End of prototype*/

void
init_squfof(void)
{
  return;
}

/* return one non-trivial divisor of n > 1 using Shanks's SQUFOF */
GEN
squfof(GEN n, long prec)
{
  GEN p1;
  GEN p = polx[fetch_user_var("p")], D = polx[fetch_user_var("D")], d = polx[fetch_user_var("d")], b = polx[fetch_user_var("b")], f = polx[fetch_user_var("f")], l = polx[fetch_user_var("l")], q = polx[fetch_user_var("q")], lq = polx[fetch_user_var("lq")], i = polx[fetch_user_var("i")], a = polx[fetch_user_var("a")], as = polx[fetch_user_var("as")], j = polx[fetch_user_var("j")], bb = polx[fetch_user_var("bb")], gs = polx[fetch_user_var("gs")];
  GEN p2;	  /* vec */
  GEN g = polx[fetch_user_var("g")], b1 = polx[fetch_user_var("b1")], p3;
  if (!gcmp0(gisprime(n, 0)) || !gcmp0(gcarrecomplet(n, &n)))
  {
    GEN p4;
    p4 = n;
    return p4;
  }
  p = gcoeff(gboundfact(n, 0), 1, 1);
  if (!gequal(p, n))
  {
    GEN p5;
    p5 = p;
    return p5;
  }
  if (gcmp1(gmodgs(n, 4)))
  {
    D = n;
    d = racine(D);
    b = mpadd(gshift(gdiventgs(gsubgs(d, 1), 2), 1), gun);
  }
  else
  {
    D = gshift(n, 2);
    d = racine(D);
    b = gshift(gdiventgs(d, 2), 1);
  }
  f = Qfb0(gun, b, gshift(gsub(gsqr(b), D), -2), NULL, prec);
  l = racine(d);
  q = cgetg(1,t_VEC);
  lq = gzero;
  i = gzero;
  while (1)
  {
    i = gaddgs(i, 1);
    f = qfbred0(f, 3, D, d, NULL);
    a = compo(f, 1);
    if ((gcmp0(gmodgs(i, 2))) && !gcmp0(gcarrecomplet(a, &as)))
    {
      j = gun;
      while (gcmp(j, lq) <= 0)
      {
        if (gequal(as, gel(q, itos(j))))
          break;
        j = gaddgs(j, 1);
      }
      if (gcmp(j, lq) > 0)
        break;
    }
    if (gcmp(gabs(a, prec), l) <= 0)
    {
      q = concat(q, gabs(a, prec));
      pariputsf("%Z\n",q);
      lq = gaddgs(lq, 1);
    }
  }
  pariputsf("i = %Z\n",i);
  pariputsf("%Z\n",f);
  bb = compo(f, 2);
  p2 = cgetg(4, t_VEC);
  gel(p2, 1) = as;
  gel(p2, 2) = bb;
  gel(p2, 3) = D;
  gs = gcd0(p2, NULL, 0);
  if (gcmpgs(gs, 1) > 0)
  {
    GEN p6;
    p6 = gs;
    return p6;
  }
  f = Qfb0(as, gneg(bb), gmul(as, compo(f, 3)), NULL, prec);
  g = qfbred0(f, 3, D, d, NULL);
  b = compo(g, 2);
  do
  {
    b1 = b;
    g = qfbred0(g, 3, D, d, NULL);
    b = compo(g, 2);
  } while(!gequal(b1, b));
  a = gabs(compo(g, 1), prec);
  if (!gcmp0(gmodgs(a, 2)))
    p3 = a;
  else
    p3 = gshift(a, -1);
  p1 = p3;
  return p1;
}




2>[H[JPARI/GP Workshop      IHP    2004/09/14      Slide 3/12
Using -g for garbage collecting
$ gp2c-run -g pari/examples/squfof.gp
         GP/PARI CALCULATOR Version 2.2.10 (development CHANGES-1.1117)
                i686 running linux (ix86 kernel) 32-bit version
         compiled: Mar  8 2005, gcc-2.95.4 20011002 (Debian prerelease)
               (readline v4.2a enabled, extended help available)

                       Copyright (C) 2003 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 = 4000000, primelimit = 500000
? squfof(122294051504814979)
[20137]
[20137, 3445]
i = 46474
Qfb(321233929, 131349818, -367273962, 0.E-28)
[0;35m%1 = [0;36m73823023[0m
? \rpari/examples/squfof.gp
? squfof(n)=if(isprime(n)||issquare(n,&n),return(n));p=factor(n,0)[1,1];if(p!=n,return(p));if(n%4==1,D=n;d=sqrtint(D);b=(((d-1)\2)<<1)+1,D=n<<2;d=sqrtint(D);b=(d\2)<<1);f=Qfb(1,b,(b^2-D)>>2);l=sqrtint(d);q=[];lq=0;i=0;while(1,i++;f=qfbred(f,3,D,d);a=component(f,1);if(!(i%2)&&issquare(a,&as),j=1;while(j<=lq,if(as==q[j],break);j++);if(j>lq,break));if(abs(a)<=l,q=concat(q,abs(a));print(q);lq++;));print("i = ",i);print(f);bb=component(f,2);gs=gcd([as,bb,D]);if(gs>1,return(gs));f=Qfb(as,-bb,as*component(f,3));g=qfbred(f,3,D,d);b=component(g,2);until(b1==b,b1=b;g=qfbred(g,3,D,d);b=component(g,2););a=abs(component(g,1));if(a%2,a,a>>1);
[0mGoodbye!
/*-*- compile-command: "/usr/bin/gcc -c -o pari/examples/squfof.gp.o -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -I/home/bill/pari/i686/include pari/examples/squfof.gp.c && /usr/bin/gcc -o pari/examples/squfof.gp.so -shared -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -Wl,-shared pari/examples/squfof.gp.o -lc -ldl -lm -L/home/bill/pari/i686/lib -lpari"; -*-*/
#include <pari/pari.h>
/*
GP;install("squfof","D0,G,p","squfof","./pari/examples/squfof.gp.so");
GP;install("init_squfof","v","init_squfof","./pari/examples/squfof.gp.so");
*/
GEN squfof(GEN n, long prec);
void init_squfof(void);
/*End of prototype*/

void
init_squfof(void)
{
  pari_sp ltop = avma;
  avma = ltop;
  return;
}

/* return one non-trivial divisor of n > 1 using Shanks's SQUFOF */
GEN
squfof(GEN n, long prec)
{
  pari_sp ltop = avma;
  GEN p1 = gzero;
  GEN p = polx[fetch_user_var("p")], D = polx[fetch_user_var("D")], d = polx[fetch_user_var("d")], b = polx[fetch_user_var("b")], f = polx[fetch_user_var("f")], l = polx[fetch_user_var("l")], q = polx[fetch_user_var("q")], lq = polx[fetch_user_var("lq")], i = polx[fetch_user_var("i")], a = polx[fetch_user_var("a")], as = polx[fetch_user_var("as")], j = polx[fetch_user_var("j")], bb = polx[fetch_user_var("bb")], gs = polx[fetch_user_var("gs")];
  GEN p2 = gzero;	  /* vec */
  GEN g = polx[fetch_user_var("g")], b1 = polx[fetch_user_var("b1")], p3 = gzero;
  if (!gcmp0(gisprime(n, 0)) || !gcmp0(gcarrecomplet(n, &n)))
  {
    GEN p4 = gzero;
    p4 = gcopy(n);
    p4 = gerepileupto(ltop, p4);
    return p4;
  }
  p = gcopy(gcoeff(gboundfact(n, 0), 1, 1));
  if (!gequal(p, n))
  {
    GEN p5 = gzero;
    p5 = p;
    p5 = gerepileupto(ltop, p5);
    return p5;
  }
  if (gcmp1(gmodgs(n, 4)))
  {
    D = gcopy(n);
    d = racine(D);
    b = mpadd(gshift(gdiventgs(gsubgs(d, 1), 2), 1), gun);
  }
  else
  {
    D = gshift(n, 2);
    d = racine(D);
    b = gshift(gdiventgs(d, 2), 1);
  }
  f = Qfb0(gun, b, gshift(gsub(gsqr(b), D), -2), NULL, prec);
  l = racine(d);
  q = cgetg(1,t_VEC);
  lq = gzero;
  i = gzero;
  {
    pari_sp btop = avma, st_lim = stack_lim(btop, 1);
    GEN *bptr[] = { &i, &f, &a, &as, &j, &q, &lq };
    while (1)
    {
      i = gaddgs(i, 1);
      f = qfbred0(f, 3, D, d, NULL);
      a = compo(f, 1);
      if ((gcmp0(gmodgs(i, 2))) && !gcmp0(gcarrecomplet(a, &as)))
      {
        j = gun;
        {
          pari_sp btop = avma;
          while (gcmp(j, lq) <= 0)
          {
            if (gequal(as, gel(q, itos(j))))
              break;
            j = gaddgs(j, 1);
            j = gerepileupto(btop, j);
          }
        }
        if (gcmp(j, lq) > 0)
          break;
      }
      if (gcmp(gabs(a, prec), l) <= 0)
      {
        q = concat(q, gabs(a, prec));
        pariputsf("%Z\n",q);
        lq = gaddgs(lq, 1);
      }
      if (low_stack(st_lim, stack_lim(btop, 1)))
        gerepilemany(btop, bptr, 7);
    }
  }
  pariputsf("i = %Z\n",i);
  pariputsf("%Z\n",f);
  bb = compo(f, 2);
  p2 = cgetg(4, t_VEC);
  gel(p2, 1) = as;
  gel(p2, 2) = bb;
  gel(p2, 3) = D;
  gs = gcd0(p2, NULL, 0);
  if (gcmpgs(gs, 1) > 0)
  {
    GEN p6 = gzero;
    p6 = gs;
    p6 = gerepileupto(ltop, p6);
    return p6;
  }
  f = Qfb0(as, gneg(bb), gmul(as, compo(f, 3)), NULL, prec);
  g = qfbred0(f, 3, D, d, NULL);
  b = compo(g, 2);
  {
    pari_sp btop = avma, st_lim = stack_lim(btop, 1);
    GEN *bptr[] = { &b1, &g, &b };
    do
    {
      b1 = b;
      g = qfbred0(g, 3, D, d, NULL);
      b = compo(g, 2);
      if (low_stack(st_lim, stack_lim(btop, 1)))
        gerepilemany(btop, bptr, 3);
    } while(!gequal(b1, b));
  }
  a = gabs(compo(g, 1), prec);
  if (!gcmp0(gmodgs(a, 2)))
    p3 = a;
  else
    p3 = gshift(a, -1);
  p1 = p3;
  p1 = gerepileupto(ltop, p1);
  return p1;
}




3>[H[JPARI/GP Workshop      IHP    2004/09/14      Slide 4/12
$ gp2c-run -g -s_c pari/examples/squfof.gp
         GP/PARI CALCULATOR Version 2.2.10 (development CHANGES-1.1117)
                i686 running linux (ix86 kernel) 32-bit version
         compiled: Mar  8 2005, gcc-2.95.4 20011002 (Debian prerelease)
               (readline v4.2a enabled, extended help available)

                       Copyright (C) 2003 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 = 4000000, primelimit = 500000
? \rpari/examples/squfof.gp
? squfof(n)=if(isprime(n)||issquare(n,&n),return(n));p=factor(n,0)[1,1];if(p!=n,return(p));if(n%4==1,D=n;d=sqrtint(D);b=(((d-1)\2)<<1)+1,D=n<<2;d=sqrtint(D);b=(d\2)<<1);f=Qfb(1,b,(b^2-D)>>2);l=sqrtint(d);q=[];lq=0;i=0;while(1,i++;f=qfbred(f,3,D,d);a=component(f,1);if(!(i%2)&&issquare(a,&as),j=1;while(j<=lq,if(as==q[j],break);j++);if(j>lq,break));if(abs(a)<=l,q=concat(q,abs(a));print(q);lq++;));print("i = ",i);print(f);bb=component(f,2);gs=gcd([as,bb,D]);if(gs>1,return(gs));f=Qfb(as,-bb,as*component(f,3));g=qfbred(f,3,D,d);b=component(g,2);until(b1==b,b1=b;g=qfbred(g,3,D,d);b=component(g,2););a=abs(component(g,1));if(a%2,a,a>>1);
? #
   timer = 1 (on)
? squfof_c(122294051504814979)
[20137]
[20137, 3445]
i = 46474
Qfb(321233929, 131349818, -367273962, 0.E-28)
time = [0;33m140 ms[0m.
[0;35m%1 = [0;36m73823023[0m
? squfof(122294051504814979)
[20137]
[20137, 3445]
i = 46474
Qfb(321233929, 131349818, -367273962, 0.E-28)
time = [0;33m860 ms[0m.
[0;35m%2 = [0;36m73823023[0m
[0mGoodbye!



4>[H[JPARI/GP Workshop      IHP    2004/09/14      Slide 5/12
How can I compile directly with gp2c ?
$ gp2c pari/examples/squfof.gp > squfof.gp.c
$ head -n 6 squfof.gp.c
/*-*- compile-command: "/usr/bin/gcc -c -o pari/examples/squfof.gp.o -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -I/home/bill/pari/i686/include pari/examples/squfof.gp.c && /usr/bin/gcc -o pari/examples/squfof.gp.so -shared -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -Wl,-shared pari/examples/squfof.gp.o -lc -ldl -lm -L/home/bill/pari/i686/lib -lpari"; -*-*/
#include <pari/pari.h>
/*
GP;install("squfof","D0,G,p","squfof","./pari/examples/squfof.gp.so");
GP;install("init_squfof","v","init_squfof","./pari/examples/squfof.gp.so");
*/



5>[H[JPARI/GP Workshop      IHP    2004/09/14      Slide 6/12
$ gp2c -s_c pari/examples/squfof.gp > squfof2.gp.c
$ head -n 6 squfof2.gp.c
/*-*- compile-command: "/usr/bin/gcc -c -o pari/examples/squfof.gp.o -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -I/home/bill/pari/i686/include pari/examples/squfof.gp.c && /usr/bin/gcc -o pari/examples/squfof.gp.so -shared -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -Wl,-shared pari/examples/squfof.gp.o -lc -ldl -lm -L/home/bill/pari/i686/lib -lpari"; -*-*/
#include <pari/pari.h>
/*
GP;install("squfof","D0,G,p","squfof","./pari/examples/squfof.gp.so");
GP;install("init_squfof","v","init_squfof","./pari/examples/squfof.gp.so");
*/
--- squfof.gp.c	Tue Mar  8 19:05:08 2005
+++ squfof2.gp.c	Tue Mar  8 19:05:09 2005
@@ -1,8 +1,8 @@
 /*-*- compile-command: "/usr/bin/gcc -c -o pari/examples/squfof.gp.o -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -I/home/bill/pari/i686/include pari/examples/squfof.gp.c && /usr/bin/gcc -o pari/examples/squfof.gp.so -shared -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -Wl,-shared pari/examples/squfof.gp.o -lc -ldl -lm -L/home/bill/pari/i686/lib -lpari"; -*-*/
 #include <pari/pari.h>
 /*
-GP;install("squfof","D0,G,p","squfof","./pari/examples/squfof.gp.so");
-GP;install("init_squfof","v","init_squfof","./pari/examples/squfof.gp.so");
+GP;install("squfof","D0,G,p","squfof_c","./pari/examples/squfof.gp.so");
+GP;install("init_squfof","v","init_squfof_c","./pari/examples/squfof.gp.so");
 */
 GEN squfof(GEN n, long prec);
 void init_squfof(void);



6>[H[JPARI/GP Workshop      IHP    2004/09/14      Slide 7/12
Using GP2C to find errors in GP scripts
$ gp2c -W pari/examples/squfof.gp > squfof.gp.c



7>[H[JPARI/GP Workshop      IHP    2004/09/14      Slide 8/12
Using compiled functions in a new program
$ cat fact.gp
install("squfof","D0,G,p","squfof","./pari/examples/squfof.gp.so");
fact_mersenne(p)=squfof(2^p-1)
gp2c-run -pmy_ -g pari/examples/squfof.gp
         GP/PARI CALCULATOR Version 2.2.10 (development CHANGES-1.1117)
                i686 running linux (ix86 kernel) 32-bit version
         compiled: Mar  8 2005, gcc-2.95.4 20011002 (Debian prerelease)
               (readline v4.2a enabled, extended help available)

                       Copyright (C) 2003 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 = 4000000, primelimit = 500000
[0mGoodbye!
$ cat fact.gp
install("my_squfof","D0,G,p","squfof","./pari/examples/squfof.gp.so");
fact_mersenne(p)=squfof(2^p-1)
gp2c-run -g fact.gp
         GP/PARI CALCULATOR Version 2.2.10 (development CHANGES-1.1117)
                i686 running linux (ix86 kernel) 32-bit version
         compiled: Mar  8 2005, gcc-2.95.4 20011002 (Debian prerelease)
               (readline v4.2a enabled, extended help available)

                       Copyright (C) 2003 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 = 4000000, primelimit = 500000
? fact_mersenne(67)
i = 2418
Qfb(10825778209, 4021505768, -13258245519, 0.E-28)
[0;35m%1 = [0;36m193707721[0m
[0mGoodbye!



8>[H[JPARI/GP Workshop      IHP    2004/09/14      Slide 9/12
Optimizing GP2c output using types
lucas(p) =
{
  local(u,q); u=4; q=1<<p - 1;
  for(k=3,p, u = (sqr(u)-2)%q);
  u == 0
}
$ gp2c lucas.gp | less
/*-*- compile-command: "/usr/bin/gcc -c -o lucas.gp.o -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -I/home/bill/pari/i686/include lucas.gp.c && /usr/bin/gcc -o lucas.gp.so -shared -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -Wl,-shared lucas.gp.o -lc -ldl -lm -L/home/bill/pari/i686/lib -lpari"; -*-*/
#include <pari/pari.h>
/*
GP;install("lucas","lD0,G,","lucas","./lucas.gp.so");
GP;install("init_lucas","v","init_lucas","./lucas.gp.so");
*/
long lucas(GEN p);
void init_lucas(void);
/*End of prototype*/

void
init_lucas(void)
{
  return;
}

long
lucas(GEN p)
{
  long l1;	  /* bool */
  GEN u, q;
  u = stoi(4);
  q = subis(shifti(gun, itos(p)), 1);
  {
    GEN k;
    for (k = stoi(3); gcmp(k, p) <= 0; k = gaddgs(k, 1))
      u = gmod(gsubgs(gsqr(u), 2), q);
  }
  l1 = gcmp0(u);
  return l1;
}

$ gp2c -g lucas.gp | less
/*-*- compile-command: "/usr/bin/gcc -c -o lucas.gp.o -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -I/home/bill/pari/i686/include lucas.gp.c && /usr/bin/gcc -o lucas.gp.so -shared -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -Wl,-shared lucas.gp.o -lc -ldl -lm -L/home/bill/pari/i686/lib -lpari"; -*-*/
#include <pari/pari.h>
/*
GP;install("lucas","lD0,G,","lucas","./lucas.gp.so");
GP;install("init_lucas","v","init_lucas","./lucas.gp.so");
*/
long lucas(GEN p);
void init_lucas(void);
/*End of prototype*/

void
init_lucas(void)
{
  pari_sp ltop = avma;
  avma = ltop;
  return;
}

long
lucas(GEN p)
{
  pari_sp ltop = avma;
  long l1;	  /* bool */
  GEN u = gzero, q = gzero;
  u = stoi(4);
  q = subis(shifti(gun, itos(p)), 1);
  {
    pari_sp btop = avma, st_lim = stack_lim(btop, 1);
    GEN k = gzero;
    GEN *bptr[] = { &k, &u };
    for (k = stoi(3); gcmp(k, p) <= 0; k = gaddgs(k, 1))
    {
      u = gmod(gsubgs(gsqr(u), 2), q);
      if (low_stack(st_lim, stack_lim(btop, 1)))
        gerepilemany(btop, bptr, 2);
    }
  }
  l1 = gcmp0(u);
  avma = ltop;
  return l1;
}




9>[H[JPARI/GP Workshop      IHP    2004/09/14      Slide 10/12
lucas(p[0;31m:small[0m) =
{
  local(u,q); u=4; q=1<<p - 1;
  for(k=3,p, u = (sqr(u)-2)%q);
  u == 0
}
$ gp2c lucas.gp | less
/*-*- compile-command: "/usr/bin/gcc -c -o lucas.gp.o -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -I/home/bill/pari/i686/include lucas.gp.c && /usr/bin/gcc -o lucas.gp.so -shared -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -Wl,-shared lucas.gp.o -lc -ldl -lm -L/home/bill/pari/i686/lib -lpari"; -*-*/
#include <pari/pari.h>
/*
GP;install("lucas","lD0,L,","lucas","./lucas.gp.so");
GP;install("init_lucas","v","init_lucas","./lucas.gp.so");
*/
long lucas(long p);
void init_lucas(void);
/*End of prototype*/

void
init_lucas(void)
{
  return;
}

long
lucas(long p)
{
  long l1;	  /* bool */
  GEN u, q;
  u = stoi(4);
  q = subis(shifti(gun, p), 1);
  {
    long k;
    for (k = 3; k <= p; ++k)
      u = gmod(gsubgs(gsqr(u), 2), q);
  }
  l1 = gcmp0(u);
  return l1;
}




10>[H[JPARI/GP Workshop      IHP    2004/09/14      Slide 11/12
lucas(p) =
{
  local(u,q); u=4; q=1<<p - 1;
  for(k[0;31m:small[0m=3,p, u = (sqr(u)-2)%q);
  u == 0
}
$ gp2c lucas.gp | less
/*-*- compile-command: "/usr/bin/gcc -c -o lucas.gp.o -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -I/home/bill/pari/i686/include lucas.gp.c && /usr/bin/gcc -o lucas.gp.so -shared -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -Wl,-shared lucas.gp.o -lc -ldl -lm -L/home/bill/pari/i686/lib -lpari"; -*-*/
#include <pari/pari.h>
/*
GP;install("lucas","lD0,G,","lucas","./lucas.gp.so");
GP;install("init_lucas","v","init_lucas","./lucas.gp.so");
*/
long lucas(GEN p);
void init_lucas(void);
/*End of prototype*/

void
init_lucas(void)
{
  return;
}

long
lucas(GEN p)
{
  long l1;	  /* bool */
  GEN u, q;
  u = stoi(4);
  q = subis(shifti(gun, itos(p)), 1);
  {
    long k;
    for (k = 3; gcmpsg(k, p) <= 0; ++k)
      u = gmod(gsubgs(gsqr(u), 2), q);
  }
  l1 = gcmp0(u);
  return l1;
}




11>[H[JPARI/GP Workshop      IHP    2004/09/14      Slide 12/12
lucas(p[0;31m:small[0m) =
{
  local(u[0;31m:int[0m,q:int); u=4; q=1<<p - 1;
  for(k=3,p, u = (sqr(u)-2)%q);
  u == 0
}
$ gp2c lucas.gp | less
/*-*- compile-command: "/usr/bin/gcc -c -o lucas.gp.o -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -I/home/bill/pari/i686/include lucas.gp.c && /usr/bin/gcc -o lucas.gp.so -shared -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -Wl,-shared lucas.gp.o -lc -ldl -lm -L/home/bill/pari/i686/lib -lpari"; -*-*/
#include <pari/pari.h>
/*
GP;install("lucas","lD0,L,","lucas","./lucas.gp.so");
GP;install("init_lucas","v","init_lucas","./lucas.gp.so");
*/
long lucas(long p);
void init_lucas(void);
/*End of prototype*/

void
init_lucas(void)
{
  return;
}

long
lucas(long p)
{
  long l1;	  /* bool */
  GEN u, q;	  /* int */
  u = stoi(4);
  q = subis(shifti(gun, p), 1);
  {
    long k;
    for (k = 3; k <= p; ++k)
      u = modii(subis(sqri(u), 2), q);
  }
  l1 = cmpis(u, 0) == 0;
  return l1;
}




12>[H[JPARI/GP Workshop      IHP    2004/09/14      Slide 12/12
lucas(p[0;31m:small[0m) =
{
  local(u[0;31m:int[0m,q:int); u=4; q=1<<p - 1;
  for(k=3,p, u = (sqr(u)-2)%q);
  u == 0
}
$ gp2c lucas.gp | less
/*-*- compile-command: "/usr/bin/gcc -c -o lucas.gp.o -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -I/home/bill/pari/i686/include lucas.gp.c && /usr/bin/gcc -o lucas.gp.so -shared -O3 -DGCC_INLINE -Wall -fomit-frame-pointer -Wl,-shared lucas.gp.o -lc -ldl -lm -L/home/bill/pari/i686/lib -lpari"; -*-*/
#include <pari/pari.h>
/*
GP;install("lucas","lD0,L,","lucas","./lucas.gp.so");
GP;install("init_lucas","v","init_lucas","./lucas.gp.so");
*/
long lucas(long p);
void init_lucas(void);
/*End of prototype*/

void
init_lucas(void)
{
  return;
}

long
lucas(long p)
{
  long l1;	  /* bool */
  GEN u, q;	  /* int */
  u = stoi(4);
  q = subis(shifti(gun, p), 1);
  {
    long k;
    for (k = 3; k <= p; ++k)
      u = modii(subis(sqri(u), 2), q);
  }
  l1 = cmpis(u, 0) == 0;
  return l1;
}

