%I
%S 6,10,12,14,15,18,20,21,22,24,26,28,33,35,36,39,40,44,52,55
%N Smallest set of distinct numbers with property that the sum of their reciprocals is 1 and each number is of the form p^i*q^j where p and q are distinct primes and i and j are positive.
%C In other words, smallest set of distinct numbers in A007774 (i.e., numbers that are divisible by exactly 2 different primes) whose reciprocals sum to 1.
%C This is the only set of 20 distinct numbers with the specified property. (The sum of the reciprocals of the first 19 numbers that are divisible by exactly 2 different primes is 229926/230945, which is less than 1, so no set of fewer than 20 such numbers can have reciprocals that sum to 1.)
%C There exist 8 such sets of 21 numbers, and 218 such sets of 22 numbers.
%C For a set of 25 distinct numbers having the specified property, see A201463.
%C The MAGMA program (as of Feb 12 2014) is very fast at n=20, but very slow if n is changed to a value of 22 or larger. The main opportunity for improving its efficiency may be in the way it handles the search for the last two terms of a set, given the first n2 terms at the current point in the search tree. Specifically, if the remainder r = 1  (1/i[1] + 1/i[2] + ... + 1/i[n2]) is close to zero, e.g., 10^8, then it tests about 10^8 values of i[n1] (from about 10^8 to about 2*10^8), instead of applying some smarter approach that makes use of the factorization of the denominator of the remainder r.
%D E. J. Barbeau, Expressing one as a sum of distinct reciprocals: comments and a bibliography, Eureka (Ottawa), 3 (1977), 178181.
%D N. Burshtein, Improving solutions of Sum_{i=1..k} 1/x_i = 1 ..., Discrete Math., 306 (2006), 14381439.
%o (MAGMA) n:=20; i:=[]; j:=1; i[j]:=0; r:=1; while true do i[j]+:=1; while #Factorization(i[j]) ne 2 do i[j]+:=1; end while; if (n(j1))/i[j] lt r then if j eq 1 then "done"; break; end if; j:=1; r+:=1/i[j]; elif j eq n1 then TestIntN:=Floor(1/(r1/i[j])); if TestIntN le i[n1] then j:=1; r+:=1/i[j]; elif (r  1/i[j] eq 1/TestIntN) and (#Factorization(TestIntN) eq 2) then i[n]:=TestIntN; i; end if; else r:=1/i[j]; j+:=1; i[j]:=i[j1]; if 1/(i[j]+1) ge r then i[j]:=Floor(1/r); end if; end if; end while; // _Jon E. Schoenfield_, Feb 12 2014
%Y Cf. A007774, A201463.
%K nonn,fini,full
%O 1,1
%A _Jon E. Schoenfield_, Feb 02 2014
%E Improved MAGMA program and comments about opportunities for further improvement from _Jon E. Schoenfield_, Feb 12 2014
