Bug or undocumented behaviour in Intersection

Question

I have two lists consisting of element pairs:

a1 = Table[{i, i^2}, {i, 10}] // N;
a2 = Table[{2 i, Sqrt[i]}, {i, 5}] // N;

I want to pull out from both lists those pairs which have the same first element. Using Intersection with the SameTest option:

In[84]:= d1 = Intersection[a1, a2, SameTest -> (First[#1] === First[#2] &)]

Out[84]= {{2., 4.}, {4., 16.}, {6., 36.}, {8., 64.}, {10., 100.}}

gives a different result as compared to

In[85]:= d2 = Intersection[a2, a1, SameTest -> (First[#1] === First[#2] &)]

Out[85]= {{2., 1.}, {4., 1.41421}, {6., 1.73205}, {8., 2.}, {10., 2.23607}}

d1 contains the element pairs of the list a1 for which there exists also a similar first element list a2. d2 contains the element pairs of a2 for which there exists also a similar first element in list a1. Thus Intersection is not commutative.

I would like to get as "Intersection" all element pairs a1 and of a2 for which there is a match in the first element. I can of course call Intersection as shown above two times and then join the results d1 and d2. But is there not a more efficient way so that I need to call Intersection only once? I ask because in my real case the lists are very large and the SameTest is more complicated than the identity test used in this example.

Answer 1

Maybe this produces does what you want (it is basically what Manuel --Moe-- G suggested).

f[a1_?MatrixQ, a2_?MatrixQ] := Thread /@ List @@@ Normal[
    Select[
     Merge[{
       AssociationThread @@ Transpose[a1],
       AssociationThread @@ Transpose[a2]
       },
      Identity
      ],
     Length[#] >= 2 &
     ]
    ]

f[a1, a2]

{{{2., 4.}, {2., 1.}}, {{4., 16.}, {4., 1.41421}}, {{6., 36.}, {6., 1.73205}}, {{8., 64.}, {8., 2.}}, {{10., 100.}, {10., 2.23607}}}

Alternate Method

One can obtain a more efficient implementation (without binning) with Nearest:

g[a1_?MatrixQ, a2_?MatrixQ] := 
 Module[{idx, rowpointers, columnindices, A, i, j},
  idx = Nearest[a2[[All, 1]] -> "Index", a1[[All, 1]], {\[Infinity], 0}];
  rowpointers = Join[{0}, Accumulate[Length /@ idx]];
  columnindices = Partition[Join @@ idx, 1];
  A = SparseArray @@ {Automatic,
     {Length[a1], Length[a2]}
     ,
     0,
     {1,
      {rowpointers, columnindices},
      ConstantArray[1, rowpointers[[-1]]]
      }
     };
  {i, j} = Transpose[A["NonzeroPositions"]];
  Union[a1[[i]], a2[[j]]]
  ]

Now

n = 1000000;
a1 = Transpose[{Range[1., 2 n], Range[1., 2 n]^2}];
a2 = Transpose[{Range[2., 2 n, 2.], Sqrt[Range[1., n]]}];

resultf = f[a1, a2]; // AbsoluteTiming // First    
resultg = g[a1, a2]; // AbsoluteTiming // First    
Union @@ resultf == resultg

14.7094

1.26785

True

You can obtain essentially the same binning (up to reordering) with GatherBy, but that takes a bit of extra time:

resultg2 = GatherBy[resultg, First]; // AbsoluteTiming // First
Sort[Sort /@ resultg2] == Sort[Sort /@ resultf]

1.26529

True

Discussion

The second method is about 10 times faster for quite a range of values of n. Since Nearest employs a space partition tree that should have complexity $O( n \, \log(n))$ (instead of $O(n^2)$ for a naive approach), this tells me that also Merge uses some clever $O( n \, \log(n))$-method to perform the lookup. I've once read somewhere here on Mathematica.StackExchange that "Mathematica's" Association data type is based on a special kind of trie, so that is entirely possible. However, in this special case (where the keys are members of a suitable metric space), a space partition tree seems to be faster.

Is Intersection's current behavior a bug?

I think it isn't. By assigning a SameTest, you effectively tell Mathematica that it should not care about which representative to choose within an equivalence class. So Mathematica chooses the first(?) representative it can find. Listing all representatives of the equivalence class is just not what you asked for.

Answer 2

After posting the question I found out that solutions to similar problems had already been posted on this website. For my problem I constructed the following solution:

[a1_?MatrixQ, a2_?MatrixQ] := Flatten[#, 1] & /@ Select[GatherBy[Flatten[{a1, a2}, 1], First], Length[#] == 2 &]

Faster than the f-function of @Henrik, but somewhat slower than the brilliant g-function of @Henrik. Works well with matrices that have different no of columns, and it is also straightforward to implement more complicated tests for "sameness" of the first elements.

Answer 3

a = N @ Table[{i, i^2}, {i, 10}];

b = N @ Table[{2 i, Sqrt[i]}, {i, 5}];

Using KeyIntersection

c = Merge[# &] @ KeyIntersection[{Rule @@@ a, Rule @@@ b}]

<|2. -> {4., 1.}, 4. -> {16., 1.41421}, 6. -> {36., 1.73205}, 8. -> {64., 2.}, 10. -> {100., 2.23607}|>

c // Dataset

Convert to list (if wanted)

Thread /@ KeyValueMap[List] @ c

{{{2., 4.}, {2., 1.}}, {{4., 16.}, {4., 1.41421}}, {{6., 36.}, {6., 1.73205}}, {{8., 64.}, {8., 2.}}, {{10., 100.}, {10., 2.23607}}}