Is there any way to instruct Complement to skip the sorting part? If the answer is no (likely), the next question would be: how can I remove from a sorted list another sorted list (in an efficient way)?

I would like the output of the following code to be {1,3} as in the usual Complement function:

universe = {1, 2, 3}; subst = {2} ; myComplement[universe, subst] 

But the output of the following code should be {3, 2, 1} because it stopped trying to remove terms after it realized that subst[[1]] > universe[[1]].

universe = {3, 2, 1}; subst = {2}; myComplement[universe, subst] 

This question is different from a simple use of Complement or DeleteCases because we must use the fact that the lists are (assumed) sorted for efficiency. In the example below, the list universe has $2^{25}-1$ elements.

The motivation for this question.

While answering this post, I tried to go from listing the classes of equivalence of 4x4 matrices (as requested by the original question) to 5x5 case (just to push it). My answer solved the 4x4 case in 4 seconds. But the same algorithm, I estimated, will take 3 hours for the 5x5 case.

toMatrix = Partition[IntegerDigits[#, 2, 25], 5] &; toInteger = FromDigits[Flatten[#], 2] &; allSyms = Module[{s1, s3, s4, s6}, {#, s1 = Reverse[#, {2}], Reverse[#], s3 = Transpose[#], s4 = Reverse[s3, {2}], Reverse[Transpose[s1], {2}], s6 = Reverse[Transpose[s4], {2}], Reverse[Transpose[s6], {2}]}] &; casesToCheck = Range[0, 2^25 - 1]; Timing[answer = {MatrixForm@toMatrix@First[#], Length[#]} & /@ (Reap[ NestWhile[ Complement[#, Sow[Union[toInteger /@ allSyms[toMatrix[#[[1]]]]]]] &, casesToCheck, (n = Length[#]) > 0 &]][[2, 1]]); Length[answer]] 

Obviously, there is a major bottleneck that is triggered now. I guess it is because Complement sorts the input. Instead of sorting a list of 65536 integers, the 5x5 case deals with a list of 33554431 integers.