Package: pari

Version: 2.6.2

Severity: wishlist

Tags: patch

the current implementation of ispowerful does not use the full strength of the trial division it does. In particular, if there are no factors below B in the unfactored part n < B^3, it can be powerful (if and) only if n = p^2 for some prime p.

This should not be uncommon, and it replaces a hard calculation (factoring a composite, up to a C18) with an easy one (check if a number is square).

I have included a simple fix below. Minor optimizations are available for checking whether the number is small enough (e.g., since tridiv_bound is always a power of 2, just shift). Perhaps better would be to check every time a factor is found, though this introduces some tradeoffs.

long

ispowerful(GEN n)

{

  pari_sp av = avma;

  GEN F;

  ulong p, bound;

  long i, l, v;

  forprime_t S;

  if ((F = check_arith_all(n, "ispowerful")))

  {

    GEN P = gel(F,1), E = gel(F,2);

    if (lg(P) == 1) return 1; /* 1 */

    if (!signe(gel(P,1))) return 1; /* 0 */

    l = lg(E);

    for (i = 1; i < l; i++)

      if (equali1(gel(E,i))) return 0;

    return 1;

  }

  if (!signe(n)) return 1;

  if (mod4(n) == 2) return 0;

  n = shifti(n, -vali(n));

  if (is_pm1(n)) return 1;

  setabssign(n);

  u_forprime_init(&S, 3, bound = tridiv_bound(n));

  while ((p = u_forprime_next_fast(&S)))

  {

    int stop;

    v = Z_lvalrem_stop(&n, p, &stop);

    if (v)

    {

      if (v == 1) { avma = av; return 0; }

      if (stop) { avma = av; return is_pm1(n); }

    }

  }

  

  l = lg(primetab);

  for (i = 1; i < l; i++)

  {

    v = Z_pvalrem(n, gel(primetab,i), &n);

    if (v)

    {

      if (v == 1) { avma = av; return 0; }

      if (is_pm1(n)) { avma = av; return 1; }

    }

  }

  

  /* no need to factor: must be p^2 or not powerful */

  if(cmpii(powuu(bound+1, 3), n) > 0) {

    long res = Z_issquare(n);

    avma = av;

    return res;

  }

  

  if (ifac_isprime(n)) { avma=av; return 0; }

  /* large composite without small factors */

  v = ifac_ispowerful(n);

  avma = av; return v;

}

Charles Greathouse
Analyst/Programmer
Case Western Reserve University