A366251 - OEIS

%I #45 Oct 27 2023 08:24:14

%S 1,3,5,7,9,11,13,17,19,21,23,25,27,29,31,35,37,39,41,43,47,49,53,55,

%T 57,59,61,63,65,67,71,73,77,79,81,83,85,89,91,93,95,97,101,103,107,

%U 109,111,113,115,117,119,121,125,127,129,131,133,137,139,143,145

%N Numbers with a coprime divisor shift.

%C A number k has a coprime divisor shift s if GCD(d + s, k) = 1 for all divisors d of k.

%C If k is in the sequence, then all divisors of k are in the sequence too.

%C From _David A. Corneth_, Oct 06 2023: (Start)

%C As a consequence of the above, if some m is not in the sequence then any multiple of m is not in the sequence either. So all terms are odd as 2 is not in the sequence.

%C To see if k is a term we can do the following:

%C - Check whether k is even; if so then k is not in the sequence.

%C - Make a list of the divisors of k. Let tau(k) be the number of divisors of k. Then for each prime p <= tau(k) if the number of residues mod p of the divisors of k is equal to p then k is not in the sequence. Otherwise by the Chinese Remainder Theorem we can find a number s such that gcd(d + s, n) = 1 for all d.

%C So every odd prime is a term. (End)

%H M. Farrokhi D. G., <a href="/A366251/b366251.txt">Table of n, a(n) for n = 1..10000</a>

%e 1 is a term since GCD(1 + 0, 1) = 1.

%e 2 is not a term since GCD(1 + s, 2) > 1 or GCD(2 + s, 2) > 1 for all nonnegative integers s.

%e 3 is a term since GCD(1 + 1, 3) = GCD(3 + 1, 3) = 1.

%p aList := proc(len) local isds, findds;

%p isds := (k, s) -> andmap(d -> igcd(d + s, k) = 1, NumberTheory:-Divisors(k));

%p findds := k -> ormap(s -> isds(k, s), [seq(1..k)]);

%p select(k -> findds(k), [seq(1..len)]) end:

%p aList(133);  # _Peter Luschny_, Oct 06 2023

%p # More efficient, after _David A. Corneth_:

%p isa := proc(n) local d, p, t; p := 3;

%p if irem(n, 2) = 0 then return false fi;

%p d := NumberTheory:-Divisors(n);

%p while p < nops(d) do

%p    {seq(irem(t, p), t = d)};

%p    if nops(%) = p then return false fi;

%p    p := nextprime(p);

%p od: true end:

%p aList := len -> select(isa, [seq(1..len)]):

%p aList(145);  # _Peter Luschny_, Oct 07 2023

%t isa[n_] := Module[{d, p, t}, p = 3; If[Mod[n, 2] == 0, Return[False]]; d = Divisors[n]; While[p < Length[d], If[Length[Union@Table[Mod[t, p], {t, d}]] == p, Return[False]]; p = NextPrime[p]]; True];

%t aList[len_] := Select[Range[len], isa];

%t aList[145] (* _Jean-François Alcover_, Oct 27 2023, after _Peter Luschny_ *)

%o (PARI)

%o isds(k,s)={fordiv(k, d, if(gcd(d+s, k)<>1, return(0))); 1}

%o findds(k)={for(s=1, k, if(isds(k,s), return(s))); 0}

%o select(k->findds(k), [1..150]) \\ _Andrew Howroyd_, Oct 05 2023

%o (PARI)

%o is(n) = {

%o 	if(!bitand(n, 1), return(0));

%o 	my(d = divisors(n));

%o 	forprime(p = 3, #d,

%o 		if(#Set(d % p) == p,

%o 			return(0)

%o 		)

%o 	); 1	

%o } \\ _David A. Corneth_, Oct 06 2023

%o (GAP)

%o CoprimeDivisorShift := function(n)

%o local shift, divisors;

%o shift := 0;

%o if not IsPrimeInt(n) and First(PrimeDivisors(n), p -> CoprimeDivisorShift(n / p) = infinity) <> fail then

%o 	shift := infinity;

%o fi;

%o if shift < infinity then

%o 	divisors := DivisorsInt(n);

%o 	while shift < n and First(divisors, d -> GcdInt(d + shift, n) > 1) <> fail do

%o 		shift := shift + 1;

%o 	od;

%o 	if shift = n then

%o 		shift := infinity;

%o 	fi;

%o fi;

%o return shift;

%o end;

%Y Cf. A366219, A366330.

%K nonn,easy

%O 1,2

%A _M. Farrokhi D. G._, Oct 05 2023