I have a tree-based method that has the right asymptotics but a very high coefficient. The upshot being, it will not compete with other methods until we get past 10^6 or so in list size. With considerable work that tree structure could be flattened so that Compile
might be brought into play.
The basic tree layout is {left subtree, node, right subtree} where a node is a 4-tuple of the form {value, size of left subtree, size of right subtree, present}. That last is a boolean flag indicating whether the node is still around or has been removed.
{leftsubtree, node, rightsubtree} = Range[3];
emptyTree = {};
leftChild[tree_] := tree[[leftsubtree]]
SetAttributes[setLeftChild, HoldFirst]
setLeftChild[tree_, left_] := tree[[leftsubtree]] = left
rightChild[tree_] := tree[[rightsubtree]]
SetAttributes[setRightChild, HoldFirst]
setRightChild[tree_, right_] := tree[[rightsubtree]] = right
value[tree_] := tree[[node, 1]]
numLeft[emptyTree] := 0
numLeft[tree_] := tree[[node, 2]]
SetAttributes[decNumLeft, HoldFirst]
decNumLeft[node_] := node[[2]]--
numRight[emptyTree] := 0
numRight[tree_] := tree[[node, 3]]
SetAttributes[decNumRight, HoldFirst]
decNumRight[node_] := node[[3]]--
isPresent[tree_] := tree[[node, 4]]
SetAttributes[setIsPresent, HoldFirst]
setIsPresent[node_, flag_] := node[[4]] = flag
makeNode[val_] := {val, 1, 0, True};
makeNode[val_, numl_, numr_] := {val, numl, numr, True}
makeTree[val_] := {emptyTree, makeNode[val], emptyTree}
makeTree[val_, left_, right_] := {left,
makeNode[val, numLeft[left] + numRight[left] + 1,
numLeft[right] + numRight[right]], right}
newTree[tlist_List] := Block[
{len, len2, left, right},
len = Length[tlist];
If[len == 0, Return[emptyTree]];
If[len == 1, Return[makeTree[tlist[[1]]]]];
len2 = Quotient[len, 2];
left = newTree[Take[tlist, len2]];
right = newTree[Drop[tlist, len2 + 1]];
makeTree[tlist[[len2 + 1]], left, right]
]
SetAttributes[removeElement, HoldFirst];
removeElement[tree_, elem_Integer] := Module[{newnode, val},
newnode = tree[[node]];
If[isPresent[tree] && numLeft[tree] == elem, val = value[tree];
setIsPresent[newnode, False];
decNumLeft[newnode];
tree[[node]] = newnode;
Return[val];];
If[numLeft[tree] < elem, decNumRight[newnode];
tree[[node]] = newnode;
removeElement[tree[[rightsubtree]], elem - numLeft[tree]]
,(*else*)
decNumLeft[newnode];
tree[[node]] = newnode;
removeElement[tree[[leftsubtree]], elem]]]
Here is an example. For the "take" list we have a max jump size of 100. To be safe we set the final 101 elements equal to 1.
n = 10^4;
maxJump = 100;
originals = Range[n];
takelist =
Join[RandomInteger[{1, maxJump}, n - maxJump - 1],
ConstantArray[1, maxJump + 1]];
AbsoluteTiming[tree = newTree[originals];]
AbsoluteTiming[
result = Table[removeElement[tree, takelist[[j]]], {j, n}];]
(* Out[66]= {0.110308, Null}
Out[67]= {2.858037, Null} *)
This is, suffice it to say, hugely slower than the other solutions presented thus far. But the complexity is O(n log n) and it does in fact cross over with the originally presented method somewhere between one and two million.
Edit 1
Since we do not need an on line (that is, dynamically changing) structure, it is faster to use a "flat" tree. We have basically the same tree as I used above but flattened into a tensor form.
Here is a reference implementation. I simplified the code somewhat by removing things I wasn't using in the first place.
makeNode[val_] := {val, 1,(*0,*)1};
makeNode[val_, numl_] := {val, numl, 1}
makeTree[val_] := {makeNode[val]}
makeTree[val_, left_, right_] := Module[{node},
node = makeNode[val, Length[left] + 1];
Join[left, {node}, right]]
newTree[tlist_List] :=
Block[{len, len2, left, right}, len = Length[tlist];
If[len == 0, Return[{}]];
If[len == 1, Return[makeTree[tlist[[1]]]]];
len2 = Quotient[len, 2];
left = newTree[Take[tlist, len2]];
right = newTree[Drop[tlist, len2 + 1]];
makeTree[tlist[[len2 + 1]], left, right]]
SetAttributes[removeElement, HoldFirst];
removeElement[tree_, elem_Integer] := Catch[Module[
{node, mid = Quotient[Length[tree], 2] + 1, start = 0,
end = Length[tree] + 1, j = 0, sum = 0},
While[True && j < 20, j++;
If[mid == start || mid == end, Throw[-1]];
node = tree[[mid]];
If[node[[3]] == 1 && sum + node[[2]] == elem,
tree[[mid, 3]] = 0;
tree[[mid, 2]]--;
Throw[node[[1]]];];
If[sum + node[[2]] < elem,
sum += node[[2]];
{start, mid} = {mid, mid + Quotient[(end - 1 - mid), 2] + 1};
,(*else*)
tree[[mid, 2]]--;
{end, mid} = {mid, start + Quotient[(mid - start - 1), 2] + 1};];
]
]]
To be continued.
Edit 2
We can now put this through Compile
as below.
newTreeC = Compile[{{tlist, _Integer, 1}},
Catch[Module[{len, len2, left, right},
len = Length[tlist];
If[len == 1, Throw[{{tlist[[1]], 1, 1}}]];
If[len == 2, Throw[{{tlist[[1]], 1, 1}, {tlist[[2]], 2, 1}}]];
len2 = Quotient[len, 2];
left = newTreeC[tlist[[1 ;; len2]]];
right = newTreeC[tlist[[len2 + 2 ;; -1]]];
Join[left, {{tlist[[len2 + 1]], Length[left] + 1, 1}},
right]]], {{newTreeC[_], _Integer, 2}},
RuntimeOptions -> "Speed"];
takeByListC = Compile[{{otree, _Integer, 2}, {tlist, _Integer, 1}},
Module[
{tree = otree, node, mid, start, end, sum, elem},
Table[elem = tlist[[k]];
mid = Quotient[Length[otree], 2] + 1;
start = 0;
end = Length[otree] + 1;
sum = 0;
Catch[While[True,
If[mid == start || mid == end, Throw[-1]];
node = tree[[mid]];
If[node[[3]] == 1 && sum + node[[2]] == elem,
tree[[mid, 3]] = 0;
tree[[mid, 2]]--;
Throw[node[[1]]],
If[sum + node[[2]] < elem,
sum += node[[2]];
{start, mid} = {mid, mid + Quotient[(end - 1 - mid), 2] + 1};
,(*else*)
tree[[mid, 2]]--;
{end, mid} = {mid,
start + Quotient[(mid - start - 1), 2] + 1}];
]]]
, {k, Length[tlist]}]], RuntimeOptions -> "Speed"];
This test, on a 6 year old laptop, will at least give an idea of performance.
I think I may finally have all the bugs out but caveat emptor.
n = 10^6;
maxJump = 100;
originals = Range[n];
takelist =
Join[RandomInteger[{1, maxJump}, n - maxJump - 1],
ConstantArray[1, maxJump + 1]];
AbsoluteTiming[tree = newTreeC[originals];]
AbsoluteTiming[result = takeByListC[tree, takelist];]
(* Out[45]= {5.290000, Null}
Out[46]= {18.920000, Null} *)
Might be faster with CompilationTarget -> "C"
. If I had a machine that knew what a C compiler looked like.
Edit 3
[I guess I'm spending too much time on this. Nice problem it has turned out to be.]
Here is an improved compiled version. Still around twice as slow as the code by @2012rcampion.
newTreeC = Compile[{{len, _Integer}}, Catch[Block[
{len2, left, right},
If[len == 1, Throw[{{1, 1}}]];
If[len == 2, Throw[{{1, 1}, {2, 1}}]];
len2 = Quotient[len, 2];
left = newTreeC[len2];
right = newTreeC[len - len2 - 1];
Join[left, {{Length[left] + 1, 1}}, right]]],
{{newTreeC[_], _Integer, 2}},
RuntimeOptions -> "Speed", CompilationTarget -> "C"];
takeByListC =
Compile[{{otree, _Integer, 2}, {tlist, _Integer, 1}},
Block[{tree = otree, node, mid, start, end, sum, elem,
saveend = Length[otree] + 1,
savemid = Floor[Length[otree]/2.] + 1},
Table[
elem = tlist[[k]];
mid = savemid;
start = 0;
end = saveend;
sum = 0;
Catch[While[True,
node = tree[[mid]];
If[node[[2]] == 1 && sum + node[[1]] == elem,
tree[[mid, 2]] = 0;
tree[[mid, 1]]--;
Throw[mid]
];
If[sum + node[[1]] < elem,
sum += node[[1]];
start = mid;
mid = mid + Floor[(end - 1 - mid)/2.] + 1;
,(*else*)
tree[[mid, 1]]--;
end = mid;
mid = start + Floor[(mid - start - 1)/2.] + 1;
];
]]
, {k, Length[tlist]}]],
RuntimeOptions -> "Speed", CompilationTarget -> "C"];
Here is code similar to the method of @halirutan, but moving result values instead of index list values. It is slightly longer than need be so as to avoid using RotateRight
. That usage would make is somewhat slower; I guess it involves calls to the run time library.
iterativeTakeByList = Compile[{{tlist, _Integer, 1}}, Block[
{len = Length[tlist], result, elem},
result = Range[len];
Do[
elem = result[[j + tlist[[j]] - 1]];
Do[result[[j + k + 1]] = result[[j + k]],
{k, tlist[[j]] - 2, 0, -1}];
result[[j]] = elem;
, {j, len}
];
result
],
RuntimeOptions -> "Speed", CompilationTarget -> "C"];
Here is a test on a list of a million elements. With a maximum jump of 1000 it is still competitive with the best of the tree-based codes.
n = 10^6;
maxJump = 1000;
takelist =
Join[RandomInteger[{1, maxJump}, n - maxJump - 1],
ConstantArray[1, maxJump + 1]];
AbsoluteTiming[result = iterativeTakeByList[takelist];]
(* Out[179]= {0.762602, Null} *)
Delete
orDrop
at each step. They both require shifting or copying elements, so the algorithms are $O(n^2)$. I think there is a method using a tree structure that's $O(n\log n)$, I'm going to sleep on it and see if my USACO training comes back to me =) $\endgroup$