cp's OEIS Frontend

This is a front-end for the Online Encyclopedia of Integer Sequences, made by Christian Perfect. The idea is to provide OEIS entries in non-ancient HTML, and then to think about how they're presented visually. The source code is on GitHub.

A073103 Number of solutions to x^4 == 1 (mod n).

Original entry on oeis.org

1, 1, 2, 2, 4, 2, 2, 4, 2, 4, 2, 4, 4, 2, 8, 8, 4, 2, 2, 8, 4, 2, 2, 8, 4, 4, 2, 4, 4, 8, 2, 8, 4, 4, 8, 4, 4, 2, 8, 16, 4, 4, 2, 4, 8, 2, 2, 16, 2, 4, 8, 8, 4, 2, 8, 8, 4, 4, 2, 16, 4, 2, 4, 8, 16, 4, 2, 8, 4, 8, 2, 8, 4, 4, 8, 4, 4, 8, 2, 32, 2, 4, 2, 8, 16, 2, 8, 8, 4, 8, 8, 4, 4, 2, 8, 16, 4, 2, 4, 8
Offset: 1

Views

Author

Benoit Cloitre, Aug 19 2002

Keywords

Comments

a(n) = 2*A060594(n) for n = 5, 10, 13, 15, 16, 17, 20, 25, 26, 29, .... This subsequence, which contains all the primes of form 4k+1, seems to be asymptotic to 2n.
Shadow transform of A123865. - Michel Marcus, Jun 06 2013

Crossrefs

Programs

  • Maple
    a:= n-> add(`if`(irem(j^4-1, n)=0, 1, 0), j=0..n-1):
    seq(a(n), n=1..120);  # Alois P. Heinz, Jun 06 2013
    # alternative
    A073103 := proc(n)
        local a,pf,p,r;
        a := 1 ;
        for pf in ifactors(n)[2] do
            p := op(1,pf);
            r := op(2,pf);
            if p = 2 then
                a := a*p^min(r-1,3) ;
            else
                if modp(p,4) = 1 then
                    a := 4*a ;
                else
                    a := 2*a ;
                end if;
            end if;
        end do:
        a ;
    end proc: # R. J. Mathar, Mar 02 2015
  • Mathematica
    a[n_] := Sum[If[Mod[j^4-1, n] == 0, 1, 0], {j, 0, n-1}]; Table[a[n], {n, 1, 120}] (* Jean-François Alcover, Jun 12 2015, after Alois P. Heinz *)
    f[2, e_] := 2^Min[e-1, 3]; f[p_, e_] := If[Mod[p, 4] == 1, 4, 2]; a[1] = 1; a[n_] := Times @@ f @@@ FactorInteger[n]; Array[a, 100] (* Amiram Eldar, Sep 17 2020 *)
  • PARI
    a(n)=sum(i=1,n,if((i^4-1)%n,0,1))
    
  • PARI
    a(n)=my(f=factor(n)); prod(i=1,#f~, if(f[i,1]==2, 2^min(f[i,2]-1,3), if(f[i,1]%4==1, 4, 2))) \\ Charles R Greathouse IV, Mar 02 2015
    
  • Python
    from math import prod
    from sympy import primefactors
    def A073103(n): return (1<>1) for p in primefactors(n>>s)) # Chai Wah Wu, Oct 26 2022

Formula

Sum_{k=1..n} a(k) seems to be asymptotic to C*n*Log(n) with C>1.4...(when Sum_{k=1..n} A060594(k) is asymptotic to C/2*n*Log(n)).
Multiplicative with a(p^e) = p^min(e-1, 3) if p = 2, 4 if p == 1 (mod 4), 2 if p == 3 (mod 4). - David W. Wilson, Jun 09 2005
In fact, Sum_{k=1..n} a(k) is asymptotic to c*n*log(n)^2 where 2*c=0.190876.... - Steven Finch, Aug 12 2009 [c = 7/(2*Pi^3) * Product_{p prime == 1 (mod 4)} (1 - 4/(p+1)^2) = 0.0954383605... (Finch et al., 2010). - Amiram Eldar, Mar 26 2021]