Partitioning a number into consecutive integers
The sum of consecutive numbers from $a$ to $b$ is
$$\frac{(a+b)(b-a+1)}{2}$$
hence simply
f[n_] := {a, b} /.
Solve[(a + b) (b - a + 1)/2 == n && 0 < a < n && 0 < b < n, {a, b}, Integers]
f[45] // AbsoluteTiming
{0.019466, {{1, 9}, {5, 10}, {7, 11}, {14, 16}, {22, 23}}}
It is straightforward and rather fast. As a test case:
f[4500] // AbsoluteTiming
{0.063403, {{23, 97}, {27, 98}, {78, 122}, {93, 132}, {168, 192}, {176, 199}, {293, 307}, {496, 504}, {559, 566}, {898, 902}, {1499, 1501}}}
This distribution is interesting:
ListPlot[tab = Table[{i, Length @ f[i]}, {i, 1000}], Filling -> Axis,
Frame -> True, PlotRange -> All,
FrameLabel -> {"n", "Length@f[n]"}]; // AbsoluteTiming
{12.3098, Null}
The biggest number of partitions (Length@f[n] == 15
) for n <= 1000
is for n = 945
:
Select[tab, #[[2]] == Max @ tab[[All,2]] &]
{{945, 15}}
while Length@f[944] == 1
and Length@f[946] == 3
.
It works also with huge numbers:
f[10^11] // AbsoluteTiming
{0.149738, {{60688, 451312}, {925363, 1027762}, {1240938, 1319062}, {4872573, 4893052}, {6392188, 6407812}, {24412015, 24416110}, {31998438, 32001562}, {159999688, 160000312}, {799999938, 800000062}, {3999999988, 4000000012}, {19999999998, 20000000002}}}
f[10^40] // AbsoluteTiming // Short
{1.55024,{{39445166261547851563,146819348661547851562},<<38>>,{<<1>>}}}
(I'm surprised it works so well, and I think that the timing is quite acceptable too.)
A sample rejection method (brute force, with just a little trick to use Ceiling[n/2]
):
n = 45;
MinMax /@
First /@ Select[
Flatten[#, 1] & @
Table[{Range[i, j], Total @ Range[i, j]}, {j, 1, Ceiling[n/2]}, {i, 1, j}],
#[[2]] == n &]; // AbsoluteTiming
{0.000803, Null}
works well for this simple case (however, on my laptop it needs 14 sec for n = 4000
) , but when
n = 4500;
it crashes the kernel.
Additionally,
g[n_] := Cases[FrobeniusSolve[Range @ n, n], {0 ..., 1 .., 0 ..}, Infinity]
is very slow:
(o = g[45]); // AbsoluteTiming
{43.871, Null}
although works:
MinMax /@ (Pick[Range@45, #, 1] & /@ o) // Sort
{{1, 9}, {5, 10}, {7, 11}, {14, 16}, {22, 23}}
I took it as a challenge to avoid using Solve
, which can be slower than a direct assault. If $a$ is the first number in the sum of consecutive positive integers, and $k$ is the count of integers summing to $n$, then $n=k*a+k(k-1)/2$. Solve this for $a=n/k-(k-1)/2$, with bounds $1 \le k \le {\rm Floor}[(\sqrt{8n+1}-1)/2]$. Consider the odd and even divisors of $n$ and $2n$, respectively, to give the following. This includes the $k=1$ case of the sum of one number, just $n$.
CoreyPartition[n_] :=
Block[{bound = Floor[(Sqrt[1 + 8 n] - 1)/2], oddk, evenk, k, e},
oddk = Pick[#,OddQ[#]] &[Pick[#, UnitStep[bound - #], 1]&[Divisors[n]]];
evenk = Pick[#, UnitStep[bound - #], 1] &[Divisors[2 n]];
e = IntegerExponent[n, 2];
evenk = Pick[evenk, IntegerExponent[evenk, 2], 1 + e];
k = Reverse[Sort[Join[oddk, evenk]]];
Transpose[{n/k - (k - 1)/2, n/k + (k - 1)/2}]]
A slight modification to f[n]
for cases with no solution, such as $n=1,2,4,8,\ldots$.
f[n_] := If[# == {}, {}, {a, b} /. #]&[
Solve[(a + b) (b - a + 1)/2 == n && 0 < a < n && 0 < b < n,
{a, b}, Integers]]
The two solutions agree, for example:
Sum[Total[Most[CoreyPartition[n]] - f[n]], {n, 1, 200}]
However, testing with large integers such as n=10^40
showed CoreyPartition[n]
was over 200 times faster.
Update
Given DivisorPairs
from Mr.Wizard:
DivisorPairs[n_] :=
Thread[{#, Reverse[#]}][[ ;; Ceiling[Length[#]/2]]] &[Divisors[n]]
there is the following one-liner which is twice as fast as CoreyPartition
.
CoreyFastPartition[n_] :=
Reverse[Pick[#, Total[Mod[#, 4], {2}], 2]] &[ DivisorPairs[8 n]] /.
{r_Integer, s_Integer} -> {(s - r + 2)/4, (s + r - 2)/4}
This code returns the trivial solution {n,n}, which may be removed with Most
.
d[n_] := With[{dv = Divisors[n]}, {#, n/#} & /@
Pick[dv, # < Sqrt[n] & /@ dv]]
f[a_, b_] := With[{p = a - 1, q = b}, {(q - p)/2, (p + q)/2}]
res[n_] := Rest@Cases[f @@@ d[2 n], {_Integer, _Integer}]
So,
ListPlot[{#, Length[res[#]]} & /@ Range[1000], Filling -> Axis]
ListPlot[DeleteCases[{#, #2 - #1 + 1 & @@ (res[#][[-1]])} & /@
Range[1000], {_, {}}], Filling -> Axis,
Epilog -> {Red, Point[{# (# + 1)/2, #} & /@ Range[45]]},
PlotLabel -> "Maximum Length of run"]
Cases[{#, res@#} & /@ Range[5000], {x_, {}} :> x]
Powers of 2 no non-trivial representation. Triangular numbers shown as red points.