%I #21 Oct 19 2024 16:17:36
%S 738,748,774,846,954,1062,1098,1206,1278,1314,1422,1494,1602,1746,
%T 1818,1854,1926,1962,2034,2286,2358,2466,2502,2682,2718,2826,2934,
%U 3006,3114,3222,3258,3438,3474,3492,3546,3582,3636,3708,3798,3852,3924,4014,4068,4086
%N Non-deficient numbers with even sigma which are not Zumkeller.
%C Numbers which are non-deficient (sigma(n) >= 2n) [A023196] such that sigma(n) [A000203] is even but which are not Zumkeller numbers [A083207], i.e., the positive factors of n cannot be partitioned into two disjoint parts so that the sums of the two parts are equal.
%H Amiram Eldar, <a href="/A171641/b171641.txt">Table of n, a(n) for n = 1..10000</a> (terms 1..1000 from Chai Wah Wu)
%H Peter Luschny, <a href="http://www.luschny.de/math/seq/ZumkellerNumbers.html"> Zumkeller Numbers</a>.
%p with(NumberTheory):
%p isA171641 := proc(n) local s, p, i, P;
%p s := SumOfDivisors(n);
%p if s::odd or s < n*2 then false else
%p P := mul(1 + x^i, i in Divisors(n));
%p 0 = coeff(P, x, s/2) fi end:
%p select(isA171641, [seq(1..4100)]); # _Peter Luschny_, Oct 19 2024
%t Reap[For[n = 2, n <= 4000, n = n+2, sigma = DivisorSigma[1, n]; If[sigma >= 2n && EvenQ[sigma] && Coefficient[ Times @@ (1 + x^Divisors[n]) // Expand, x, sigma/2] == 0, Print[n]; Sow[n]]]][[2, 1]] (* _Jean-François Alcover_, Jul 26 2013 *)
%o (Python)
%o from sympy import divisors
%o import numpy as np
%o A171641 = []
%o for n in range(2,10**6):
%o d = divisors(n)
%o s = sum(d)
%o if not s % 2 and 2*n <= s:
%o d.remove(n)
%o s2, ld = int(s/2-n), len(d)
%o z = np.zeros((ld+1,s2+1),dtype=int)
%o for i in range(1,ld+1):
%o y = min(d[i-1],s2+1)
%o z[i,range(y)] = z[i-1,range(y)]
%o z[i,range(y,s2+1)] = np.maximum(z[i-1,range(y,s2+1)],z[i-1,range(0,s2+1-y)]+y)
%o if z[ld,s2] != s2:
%o A171641.append(n)
%o # _Chai Wah Wu_, Aug 19 2014
%Y Cf. A083207, A023196.
%K nonn,changed
%O 1,1
%A _Peter Luschny_, Dec 14 2009