11
$\begingroup$

This is less a question and more asking if someone has implemented this already, with more skill.

I need to perform the Outer-like generalized outer product of a list of lists (also a form of Tuples). I need to do it in a lazy way because the lists will become very large (many tens of thousands of elements).

I am using a Module to hold state as iterate over the products

loopOver[list0_]:=Module[
  {len0,index,get,increment,isDoneQ,doneB},
  len0=Length[list0];
  Do[index[j]=1,{j,len0}];
  doneB=False;
  get[]:=Table[list0[[j,index[j]]],{j,len0}];
  increment[]:=Block[
    {stack,p},
    If[doneB,Return[]];
    stack={len0};
    While[Length[stack]>0,
     p=First@stack;
     stack=Rest@stack;
     index[p]=index[p]+1;
     If[index[p]>Length[list0[[p]]],
      (
       If[p===1,doneB=True;Return[]];
       index[p]=1;stack=Append[stack,p-1]
       )
      ]
     ]
    ];
  isDoneQ[]:=doneB;
  {
   "get"->get,
   "increment"->increment,
   "isDoneQ"->isDoneQ
   }
  ]

And you use it as such:

Block[
 {a, get, increment, isDoneQ},
 a = {
   {"11", "12", "13"},
   {"21", "22"},
   {"31"},
   {"41", "42"}
   };
 {get,increment,isDoneQ}=loopOver[a][[All,2]];
 While[!isDoneQ[],
  Print[get[]];
  increment[]
  ]
 ]

outputting the expected outer of 3*2*1*2=12 products:

{11,21,31,41}
{11,21,31,42}
{11,22,31,41}
{11,22,31,42}
{12,21,31,41}
{12,21,31,42}
{12,22,31,41}
{12,22,31,42}
{13,21,31,41}
{13,21,31,42}
{13,22,31,41}
{13,22,31,42}

Code review would be appreciated too. I hope my code is high on self documentation where it is low on performance and use of MMA functional coding style.

$\endgroup$
  • 1
    $\begingroup$ You may find this answer of mine relevant. There, I have illustrated lazy Tuples, from which lazy Outer is pretty easy to get. Unfortunately, that answer relies on undocumented functionality, so I can't at the moment recommend it for anything more than an illustration, since it is not based on officially supported features. $\endgroup$ – Leonid Shifrin Mar 9 '16 at 23:51
  • $\begingroup$ @LeonidShifrin that question is already favorited, but not completely mentally digested ;) $\endgroup$ – Manuel --Moe-- G Mar 9 '16 at 23:55
  • $\begingroup$ Related: (9554) $\endgroup$ – Mr.Wizard Aug 12 '17 at 2:23
13
$\begingroup$

The implementation of lazy tuples here pretty much contains the solution to the lazy Outer problem. I will take the relevant parts from that code.

The following code constructs a function take, which would, given the start and end positions in the flat list of the resulting combinations, extract the elements:

