%I #21 Sep 08 2022 08:45:04
%S 3,5,6,10,11,17,18,20,26,29,30,38,41,44,45,46,51,55,56,59,71,80,85,88,
%T 90,98,101,105,107,114,116,118,126,132,137,140,141,145,149,150,152,
%U 153,155,158,160,161,177,178,179,185,188,191,197,203,206,207,209,212,227
%N Numbers k such that usigma(k) +-1 are twin primes, where usigma(k) is the sum of the unitary divisors of k (A034448).
%H Amiram Eldar, <a href="/A065873/b065873.txt">Table of n, a(n) for n = 1..10000</a> (terms 1..1000 from Harry J. Smith)
%e 3 is a term since usigma(3) = 4 and (4-1, 4+1) = (3, 5) are twin primes.
%t f[n_] := Block[ {a = FactorInteger[n], k = l = s = 1}, l = Length[a]; While[k <= l, s = s * (a[[k, 1]]^a[[k, 2]] + 1); k++ ]; Return[s]]; Select[ Range[250], PrimeQ[ f[ # ] + 1] && PrimeQ[ f[ # ] - 1] & ]
%o (PARI) usigma(n)= { local(f,s=1); f=factor(n); for(i=1, matsize(f)[1], s*=1 + f[i, 1]^f[i, 2]); return(s) }
%o { n=0; for (m=1, 10^9, u=usigma(m); if (isprime(u - 1) && isprime(u + 1), write("b065873.txt", n++, " ", m); if (n==1000, return)) ) } \\ _Harry J. Smith_, Nov 02 2009
%o (Magma) usigma:=func<n|&+[d:d in Divisors(n)| Gcd(d,n div d) eq 1]>; [k:k in [1..230]| IsPrime(a) and (NextPrime(a)-2 eq a) where a is usigma(k)-1]; // _Marius A. Burtea_, Jan 16 2020
%Y Cf. A001097, A001359, A014574, A006512, A034448.
%K nonn
%O 1,1
%A _Robert G. Wilson v_, Dec 07 2001
|