login

Reminder: The OEIS is hiring a new managing editor, and the application deadline is January 26.

Let p = prime(n); a(n) = number of pairs of primes (i,j) with i < j < p such that i*j == 1 (mod p).
6

%I #57 Feb 28 2020 01:43:47

%S 0,0,1,1,0,1,1,2,1,0,1,2,2,2,2,2,2,3,4,3,2,3,1,2,2,4,2,3,3,4,5,1,3,4,

%T 2,3,4,4,5,1,4,3,6,4,6,4,5,6,4,4,5,4,5,4,8,7,3,5,5,7,7,6,4,10,9,9,8,6,

%U 4,9,9,7,8,6,8,9,9,5,9,6,9,9,6,4,10

%N Let p = prime(n); a(n) = number of pairs of primes (i,j) with i < j < p such that i*j == 1 (mod p).

%C _David James Sycamore_ and Karl Love computed this sequence in early October 2018 and studied its asymptotic growth.

%C See link for conjectured asymptotic behavior of this sequence, showing (red line) first asymptotic term and residual bounds (blue lines). It appears that a(n) ~ Li(p)^2/(2*p), based upon the assumption that the probabilities that the first and second numbers in a pair are primes are independent. The residuals (with minor exceptions for some small p) seem to be bounded by (2/Pi^2)*sqrt(p). These observations were contributed by Karl Love (via Mapleprimes). - _David James Sycamore_, Nov 05 2018

%C If i < j < p in the definition is replaced by i <= j < p, then the only change is that a(2) becomes 1. For if i^2 == 1 (mod p) then p divides i^2-1 = (i-1)*(i+1), which implies p=i+1, so p must be 3 and i must be 2. - _N. J. A. Sloane_, Nov 06 2018

%H Chai Wah Wu, <a href="/A321005/b321005.txt">Table of n, a(n) for n = 1..10000</a>

%H Karl Love, <a href="/A321005/a321005.txt">Maple programs for A321005</a>

%H David James Sycamore, <a href="/A321005/a321005_1.jpg">Asymptotic Behaviour of A321005</a>

%F An alternative heuristic argument leads the conjecture that a(n) ~ (n-4)/(2*log(n)) ~ Pi(n)/2. (Essentially the same as Karl Love's conjecture in the Comments section, but simpler to calculate.) - _N. J. A. Sloane_, Nov 06 2018

%p with(numtheory);

%p f:=proc(n) local m,i,j,c; global pi;

%p if n < 4 then return(0); fi;

%p c:=0; m:=NumberTheory[pi](n); if isprime(n) then m:=m-1; fi;

%p for i from 1 to m-1 do p:=ithprime(i);

%p for j from i+1 to m do

%p if ithprime(j)*p mod n = 1 then c:=c+1; fi; od: od: c; end;

%p [seq(f(ithprime(n)),n=1..120)];

%t a[n_] := Module[{p=Prime[n]}, c=0; Do[Do[If[Mod[Prime[i]*Prime[j],p]==1, c++],{i,1,j-1}], {j,1,n-1} ]; c]; Array[a, 100] (* _Amiram Eldar_, Jan 11 2019 *)

%o (Python)

%o from sympy import nextprime

%o A321005_list, plist = [], [2]

%o for n in range(10000):

%o c, p = 0, plist[-1]

%o for j in range(n):

%o pj = plist[j]

%o for i in range(j):

%o if (plist[i]*pj) % p == 1:

%o c += 1

%o A321005_list.append(c)

%o plist.append(nextprime(p)) # _Chai Wah Wu_, Nov 02 2018

%o (PARI) a(n) = my(vp = primes(n)); sum(i=1, n-2, sum(j=i+1, n-1, (vp[i]*vp[j] % vp[n]) == 1)); \\ _Michel Marcus_, Jan 11 2019

%Y A subsequence of A321004.

%Y The record values produce A321006.

%Y A320276 is an inverse.

%K nonn

%O 1,8

%A _N. J. A. Sloane_, Nov 02 2018, following a suggestion from _David James Sycamore_