ClearAll[next];
next[{left_, _}, dim_] := 
  {left - dim*(# - 1), #} &[IntegerPart[(left - 1)/dim] + 1];

ClearAll[multiDims];
multiDims[dims_] := Rest @ Reverse @ FoldList[Times, 1, Reverse @ dims];

ClearAll[multiIndex];
multiIndex[pos_, dims : {__Integer}] :=
   Rest@FoldList[next, {pos, 0}, multiDims@dims][[All, 2]]

ClearAll[take];
take[lists : {__List}, {start_, end_}] :=
  With[{rend = Min[end, Times @@ Map[Length, lists]]},
    Transpose @ MapThread[
        Part, 
        {lists, multiIndex[Range[start, rend], Length /@ lists]}
    ]
  ];

For example,

take[{{1, 2, 3}, {4, 5, 6}}, {3, 7}] == Tuples[{{1, 2, 3}, {4, 5, 6}}][[3 ;; 7]]

(* True *)

The difference is of course, that take only computes those elements that have been requested, so can be used as a basis for a lazy implementation.

Here is then an implementation of an iterator, that would return consecutive combinations in chunks of specified length:

ClearAll[makeTupleIterator];
makeTupleIterator[lists:{__List}, chunkSize_Integer?Positive]:=
  With[{len=Times @@ Length /@ lists},
    Module[{ctr = 0},          
      If[ctr >= len,
        {},
        (* else *)
        With[{taken = take[lists,{ctr+1, Min[ctr+chunkSize,len]}]},
          ctr += Length[taken];
          taken
        ]
      ]&
    ]
  ];

Here is an example: we construct an iterator with the chunk size of 10 elements:

iter = makeTupleIterator[{{"11", "12", "13"}, {"21", "22"}, {"31"}, {"41", "42"}}, 10];

Now we use it:

iter[]

(*
{        
    {"11","21","31","41"},
    {"11","21","31","42"},
    {"11","22","31","41"},
    {"11","22","31","42"},
    {"12","21","31","41"},
    {"12","21","31","42"},
    {"12","22","31","41"},
    {"12","22","31","42"},
    {"13","21","31","41"},
    {"13","21","31","42"}
}
*)

iter[]

(* {{"13", "22", "31", "41"}, {"13", "22", "31", "42"}} *)

iter[]

(* {} *)

When we get an empty list, this tells us that the iterator has been exhausted.

This basically implements lazy tuples, and therefore also lazy Outer, more or less. You gain efficiency by picking large enough chunks, since chunk extraction (take function) is pretty fast, compared to the top-level iteration that would be needed to extract element by element.

$\endgroup$
  • $\begingroup$ Hi, Leonid. Your lazy tuple is so useful. The only problem is the performace is somewhat too slow. for example, let list=Table[Range[100], 3], Tuples[list] only takes 0.03 sec. However, to generate full tuple using iter[] takes totally 3 sec $\endgroup$ – matheorem Aug 28 '16 at 15:10
  • $\begingroup$ Hi @matheorem. Unfortunately I don't have the time to test right now, but the speed depends also on the chunk size. For a better speed, the chunk size should be large enough. Try chunk size of about 10000 or more, and see if that helps. That said, the built-in version will always be much faster. The point of lazy version is that it does not need to realize all tuples at once. $\endgroup$ – Leonid Shifrin Aug 28 '16 at 16:14
  • $\begingroup$ I have tested it, I can see there is a limit of chunk size after which further increase chunk size will not reduce timing any more. Though lazy version is destined to be slower than built-in version, but more than hundreds times slower than Tuples make it not unpractical to be severed as a general approach to circumvent memory constraints if we have to traverse all the tuples. Specific method have to be used, for example, my post mathematica.stackexchange.com/q/124984/4742 wuyingddg provides a solution which is so hard to beat. But this kind of specific solution needs much more wit $\endgroup$ – matheorem Aug 29 '16 at 0:13
  • $\begingroup$ @matheorem Perhaps you are right. I will have to look into it myself to see. There may be faster methods than the one I suggested above. I will look into it when I get a moment, I just don't know when that would be. $\endgroup$ – Leonid Shifrin Aug 29 '16 at 0:58
9
$\begingroup$

Overview

Here is a refinement of @Leonid's approach that is a bit faster. The basic idea is to create a TuplesFunction that encapsulates the tuples information, and can be applied to Part type specs. That is, instead of using:

Tuples[lists][[part]]

you would use:

tf = TuplesFunction[lists];

tf[part]

The goal is to make TuplesFunction[__] as small as possible, and TuplesFunction[__][part] as fast as possible.

Decompose

As described in other answers, it is possible to make a function that takes an integer index and the lengths of the lists, and returns the indices of the lists needed to construct the corresponding tuple. The built-in functions that can do this are NumberDecompose and IntegerDigits with MixedRadix. However, since this is a function that accepts an integer (the index) and a list of integers (the lengths of the lists), and returns a list of integers (the indices of each list needed to construct the tuple), we can use Compile:

decompose = Compile[{{n, _Integer}, {d, _Integer, 1}},
    Module[{c=n, q},
        Table[
            q = Quotient[c, i];
            c = Mod[c, i];
            q,
            {i, d}
        ]
    ],
    RuntimeAttributes->{Listable}
];

The function decompose is modeled after NumberDecompose hence the second argument should be a basis, as described in the NumberDecompose documentation. Here is a check between decompose and NumberDecompose:

decompose[103, {49, 7, 1}]
NumberDecompose[103, {49, 7, 1}]

{2, 0, 5}

{2, 0, 5}

and a little speed comparison:

r1 = decompose[Range[10,100], {49, 7, 1}];//RepeatedTiming
r2 = NumberDecompose[#, {49, 7, 1}]& /@ Range[10, 100]; //RepeatedTiming

r1 === r2

{0.000029, Null}

{0.0019, Null}

True

So, compilation definitely helps here. Now, rather than relying on the Listable attribute to get decompose to work with lists of integers, it is possible to create a compiled function that accepts lists:

decomposeList = Compile[{{n, _Integer, 1}, {d, _Integer, 1}},
    Module[{c=n, q},
        Table[
            q = Quotient[c, i];
            c = Mod[c, i];
            q,
            {i, d}
        ]
    ]
];

Let's compare the two compiled functions:

r1 = decompose[Range[1000], {10000, 100, 1}]; //RepeatedTiming
r2 = decomposeList[Range[1000], {10000, 100, 1}]; //RepeatedTiming

r1 === Transpose[r2]

{0.000089, Null}

{0.000057, Null}

True

Even faster, although the returned result is transposed. This actually turns out to be a good thing.

Index decomposition to tuple

Next, we need to convert the list of indices into a tuple. In the following examples I use the lists {Range[10], Range[5], Range[7]} so that tuples extraction is obvious. For a single list of indices, we can use MapThread:

MapThread[Part, {{Range[10], Range[5], Range[7]}, {2, 4, 3}}]

{2, 4, 3}

If we have multiple lists of indices, we can use MapThread again, but this time we need to transpose first:

Transpose @ MapThread[
    Part,
    {
        {Range[10], Range[5], Range[7]},
        Transpose[{{2,4,3}, {5,2,3}, {7,1,2}}]
    }
]

{{2, 4, 3}, {5, 2, 3}, {7, 1, 2}}

Span support

One final enhancement. It would be nice to use Span in the definition of TuplesFunction. To do this, we need a way to convert a Span specification to a list of indices. Here is a function to do this:

toList[Span[a_, b_, c_:1], max_] := With[
    {
    x = Replace[a, {All->1, UpTo[x_]:>Min[x,max]}],
    y = Replace[b, {All->max, -1->max, UpTo[x_]:>Min[x,max]}]
    },

    Range[x, y, Replace[c, {All -> If[x<=y, 1, -1], Except[_Integer]->1}]]
]

A few examples:

toList[1 ;; ;; 2, 10]
toList[UpTo[13] ;; -1, 20]
toList[UpTo[23] ;; UpTo[8], 20]
toList[UpTo[23] ;; UpTo[8] ;; All, 20]

{1, 3, 5, 7, 9}

{13, 14, 15, 16, 17, 18, 19, 20}

{}

{20, 19, 18, 17, 16, 15, 14, 13, 12, 11, 10, 9, 8}

TuplesFunction

Those are the pieces we need to use in the TuplesFunction definition. The following code block is complete. It has all the definitions for TuplesFunction as well as a summary box format and the earlier definitions for decompose, decomposeList and toList;

TuplesFunction[lists_] := With[{lens = Length /@ lists},
    TuplesFunction[
        lists,
        Reverse @ FoldList[Times, 1, Reverse @ Rest @ lens]
    ]
]

TuplesFunction[lists_, basis_][index_Integer] := With[
    {decomp = 1 + decompose[index-1, basis]},

    MapThread[Part, {lists, decomp}]
]

TuplesFunction[lists_, basis_][indices:{__Integer}] := With[
    {decomp = 1 + decomposeList[indices-1, basis]},

    Transpose @ MapThread[Part, {lists, decomp}]
]

TuplesFunction[lists_, basis_][span_Span] := With[
    {r = toList[span, Times @@ Length /@ lists]},

    TuplesFunction[lists, basis][r]
]

MakeBoxes[t:TuplesFunction[lists_List, ___], StandardForm] ^:= Module[
    {lens=Length/@lists},

    BoxForm`ArrangeSummaryBox[
        TuplesFunction,
        t,
        BarChart[lens,ImageSize->30,Axes->False],
        {
            BoxForm`MakeSummaryItem[{"Count: ",Times@@lens}, StandardForm],
            BoxForm`MakeSummaryItem[{"Length: ", lens}, StandardForm]
        },
        {},
        StandardForm,
        "Interpretable"->True
    ]
]

decompose = Compile[{{n, _Integer}, {d, _Integer, 1}},
    Module[{c=n, q},
        Table[
            q = Quotient[c, i];
            c = Mod[c, i];
            q,
            {i, d}
        ]
    ],
    RuntimeAttributes->{Listable}
];

decomposeList = Compile[{{n, _Integer, 1}, {d, _Integer, 1}},
    Module[{c=n, q},
        Table[
            q = Quotient[c, i];
            c = Mod[c, i];
            q,
            {i, d}
        ]
    ]
];

toList[Span[a_, b_, c_:1], max_] := With[
    {
    x = Replace[a, {All->1, UpTo[x_]:>Min[x,max]}],
    y = Replace[b, {All->max, -1->max, UpTo[x_]:>Min[x,max]}]
    },

    Range[x, y, Replace[c, {All -> If[x<=y, 1, -1], Except[_Integer]->1}]]
]

Here is an example:

tf = TuplesFunction[{Range[4], Range[2], Range[3]}]

r1 = Tuples[{Range[4], Range[2], Range[3]}][[5 ;; 15]];
r2 = tf[5 ;; 15];

r1 === r2

TuplesFunction[{{1, 2, 3, 4}, {1, 2}, {1, 2, 3}}, {6, 3, 1}]

True

Here is a timing comparison between using Tuples and TuplesFunction:

r1 = Tuples[{Range[100], Range[100], Range[100]}]; //RepeatedTiming

tf = TuplesFunction[{Range[100], Range[100], Range[100]}]
r2 = tf /@ Partition[Range[10^6], 1000]; //RepeatedTiming

r1 === Flatten[r2, 1]

{0.0066, Null}

TuplesFunction[CompressedData[" 1:eJzt0bVCQgEAQNFnYwN2oIKY2N1gYyt2Kzrr/2+cxdkfeMMZ7nyTxd/CT3kQ BBV8U8ZfV1JFNTVEqKWOehpopIlmosSI00IrbbTTQSdddNNDLwn66GeAJCkG STPEMCOMMsY4GSaYZIppZphljnkWWGSJZVZYZY11Nthkiyw5ttlhlz32OeCQ PEccc8IpZ5xzwSVXFLjmhlvuuOeBR5545oVX3njng0+KfIU/wh/hj39/lADu ED1K "], {10000, 100, 1}]

{0.11, Null}

True

A little over an order of magnitude slower than just using Tuples, but the memory footprint is far less. Here is a memory comparison:

r1 = Total @ Tuples[{Range[100], Range[100], Range[100]}]; //MaxMemoryUsed

r2 = Sum[Total @ tf[ Span[1000l + 1, 1000(l+1)] ], {l, 0, 999}]; //MaxMemoryUsed

r1 === r2

24003312

253096

True

One final example:

tf = TuplesFunction[{Range[10^4], Range[10^4], Range[10^4]}];

tf[10^9 -10 ;; 10^9 + 10] //AbsoluteTiming

{0.000086, {{10, 10000, 9990}, {10, 10000, 9991}, {10, 10000, 9992}, {10, 10000, 9993}, {10, 10000, 9994}, {10, 10000, 9995}, {10, 10000, 9996}, {10, 10000, 9997}, {10, 10000, 9998}, {10, 10000, 9999}, {10, 10000, 10000}, {11, 1, 1}, {11, 1, 2}, {11, 1, 3}, {11, 1, 4}, {11, 1, 5}, {11, 1, 6}, {11, 1, 7}, {11, 1, 8}, {11, 1, 9}, {11, 1, 10}}}

$\endgroup$
2
$\begingroup$

Also you can use this simple "trick" :

given you example list :

a = {{"11", "12", "13"}, {"21", "22"}, {"31"}, {"41", "42"}};

the corresponding 12 (3x2x1x2) unique combinations correspond also to the unique 12combinations of the elements position in their respective list:

enter image description here

(for example here above, for the entry #9, 13 is at position 3in the {"11", "12", "13"} list, 21 is at position 1 in its list, ...).

and, these combinations of the positions is another way to write the integers [1-12] in the base form {3,2,1,2} (* Length/@a *).

For example for entry #9:

IntegerDigits[9 - 1, MixedRadix[Length /@ a]] + 1

{3, 1, 1, 1}

In other words, given the unique integer entry, you get directly the corresponding combination of the elements.

Then we can write the function you need,

take[alist_, {start_, end_}] := 
 Table[PadLeft[IntegerDigits[n, MixedRadix[Length /@ alist]], 
     Length@alist], {n, start - 1, end - 1}] // 
   Map[Thread[List[Range@Length@alist, 1 + #]] &] // 
  Apply[Part[alist, ##] &, #, {2}] &

take[alist_] := take[alist, {1, Times @@ Length /@ alist}]

and for example :

take[a, {4, 6}] // TableForm
{
 {"11", "22", "31", "42"},
 {"12", "21", "31", "41"},
 {"12", "21", "31", "42"}
}
take[{{1, 2, 3}, {4, 5, 6}}, {3, 7}] == Tuples[{{1, 2, 3}, {4, 5, 6}}][[3 ;; 7]]

True

$\endgroup$
0
$\begingroup$

Outer itself seems to be able to do this:

Flatten[Outer[List, {"11", "12", "13"}, {"21", "22"}, {"31"}, {"41", "42"}], 3]

{{"11", "21", "31", "41"}, {"11", "21", "31", "42"}, {"11", "22", "31", "41"}, 
 {"11", "22", "31", "42"}, {"12", "21", "31", "41"}, {"12", "21", "31", "42"}, 
 {"12", "22", "31", "41"}, {"12", "22", "31", "42"}, {"13", "21", "31", "41"}, 
 {"13", "21", "31", "42"}, {"13", "22", "31", "41"}, {"13", "22", "31", "42"}}
$\endgroup$
  • $\begingroup$ Yes, I am after a version of Outer where the whole list of products doesn't have to be in memory all at once, so it is done lazily. This question mathematica.stackexchange.com/questions/85278/… does a good job of describing the problem. $\endgroup$ – Manuel --Moe-- G Mar 9 '16 at 23:57
  • 2
    $\begingroup$ Tuples[] can alternatively be used in the non-lazy case. $\endgroup$ – J. M. is away Mar 10 '16 at 0:04

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.