cp's OEIS Frontend

A278638 Numbers n such that 1/n is a difference of Egyptian fractions with all denominators < n.

%I A278638 #40 Jul 29 2020 05:08:58
%S A278638 6,12,15,18,20,21,24,28,30,33,35,36,40,42,44,45,48,52,54,55,56,60,63,
%T A278638 65,66,68,70,72,75,76,77,78,80,84,85,88,90,91,95,96,99,100,102,104,
%U A278638 105,108,110,112,114,115,117,119,120,126,130,132,133,135,136,138,140,143,144,145,147,150,152,153
%N A278638 Numbers n such that 1/n is a difference of Egyptian fractions with all denominators < n.
%C A278638 Numbers n such that we can write 1/n = Sum_{1<=k<n} c(k)/k with all c(k) in {-1,0,1}.
%C A278638 Numbers n such that A072207(n) < 2*A072207(n-1).
%C A278638 If n is in the sequence, so is k*n for all k>1 (cf. A328226).
%C A278638 Contains A001284, because 1/(m*k) = 1/(m*(k-m))-1/(k*(k-m)).
%C A278638 Disjoint from A000961.
%C A278638 2*p^k with p prime is in the sequence if and only if p=3.
%C A278638 3*p^k with p prime is in the sequence if and only if p=2,5,7 or 11.
%C A278638 4*p^k with p prime is in the sequence if and only if p=3,5,7,11,13,17 or 19.
%C A278638 For each m that is not a term, there are only finitely many primes p such that some m*p^k is a term. [Corrected by _Max Alekseyev_, Oct 08 2019]
%H A278638 Robert Israel, <a href="/A278638/b278638.txt">Table of n, a(n) for n = 1..431</a>
%H A278638 Robert Israel, <a href="/A278638/a278638_1.txt">Examples for n = 1..431</a>
%H A278638 Robert Israel, <a href="https://math.stackexchange.com/q/2029163">1/n as a difference of Egyptian fractions with all denominators < n</a>, Math StackExchange, 2017.
%e A278638 44 is in the sequence because 1/44 = (1/12 + 1/33) - 1/11.
%e A278638 4 is not in the sequence because 1/4 can't be written as the difference of sums of two subsets of {1, 1/2, 1/3}.
%p A278638 N:= 200: # to get all terms <= N
%p A278638 V:= Vector(N):
%p A278638 f:= proc(n)  option remember;
%p A278638 local F,E,p,e,k,m,L,L1,i,s,t,sg,Maybe;
%p A278638 global Rep;
%p A278638   F:= numtheory:-factorset(n);
%p A278638   if nops(F) = 1 then return false fi;
%p A278638   if ormap(m -> n < m^2 and m^2 < 2*n, numtheory:-divisors(n)) then
%p A278638     for m in numtheory:-divisors(n) do
%p A278638       if n < m^2 and m^2 < 2*n then
%p A278638         k:= n/m; Rep[n]:= [m*(k-m),-k*(k-m)]; return true
%p A278638       fi
%p A278638     od
%p A278638   fi;
%p A278638   F:= convert(F,list);
%p A278638   E:= map(p -> padic:-ordp(n,p), F);
%p A278638   i:= max[index](zip(`^`,F,E));
%p A278638   p:= F[i];
%p A278638   e:= E[i];
%p A278638   k:= n/p^e;
%p A278638   Maybe:= false;
%p A278638   for i from 3^(k-1) to 2*3^(k-1)-1 do
%p A278638     L:= (-1) +~ convert(i,base,3);
%p A278638     s:= 1/k - add(L[i]/i,i=1..k-1);
%p A278638     if numer(s) mod p = 0 then
%p A278638     Maybe:= true;
%p A278638       t:= abs(s/p^e); sg:= signum(s);
%p A278638       if  (numer(t) <= 1 and (denom(t) < n or (denom(t) < N and V[denom(t)] = 1))) or (numer(t) = 2 and denom(t) < N and V[denom(t)] = 1) then
%p A278638          L1:= subs(0=NULL, [seq(L[i]*i*p^e,i=1..k-1)]);
%p A278638          if t = 0 then ;
%p A278638          elif numer(t) = 1 and denom(t) < n then L1:= [op(L1),sg/t]
%p A278638          elif procname(2/t) then
%p A278638             L1:= ([op(L1), 2*sg/t, op(expand(sg*Rep[2/t]))])
%p A278638          else next
%p A278638          fi;
%p A278638          if max(abs~(L1)) < n then Rep[n]:= L1; return true fi;
%p A278638       fi;
%p A278638     fi
%p A278638   od:
%p A278638   if Maybe then printf("Warning: %d is uncertain\n",n)
%p A278638 else false
%p A278638 fi;
%p A278638 end proc:
%p A278638 for n from 6 to N do
%p A278638   if V[n] = 0 and f(n) then
%p A278638     V[n] := 1;
%p A278638     for j from 2*n to N by n do
%p A278638       if not assigned(Rep[j]) then
%p A278638         V[j]:= 1;
%p A278638         Rep[j] := map(`*`,Rep[n],j/n);
%p A278638         f(j):= true;
%p A278638       fi
%p A278638     od;
%p A278638   fi;
%p A278638 od:
%p A278638 select(t -> V[t]=1,[$6..N]);
%t A278638 sol[n_] := Module[{c, cc}, cc = Array[c, n-1]; FindInstance[AllTrue[cc, -1 <= # <= 1&] && 1/n == Total[cc/Range[n-1]], cc, Integers, 1]];
%t A278638 Reap[For[n = 6, n <= 200, n++, If[sol[n] != {}, Print[n]; Sow[n]]]][[2, 1]] (* _Jean-François Alcover_, Jul 29 2020 *)
%Y A278638 Cf. A000961, A001284, A072207.
%Y A278638 Contains A005279. - _Robert G. Wilson v_, Nov 27 2016
%K A278638 nonn
%O A278638 1,1
%A A278638 _Robert Israel_, Nov 24 2016