| Ilya Zakharevich on Tue, 2 Jul 2002 15:18:11 -0400 |
[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]
| Re: [PATCH] Hiccups in diffptr |
[quoting a very old message...]
On Fri, Feb 05, 1999 at 04:24:31AM -0500, I wrote:
> This is a first shot to enable yet bigger prime tables in PARI.
>
> We introduce hiccups into diffptr table in a backward-compatible way. The
> "real" end of the table is marked by "\0\0", the differences which are
> larger than 256 are introduced as
This is a much simpler implementation. First of all, it wraps all the
zillions of different ways to run through diffptr into 3 macros:
NEXT_PRIME_VIADIFF_PRECHECK, NEXT_PRIME_VIADIFF_POSTCHECK, NEXT_PRIME_VIADIFF
(I suspect all the usages of NEXT_PRIME_VIADIFF_POSTCHECK are bugs,
but I'm not sure). Second, we implement the large gaps in the
following ways: the gap 255 is fake, in the sense that it does not
correpond to a prime number. I.e, the sequence of gaps 3, 290, 5 is
coded into diffptr as
5, 255, 35, 5.
But this particular implementation is easily switchable to a different
one (or off) due to macroization.
Enjoy,
Ilya
--- ./src/modules/elliptic.c~ Sat Jun 8 16:02:31 2002
+++ ./src/modules/elliptic.c Tue Jul 2 19:50:10 2002
@@ -3055,7 +3055,7 @@ torsellnagelllutz(GEN e)
static long
torsbound(GEN e)
{
- long m, b, c, d, prime = 2;
+ long m, b, c, prime = 2;
gpmem_t av = avma;
byteptr p = diffptr;
GEN D = (GEN)e[12];
@@ -3064,8 +3064,7 @@ torsbound(GEN e)
b = c = m = 0; p++;
while (m<n)
{
- d = *p++; if (!d) err(primer1);
- prime += d;
+ NEXT_PRIME_VIADIFF_PRECHECK(prime,p);
if (smodis(D, prime))
{
b = cgcd(b, prime+1 - itos(apell0(e,prime)));
--- ./src/modules/aprcl.c~ Wed May 15 15:09:49 2002
+++ ./src/modules/aprcl.c Tue Jul 2 21:03:52 2002
@@ -914,9 +914,9 @@ step5(GEN N, int p, GEN et)
const ulong ltab = (NBITSN/kglob)+2;
av = avma;
- for (q = 3; *d; q += *d++)
+ for (q = 3; *d; )
{
- if (q%p != 1 || smodis(et,q) == 0) continue;
+ if (q%p != 1 || smodis(et,q) == 0) goto repeat;
if (smodis(N,q) == 0) return -1;
k = u_val(q-1, p);
@@ -931,6 +931,8 @@ step5(GEN N, int p, GEN et)
if (fl == 1) {ctglob = max(ctglob,ct); return 1;}
ct++;
avma = av;
+ repeat:
+ NEXT_PRIME_VIADIFF(q,d);
}
return 0;
}
--- ./src/modules/nffactor.c~ Mon Jun 10 21:43:22 2002
+++ ./src/modules/nffactor.c Tue Jul 2 19:58:56 2002
@@ -179,7 +179,7 @@ choose_prime(GEN nf, GEN bad, GEN *p, by
gpmem_t av = avma;
for (;;)
{
- pp += *pt++; if (!*pt) err(primer1);
+ NEXT_PRIME_VIADIFF_POSTCHECK(pp, pt);
if (! smodis(bad,pp)) continue;
affui(pp, q);
r = rootmod(x, q); if (lg(r) > 1) break;
--- ./src/modules/galois.c~ Fri May 10 17:10:05 2002
+++ ./src/modules/galois.c Tue Jul 2 19:51:46 2002
@@ -226,9 +226,9 @@ galmodp(GEN pol, GEN dpol, GEN TYP, long
dtyp = new_chunk(NMAX+1);
k = gr[0]; for (i=1; i<k; i++) gr[i]=1;
- for (k=1; k<15; k++, d++)
+ for (k=1; k<15; k++)
{
- p += *d; if (!*d) err(primer1);
+ NEXT_PRIME_VIADIFF_PRECHECK(p,d);
if (!smodis(dpol,p)) continue; /* p divides dpol */
p1 = simplefactmod(pol,stoi(p));
--- ./src/modules/mpqs.c~ Thu Apr 18 00:10:58 2002
+++ ./src/modules/mpqs.c Tue Jul 2 20:12:27 2002
@@ -653,8 +653,9 @@ static byteptr
mpqs_iterate_primes(long *p, byteptr primes_ptr)
{
long prime = *p;
- if (*primes_ptr)
- prime += *primes_ptr++;
+ if (*primes_ptr) {
+ NEXT_PRIME_VIADIFF(prime,primes_ptr)
+ }
else
{
gpmem_t av = avma;
@@ -749,9 +750,12 @@ static long
mpqs_count_primes(void)
{
byteptr p = mpqs_diffptr;
+ long gaps = 0;
- for ( ; *p; p++) /* empty */;
- return (p - mpqs_diffptr);
+ for ( ; *p; p++)
+ if (*p == DIFFPTR_SKIP)
+ gaps++;
+ return (p - mpqs_diffptr - gaps);
}
/**
--- ./src/modules/stark.c~ Sat Jun 8 16:37:45 2002
+++ ./src/modules/stark.c Tue Jul 2 20:01:35 2002
@@ -1541,7 +1541,7 @@ InitPrimesQuad(GEN bnr, long nmax, LISTr
R->L1 = _alloc(l, t_VECSMALL); R->L1ray = (GEN*)_alloc(l, t_VEC);
R->L11 = _alloc(l, t_VECSMALL); R->L11ray = (GEN*)_alloc(l, t_VEC);
prime = stoi(2);
- for (p = 2; p <= nmax; p += *d++, prime[2] = p)
+ for (p = 2; p <= nmax; prime[2] = p) {
switch (krogs(dk, p))
{
case -1: /* inert */
@@ -1566,6 +1566,8 @@ InitPrimesQuad(GEN bnr, long nmax, LISTr
}
break;
}
+ NEXT_PRIME_VIADIFF(p,d);
+ }
/* precompute isprincipalray(x), x in Z */
R->rayZ = (GEN*)cgetg(condZ, t_VEC);
for (i=1; i<condZ; i++)
@@ -1592,7 +1594,7 @@ InitPrimes(GEN bnr, long nmax, LISTray *
R->L1 = _alloc(nmax, t_VECSMALL);
R->L1ray = (GEN*)_alloc(nmax, t_VEC);
prime = stoi(2);
- for (p = 2; p <= nmax; p += *d++, prime[2] = p)
+ for (p = 2; p <= nmax; prime[2] = p)
{
tabpr = primedec(nf, prime);
for (j = 1; j < lg(tabpr); j++)
@@ -1605,7 +1607,7 @@ InitPrimes(GEN bnr, long nmax, LISTray *
_append(R->L1, (GEN)np);
_append((GEN)R->L1ray, isprincipalray(bnr, pr));
}
-
+ NEXT_PRIME_VIADIFF(p,d);
}
gptr[0] = &(R->L1); gptr[1] = (GEN*)&(R->L1ray);
gerepilemany(av,gptr,2);
--- ./src/modules/subfield.c~ Sun Jun 9 22:18:18 2002
+++ ./src/modules/subfield.c Tue Jul 2 21:06:33 2002
@@ -610,25 +610,27 @@ choose_prime(GEN pol,GEN dpol,long d,GEN
minp = N*(m-1)/2;
if (DEBUGLEVEL) (void)timer2();
di++; p = stoi(2);
- while (p[2]<=minp) p[2] += *di++;
+ while (p[2]<=minp)
+ NEXT_PRIME_VIADIFF(p[2], di);
oldllist = 100000;
oldlcm = 0;
oldlistpotbl = oldff = oldn = NULL; pp = 0; /* gcc -Wall */
av = avma;
- for(k=1; k<11 || !oldn; k++,p[2]+= *di++,avma = av)
+ for(k=1; k<11 || !oldn; k++,avma = av)
{
- while (!smodis(dpol,p[2])) p[2] += *di++;
+ while (!smodis(dpol,p[2]))
+ NEXT_PRIME_VIADIFF(p[2], di);
if (k > 50) err(talker,"sorry, too many block systems in nfsubfields");
ff=(GEN)factmod0(pol,p)[1]; r=lg(ff)-1;
- if (r == 1 || r == N) continue;
+ if (r == 1 || r == N) goto repeat;
n = cgetg(r+1, t_VECSMALL);
lcm = n[1] = degpol(ff[1]);
for (j=2; j<=r; j++) { n[j] = degpol(ff[j]); lcm = clcm(lcm,n[j]); }
- if (lcm < oldlcm) continue; /* false when oldlcm = 0 */
- if (r >= BIL) { err(warner,"subfields: overflow in calc_block"); continue; }
+ if (lcm < oldlcm) goto repeat; /* false when oldlcm = 0 */
+ if (r >= BIL) { err(warner,"subfields: overflow in calc_block"); goto repeat; }
if (DEBUGLEVEL) fprintferr("p = %ld,\tlcm = %ld,\torbits: %Z\n",p[2],lcm,n);
- if (oldn && gegal(n,oldn)) continue;
+ if (oldn && gegal(n,oldn)) goto repeat;
listpotbl = potential_block_systems(N,d,n, oldllist);
if (!listpotbl) { oldlistpotbl = NULL; pp = p[2]; break; }
@@ -637,12 +639,14 @@ choose_prime(GEN pol,GEN dpol,long d,GEN
{
if (DEBUGLEVEL) msgtimer("#pbs >= %ld [aborted]",oldllist);
for (j=1; j<llist; j++) free((void*)listpotbl[j]);
- free((void*)(listpotbl-1)); continue;
+ free((void*)(listpotbl-1)); goto repeat;
}
oldllist = llist; oldlistpotbl = listpotbl;
pp = p[2]; oldff = ff; oldlcm = lcm; oldn = n;
if (DEBUGLEVEL) msgtimer("#pbs = %ld",llist);
av = avma;
+ repeat:
+ NEXT_PRIME_VIADIFF(p[2], di);
}
if (DEBUGLEVEL) fprintferr("Chosen prime: p = %ld\n",pp);
*ptlistpotbl=oldlistpotbl; *ptff=oldff; *ptlcm=oldlcm; return stoi(pp);
--- ./src/language/sumiter.c~ Tue Jul 2 18:06:40 2002
+++ ./src/language/sumiter.c Tue Jul 2 19:46:11 2002
@@ -104,7 +104,8 @@ sinitp(long a, long c, byteptr *ptr)
byteptr p = *ptr;
if (a <= 0) a = 2;
if (maxprime() < (ulong)a) err(primer1);
- while (a > c) c += *p++;
+ while (a > c)
+ NEXT_PRIME_VIADIFF(c,p);
*ptr = p; return c;
}
--- ./src/headers/paricom.h~ Thu Apr 11 01:29:53 2002
+++ ./src/headers/paricom.h Tue Jul 2 19:15:39 2002
@@ -315,3 +315,11 @@ extern void* global_err_data;
#define cmp_LEX 2
#define cmp_REV 4
#define cmp_C 8
+
+#define DIFFPTR_SKIP 255 /* Skip these entries */
+#define NEXT_PRIME_VIADIFF(p,d) \
+ { while (*(d) == DIFFPTR_SKIP) (p) += *(d)++; (p) += *(d)++; }
+#define NEXT_PRIME_VIADIFF_PRECHECK(p,d) \
+ { if (!*(d)) err(primer1); NEXT_PRIME_VIADIFF(p,d) }
+#define NEXT_PRIME_VIADIFF_POSTCHECK(p,d) \
+ { NEXT_PRIME_VIADIFF(p,d); if (!*(d)) err(primer1); }
--- ./src/basemath/arith2.c~ Tue Jul 2 18:06:39 2002
+++ ./src/basemath/arith2.c Tue Jul 2 20:05:57 2002
@@ -39,13 +39,12 @@ GEN
prime(long n)
{
byteptr p = diffptr;
- long c, prime = 0;
+ long prime = 0;
if (n <= 0) err(talker, "n-th prime meaningless if n = %ld",n);
while (n--)
{
- c = *p++; if (!c) err(primer1);
- prime += c;
+ NEXT_PRIME_VIADIFF_PRECHECK(prime,p);
}
return stoi(prime);
}
@@ -58,7 +57,10 @@ pith(long n)
if (n <= 0) err(talker, "pith meaningless if n = %ld",n);
if (maxprime() <= (ulong)n) err(primer1);
- while (prime<=(ulong)n) { prime += *p++; res++; }
+ while (prime<=(ulong)n) {
+ NEXT_PRIME_VIADIFF(prime,p);
+ res++;
+ }
return stoi(res-1);
}
@@ -66,15 +68,15 @@ GEN
primes(long n)
{
byteptr p = diffptr;
- long c, prime = 0;
+ long prime = 0;
GEN y,z;
if (n < 0) n = 0;
z = y = cgetg(n+1,t_VEC);
while (n--)
{
- c = *p++; if (!c) err(primer1);
- prime += c; *++z = lstoi(prime);
+ NEXT_PRIME_VIADIFF_PRECHECK(prime,p);
+ *++z = lstoi(prime);
}
return y;
}
@@ -208,9 +210,14 @@ initprimes0(ulong maxnum, long *lenp, lo
/* now q runs over addresses corresponding to primes */
for (q = p; ; plast = q++)
{
+ long d;
+
while (*q) q++; /* use 0-sentinel at end */
if (q >= fin) break;
- *curdiff++ = (q - plast) << 1;
+ d = (q - plast) << 1;
+ while (d >= DIFFPTR_SKIP)
+ *curdiff++ = DIFFPTR_SKIP, d -= DIFFPTR_SKIP;
+ *curdiff++ = d;
}
plast -= asize - 1;
remains -= asize - 1;
@@ -247,9 +254,6 @@ initprimes(long maxnum)
/* switch to unsigned: adding the 257 _could_ wrap to a negative number. */
ulong maxnum1 = ((maxnum<65302?65302:maxnum)+257ul);
- if (maxnum1 > 436273000)
- err(talker,"impossible to have prestored primes > 436273009");
-
p = initprimes0(maxnum1, &len, &last);
#if 0 /* not yet... GN */
static int build_the_tables = 1;
@@ -438,7 +442,7 @@ auxdecomp1(GEN n, long (*ifac_break)(GEN
p = 2;
while (*d && p <= lim1)
{
- p += *d++;
+ NEXT_PRIME_VIADIFF(p,d);
if (mpdivisis(n,p,n))
{
nb++; k=1; while (mpdivisis(n,p,n)) k++;
@@ -693,7 +697,7 @@ mu(GEN n)
while (*d && p < lim1)
{
- p += *d++;
+ NEXT_PRIME_VIADIFF(p,d);
if (mpdivisis(n,p,n))
{
if (smodis(n,p) == 0) { avma=av; return 0; }
@@ -739,7 +743,7 @@ issquarefree(GEN x)
p = 2;
while (*d && p < lim1)
{
- p += *d++;
+ NEXT_PRIME_VIADIFF(p,d);
if (mpdivisis(x,p,x))
{
if (smodis(x,p) == 0) { avma = av; return 0; }
@@ -779,7 +783,7 @@ omega(GEN n)
while (*d && p < lim1)
{
- p += *d++;
+ NEXT_PRIME_VIADIFF(p,d);
if (mpdivisis(n,p,n))
{
nb++; while (mpdivisis(n,p,n)); /* empty */
@@ -816,7 +820,7 @@ bigomega(GEN n)
while (*d && p < lim1)
{
- p += *d++;
+ NEXT_PRIME_VIADIFF(p,d);
if (mpdivisis(n,p,n))
{
do nb++; while (mpdivisis(n,p,n));
@@ -854,7 +858,7 @@ phi(GEN n)
while (*d && p < lim1)
{
- p += *d++;
+ NEXT_PRIME_VIADIFF(p,d);
if (mpdivisis(n,p,n))
{
m = mulis(m, p-1);
@@ -897,7 +901,7 @@ numbdiv(GEN n)
while (*d && p < lim1)
{
- p += *d++;
+ NEXT_PRIME_VIADIFF(p,d);
l=1; while (mpdivisis(n,p,n)) l++;
m = mulsi(l,m); if (is_pm1(n)) { return gerepileupto(av,m); }
}
@@ -936,7 +940,7 @@ sumdiv(GEN n)
while (*d && p < lim1)
{
- p += *d++;
+ NEXT_PRIME_VIADIFF(p,d);
if (mpdivisis(n,p,n))
{
m1 = utoi(p+1);
@@ -985,7 +989,7 @@ sumdivk(GEN n, long k)
while (*d && p < lim1)
{
- p += *d++;
+ NEXT_PRIME_VIADIFF(p,d);
if (mpdivisis(n,p,n))
{
pk = gpowgs(utoi(p),k); m1 = addsi(1,pk);
--- ./src/basemath/ifactor1.c~ Tue Jul 2 18:06:40 2002
+++ ./src/basemath/ifactor1.c Tue Jul 2 19:36:10 2002
@@ -508,10 +508,14 @@ snextpr(ulong p, byteptr *d, long *rcn,
{ evaltyp(t_INT)|_evallg(3), evalsigne(1)|evallgefint(3), 0 };
static ulong *pp2 = pp + 2;
static GEN gp = (GEN)pp;
- long d1 = **d, rcn0;
+ long rcn0;
- if (d1)
+ if (**d)
{
+ byteptr dd = *d;
+ long d1 = 0;
+
+ NEXT_PRIME_VIADIFF(d1,dd);
if (*rcn != NPRC)
{
rcn0 = *rcn;
@@ -527,7 +531,8 @@ snextpr(ulong p, byteptr *d, long *rcn,
err(bugparier, "[caller of] snextpr");
}
}
- return p + *(*d)++;
+ NEXT_PRIME_VIADIFF(p,*d);
+ return p;
}
/* we are beyond the diffptr table */
if (*rcn == NPRC) /* we need to initialize this now */
@@ -1313,7 +1318,7 @@ ellfacteur(GEN n, int insist)
fprintferr("\tB2 = %6lu,\tgss = %4ld*420\n", B2, gss);
flusherr();
}
- p = *d++;
+ NEXT_PRIME_VIADIFF(p,d);
/* ---B1 PHASE--- */
/* treat p=2 separately */
--- ./src/basemath/arith1.c~ Sat Jun 8 16:02:28 2002
+++ ./src/basemath/arith1.c Tue Jul 2 19:14:31 2002
@@ -2618,7 +2618,8 @@ classno(GEN x)
c=lforms=0;
while (c <= s && *p)
{
- c += *p++; k = krogs(D,c);
+ NEXT_PRIME_VIADIFF(c,p);
+ k = krogs(D,c);
if (!k) continue;
av2 = avma;
--- ./src/basemath/alglin1.c~ Mon May 13 20:25:46 2002
+++ ./src/basemath/alglin1.c Tue Jul 2 21:08:11 2002
@@ -1529,11 +1529,10 @@ ZM_inv(GEN M, GEN dM)
av2 = avma;
H = NULL;
d += 3000; /* 27449 = prime(3000) */
- for(p = 27449; ; p+= *d++)
+ for(p = 27449; ; )
{
- if (!*d) err(primer1);
dMp = umodiu(dM,p);
- if (!dMp) continue;
+ if (!dMp) goto repeat;
Hp = u_FpM_inv_sp(u_Fp_FpM(M,p), p);
if (dMp != 1) Hp = u_FpM_Fp_mul_ip(Hp, dMp, p);
@@ -1557,6 +1556,8 @@ ZM_inv(GEN M, GEN dM)
if (DEBUGMEM>1) err(warnmem,"ZM_inv");
gerepilemany(av2,gptr, 2);
}
+ repeat:
+ NEXT_PRIME_VIADIFF_PRECHECK(p,d);
}
if (DEBUGLEVEL>5) msgtimer("ZM_inv done");
return gerepilecopy(av, H);
--- ./src/basemath/base3.c~ Mon May 20 14:19:32 2002
+++ ./src/basemath/base3.c Tue Jul 2 20:06:59 2002
@@ -2003,7 +2003,7 @@ ideallistzstarall(GEN bnf,long bound,lon
if (bound > (long)maxprime()) err(primer1);
for (p[2]=0; p[2]<=bound; )
{
- p[2] += *ptdif++;
+ NEXT_PRIME_VIADIFF(p[2], ptdif);
if (DEBUGLEVEL>1) { fprintferr("%ld ",p[2]); flusherr(); }
fa = primedec(nf,p);
for (j=1; j<lg(fa); j++)
--- ./src/basemath/base1.c~ Fri Jun 7 19:41:17 2002
+++ ./src/basemath/base1.c Tue Jul 2 19:16:11 2002
@@ -2078,7 +2078,7 @@ dirzetak0(GEN nf, long N0)
while (court[2]<=N0)
{
- court[2] += *d++; if (! *d) err(primer1);
+ NEXT_PRIME_VIADIFF_POSTCHECK(court[2], d);
if (smodis(disc,court[2])) /* court does not divide index */
{ vect = (GEN) simplefactmod(pol,court)[1]; lx=lg(vect); }
else
--- ./src/basemath/buch1.c~ Sat Jun 8 16:02:29 2002
+++ ./src/basemath/buch1.c Tue Jul 2 19:19:49 2002
@@ -148,7 +148,7 @@ get_pq(GEN D, GEN z, GEN flag, GEN *ptp,
ell = 3;
while (l < 3 || ell<=MAXL)
{
- ell += *diffell++; if (!*diffell) err(primer1);
+ NEXT_PRIME_VIADIFF_POSTCHECK(ell, diffell);
if (smodis(z,ell) && kross(d,ell) > 0)
{
court[2]=ell; form = redimag(primeform(D,court,0));
@@ -979,7 +979,7 @@ factorbasequad(GEN Disc, long n2, long n
default: /* split */
i++; numfactorbase[p]=i; factorbase[i] = p;
}
- p += *d++; if (!*d) err(primer1);
+ NEXT_PRIME_VIADIFF_POSTCHECK(p, d);
if (KC == 0 && p>n) KC = i;
}
if (!KC) { free(factorbase); free(numfactorbase); return; }
--- ./src/basemath/buch2.c~ Fri Jun 7 23:13:26 2002
+++ ./src/basemath/buch2.c Tue Jul 2 19:20:29 2002
@@ -272,8 +272,7 @@ FBgen(GEN nf,long n2,long n)
{ setlg(p1,k); p1 = gerepilecopy(av,p1); }
idealbase[i] = p1;
}
- if (!*delta) err(primer1);
- p += *delta++;
+ NEXT_PRIME_VIADIFF_PRECHECK(p, delta);
if (KC == 0 && p>n) { KCZ=i; KC=ip; }
}
if (!KC) return NULL;
--- ./src/basemath/buch3.c~ Fri Jun 7 23:13:27 2002
+++ ./src/basemath/buch3.c Tue Jul 2 19:24:23 2002
@@ -763,7 +763,7 @@ testprime(GEN bnf, long minkowski)
if ((ulong)minkowski > maxprime()) err(primer1);
while (pp < minkowski)
{
- pp += *delta++;
+ NEXT_PRIME_VIADIFF(pp, delta);
if (DEBUGLEVEL>1) fprintferr("*** p = %ld\n",pp);
vectpp=primedec(bnf,stoi(pp)); nbideal=lg(vectpp)-1;
/* loop through all P | p if ramified, all but one otherwise */
@@ -1318,8 +1318,10 @@ certifybuchall(GEN bnf)
rootsofone = dummycopy(rootsofone);
rootsofone[2] = (long)algtobasis(nf, (GEN)rootsofone[2]);
- for (p = *delta++; p <= bound; p += *delta++)
+ for (p = *delta++; p <= bound; ) {
check_prime(p,bnf,cyc,cycgen,funits,rootsofone,big);
+ NEXT_PRIME_VIADIFF(p, delta);
+ }
if (nbgen)
{
@@ -1548,7 +1550,7 @@ rnfnormgroup(GEN bnr, GEN polrel)
* and including last pr^f or p^f is the same, but the last isprincipal
* is much easier! oldf is used to track this */
- p += *d++; if (!*d) err(primer1);
+ NEXT_PRIME_VIADIFF_POSTCHECK(p,d);
if (!smodis(index, p)) continue; /* can't be treated efficiently */
fa = primedec(nf, stoi(p)); lfa = lg(fa)-1;
@@ -2280,7 +2282,7 @@ discrayabslistarchintern(GEN bnf, GEN ar
if (bound > (long)maxprime()) err(primer1);
for (p[2]=0; p[2]<=bound; )
{
- p[2] += *ptdif++;
+ NEXT_PRIME_VIADIFF(p[2], ptdif);
if (!flbou && p[2]>sqbou)
{
if (DEBUGLEVEL>1) fprintferr("\nStarting rayclassno computations\n");
--- ./src/basemath/polarit2.c~ Sun Jun 9 17:29:09 2002
+++ ./src/basemath/polarit2.c Tue Jul 2 20:50:38 2002
@@ -1344,7 +1344,7 @@ DDF(GEN a, long hint)
{
gpmem_t av0 = avma;
- p += *pt++; if (!*pt) err(primer1);
+ NEXT_PRIME_VIADIFF_POSTCHECK(p,pt);
if (lead && !smodis(lead,p)) continue;
z = u_Fp_FpX(a,0, p);
if (!u_FpX_is_squarefree(z, p)) { avma = av0; continue ; }
@@ -4165,21 +4165,21 @@ nfgcd(GEN P, GEN Q, GEN nf, GEN den)
GEN M, dsol;
GEN R, ax, bo;
byteptr primepointer;
- for (p = 27449, primepointer = diffptr + 3000; ; p += *(primepointer++))
+ for (p = 27449, primepointer = diffptr + 3000; ; )
{
if (!*primepointer) err(primer1);
if (!smodis(den, p))
- continue;/*Discard primes dividing disc(T) or leadingcoeff(PQ) */
+ goto repeat;/*Discard primes dividing disc(T) or leadingcoeff(PQ) */
if (DEBUGLEVEL>5) fprintferr("nfgcd: p=%d\n",p);
if ((R = FpXQX_safegcd(P, Q, nf, stoi(p))) == NULL)
- continue;/*Discard primes when modular gcd does not exist*/
+ goto repeat;/*Discard primes when modular gcd does not exist*/
dR = degpol(R);
if (dR == 0) return scalarpol(gun, x);
- if (mod && dR > dM) continue; /* p divides Res(P/gcd, Q/gcd). Discard. */
+ if (mod && dR > dM) goto repeat; /* p divides Res(P/gcd, Q/gcd). Discard. */
R = polpol_to_mat(R, d);
/* previous primes divided Res(P/gcd, Q/gcd)? Discard them. */
- if (!mod || dR < dM) { M = R; mod = stoi(p); dM = dR; continue; }
+ if (!mod || dR < dM) { M = R; mod = stoi(p); dM = dR; goto repeat; }
if (low_stack(st_lim, stack_lim(btop, 1)))
{
if (DEBUGMEM>1) err(warnmem,"nfgcd");
@@ -4193,11 +4193,13 @@ nfgcd(GEN P, GEN Q, GEN nf, GEN den)
/* I suspect it must be better to take amax > bmax*/
bo = racine(shifti(mod, -1));
if ((sol = matratlift(M, mod, bo, bo, den)) == NULL)
- continue;
+ goto repeat;
sol = mat_to_polpol(sol,x,y);
dsol = primpart(sol);
if (gcmp0(pseudorem_i(P, dsol, nf))
&& gcmp0(pseudorem_i(Q, dsol, nf))) break;
+ repeat:
+ NEXT_PRIME_VIADIFF_POSTCHECK(p, primepointer);
}
}
return gerepilecopy(ltop, sol);
--- ./src/basemath/galconj.c~ Mon Jun 10 16:20:58 2002
+++ ./src/basemath/galconj.c Tue Jul 2 20:10:05 2002
@@ -1653,13 +1653,10 @@ galoisanalysis(GEN T, struct galois_anal
&& (nbtest < 3 * nbmax || !(group&ga_non_wss)) ;)
{
gpmem_t av;
- long prime_incr;
GEN ip,FS,p1;
long o,norm_o=1;
- prime_incr = *primepointer++;
- if (!prime_incr)
- err(primer1);
- p += prime_incr;
+
+ NEXT_PRIME_VIADIFF_PRECHECK(p,primepointer);
/*discard small primes*/
if (p <= min_prime)
continue;
@@ -1745,17 +1742,14 @@ galoisanalysis(GEN T, struct galois_anal
if (calcul_l && O[1]<=linf)
{
gpmem_t av;
- long prime_incr;
long l=0;
/*we need a totally splited prime l*/
av = avma;
while (l == 0)
{
long nb;
- prime_incr = *primepointer++;
- if (!prime_incr)
- err(primer1);
- p += prime_incr;
+
+ NEXT_PRIME_VIADIFF_PRECHECK(p,primepointer);
if (p<=linf) continue;
nb=FpX_nbroots(T,stoi(p));
if (nb == n)
@@ -2628,7 +2622,7 @@ galoisfindfrobenius(GEN T, GEN L, GEN de
{
gpmem_t av = avma;
long isram;
- long i,c;
+ long i;
GEN ip,Tmod;
ip = stoi(gf->p);
Tmod = lift((GEN) factmod(T, ip));
@@ -2673,10 +2667,7 @@ galoisfindfrobenius(GEN T, GEN L, GEN de
err(warner, "galoisconj _may_ hang up for this polynomial");
}
}
- c = *primepointer++;
- if (!c)
- err(primer1);
- gf->p += c;
+ NEXT_PRIME_VIADIFF_PRECHECK(gf->p, primepointer);
if (DEBUGLEVEL >= 4)
fprintferr("GaloisConj:next p=%ld\n", gf->p);
avma = av;
@@ -3057,13 +3048,11 @@ numberofconjugates(GEN T, long pdepart)
ltop2 = avma;
for (p = 0, primepointer = diffptr; nbtest < nbmax && card > 1;)
{
- long s, c;
+ long s;
long isram;
GEN S;
- c = *primepointer++;
- if (!c)
- err(primer1);
- p += c;
+
+ NEXT_PRIME_VIADIFF_PRECHECK(p,primepointer);
if (p < pdepart)
continue;
S = (GEN) simplefactmod(T, stoi(p));
--- ./src/basemath/polarit3.c~ Sun Jun 9 22:10:11 2002
+++ ./src/basemath/polarit3.c Tue Jul 2 21:11:00 2002
@@ -4224,7 +4224,7 @@ INIT:
}
/* make sure p large enough */
- while (p < (dres<<1)) p += *d++;
+ while (p < (dres<<1)) NEXT_PRIME_VIADIFF(p,d);
H = H0 = H1 = NULL;
lb = lgef(B); b = u_allocpol(degpol(B), 0);
@@ -4234,7 +4234,7 @@ INIT:
for(;;)
{
- p += *d++; if (!*d) err(primer1);
+ NEXT_PRIME_VIADIFF_POSTCHECK(p,d);
a = u_Fp_FpX(A, 0, p);
for (i=2; i<lb; i++)
@@ -4455,7 +4455,7 @@ ZX_resultant_all(GEN A, GEN B, GEN dB, u
dp = 0; /* gcc -Wall */
for(;;)
{
- p += *d++; if (!*d) err(primer1);
+ NEXT_PRIME_VIADIFF_POSTCHECK(p,d);
if (dB) { dp = smodis(dB, p); if (!dp) continue; }
a = u_Fp_FpX(A, 0, p);
@@ -4556,16 +4556,16 @@ modulargcd(GEN A0, GEN B0)
av2 = avma; H = NULL;
d += 3000; /* 27449 = prime(3000) */
- for(p = 27449; ; p+= *d++)
+ for(p = 27449; ;)
{
if (!*d) err(primer1);
- if (!umodiu(g,p)) continue;
+ if (!umodiu(g,p)) goto repeat;
a = u_Fp_FpX(A, 0, p);
b = u_Fp_FpX(B, 0, p); Hp = u_FpX_gcd_i(a,b, p);
m = degpol(Hp);
if (m == 0) { H = polun[varn(A0)]; break; } /* coprime. DONE */
- if (m > n) continue; /* p | Res(A/G, B/G). Discard */
+ if (m > n) goto repeat; /* p | Res(A/G, B/G). Discard */
if (is_pm1(g)) /* make sure lead(H) = g mod p */
Hp = u_FpX_normalize(Hp, p);
@@ -4577,7 +4577,7 @@ modulargcd(GEN A0, GEN B0)
if (m < n)
{ /* First time or degree drop [all previous p were as above; restart]. */
H = ZX_init_CRT(Hp,p,varn(A0));
- q = utoi(p); n = m; continue;
+ q = utoi(p); n = m; goto repeat;
}
qp = muliu(q,p);
@@ -4595,6 +4595,8 @@ modulargcd(GEN A0, GEN B0)
if (DEBUGMEM>1) err(warnmem,"modulargcd");
gerepilemany(av2,gptr,2);
}
+ repeat:
+ NEXT_PRIME_VIADIFF_POSTCHECK(p,d);
}
return gerepileupto(av, gmul(D,H));
}
@@ -4621,19 +4623,19 @@ QX_invmod(GEN A0, GEN B0)
/* A, B in Z[X] */
av2 = avma; U = NULL;
d += 3000; /* 27449 = prime(3000) */
- for(p = 27449; ; p+= *d++)
+ for(p = 27449; ; )
{
if (!*d) err(primer1);
a = u_Fp_FpX(A, 0, p);
b = u_Fp_FpX(B, 0, p);
/* if p | Res(A/G, B/G), discard */
- if (!u_FpX_extresultant(b,a,p, &Vp,&Up)) continue;
+ if (!u_FpX_extresultant(b,a,p, &Vp,&Up)) goto repeat;
if (!U)
{ /* First time */
U = ZX_init_CRT(Up,p,varn(A0));
V = ZX_init_CRT(Vp,p,varn(A0));
- q = utoi(p); continue;
+ q = utoi(p); goto repeat;
}
if (DEBUGLEVEL>5) msgtimer("QX_invmod: mod %ld (bound 2^%ld)", p,expi(q));
qp = muliu(q,p);
@@ -4652,6 +4654,8 @@ QX_invmod(GEN A0, GEN B0)
if (DEBUGMEM>1) err(warnmem,"QX_invmod");
gerepilemany(av2,gptr,3);
}
+ repeat:
+ NEXT_PRIME_VIADIFF_POSTCHECK(p,d);
}
D = D? gmul(D,res): res;
return gerepileupto(av, gdiv(U,D));
--- ./src/basemath/trans3.c~ Sat Jun 8 16:48:36 2002
+++ ./src/basemath/trans3.c Tue Jul 2 19:44:27 2002
@@ -1362,7 +1362,7 @@ n_s(ulong n, GEN *tab)
GEN x = NULL;
long p, e;
- for (p = 3; n > 1; p += *d++)
+ for (p = 3; n > 1; )
{
e = svaluation(n, p, &n);
if (e)
@@ -1370,6 +1370,7 @@ n_s(ulong n, GEN *tab)
GEN y = tab[pows(p,e)];
if (!x) x = y; else x = gmul(x,y);
}
+ NEXT_PRIME_VIADIFF_POSTCHECK(p,d);
}
return x;
}
@@ -1468,27 +1469,31 @@ czeta(GEN s0, long prec)
d = diffptr + 1;
if (typ(s0) == t_INT)
{ /* no explog for 1/p^s */
- for (p=2; p < nn; p += *d++)
+ for (p=2; p < nn;) {
tab[p] = divrr(unr, rpowsi(p, s0, prec));
+ NEXT_PRIME_VIADIFF_POSTCHECK(p,d);
+ }
a = divrr(unr, rpowsi(nn, s0, prec));
}
else
{ /* general case */
ms = gneg(s); p1 = cgetr(prec);
- for (p=2; p < nn; p += *d++)
+ for (p=2; p < nn;)
{
affsr(p, p1);
tab[p] = gexp(gmul(ms, mplog(p1)), prec);
+ NEXT_PRIME_VIADIFF_POSTCHECK(p,d);
}
affsr(nn,p1);
a = gexp(gmul(ms, mplog(p1)), prec);
}
sqn = (long)sqrt(nn-1.);
d = diffptr + 2; /* fill in odd prime powers */
- for (p=3; p <= sqn; p += *d++)
+ for (p=3; p <= sqn; )
{
ulong oldq = p, q = p*p;
while (q<(ulong)nn) { tab[q] = gmul(tab[p], tab[oldq]); oldq = q; q *= p; }
+ NEXT_PRIME_VIADIFF_POSTCHECK(p,d);
}
if (DEBUGLEVEL>2) msgtimer("tab[q^-s] from 1 to N-1");