Cartesian product of x and y array points into single array of 2D points
>>> numpy.transpose([numpy.tile(x, len(y)), numpy.repeat(y, len(x))])
array([[1, 4],
[2, 4],
[3, 4],
[1, 5],
[2, 5],
[3, 5]])
See Using numpy to build an array of all combinations of two arrays for a general solution for computing the Cartesian product of N arrays.
Finding the cartesian product of 2d arrays of different shapes and horizontally concatenating them row by row into 1 array
How about using a mix of repeat
and tile
(which itself uses repeat
):
In [75]: >>> ident2 = np.identity(2)
...: >>> ident3 = np.identity(3)
In [76]: np.repeat(ident3,repeats=2,axis=0)
Out[76]:
array([[1., 0., 0.],
[1., 0., 0.],
[0., 1., 0.],
[0., 1., 0.],
[0., 0., 1.],
[0., 0., 1.]])
In [77]: np.tile(ident2,(3,1))
Out[77]:
array([[1., 0.],
[0., 1.],
[1., 0.],
[0., 1.],
[1., 0.],
[0., 1.]])
In [78]: np.hstack((__,_))
Out[78]:
array([[1., 0., 0., 1., 0.],
[1., 0., 0., 0., 1.],
[0., 1., 0., 1., 0.],
[0., 1., 0., 0., 1.],
[0., 0., 1., 1., 0.],
[0., 0., 1., 0., 1.]])
Combinatoric / cartesian product of Numpy arrays without iterators and/or loop(s)
Install Scikit-Learn http://scikit-learn.org/
from sklearn.utils.extmath import cartesian
print cartesian([a_p1, a_p2])
It should output
[[ 0 20]
[ 0 21]
[ 0 22]
[ 0 23]
[ 0 24]
[ 1 20]
[ 1 21]
[ 1 22]
[ 1 23]
[ 1 24]
[ 2 20]
[ 2 21]
[ 2 22]
[ 2 23]
[ 2 24]
[ 3 20]
[ 3 21]
[ 3 22]
[ 3 23]
[ 3 24]]
This solution was taken from a similar question:
Using numpy to build an array of all combinations of two arrays
Cartesian product of arbitrary-dimensional coordinates in arbitrary-dimensional arrays
A memory efficient way is broadcasted assignment:
def cartesian_product(x, y):
if x.ndim < 2:
x = np.atleast_2d(x).T
if y.ndim < 2:
y = np.atleast_2d(y).T
sx, sy = x.shape, y.shape
sz = sy[:-1] + sx[:-1] + (sy[-1] + sx[-1],)
z = np.empty(sz, np.result_type(x, y))
# Broadcasted assignment
z[...,:sx[-1]] = x
z[...,sx[-1]:] = y.reshape(sy[:-1] + (x.ndim-1)*(1,) + (sy[-1],))
return z
In case you need the details on broadcasting, this page has you covered.
Why is numpy cartesian product slower than pure python version?
Why current implementations are slow
While the first solution is faster than the second one, it is quite inefficient since it creates a lot of temporary CPython objects (at least 6 per item of itertools.product
). Creating a lot of objects is expensive because they are dynamically allocated and reference-counted by CPython. The Numpy function cartesian_product
is pretty fast but the iteration over the resulting array is very slow because it creates a lot of Numpy views and operates on numpy.uint8
instead of CPython int
. Numpy types and functions introduce a huge overhead for very small arrays.
Numpy can be used to speed up this operation as shown by @AlainT but this is not trivial to do and Numpy does not shine to solve such problems.
How to improve performance
One solution is to use Numba to do the job yourself more efficiently and let the Numba's JIT compiler optimizes loops. You can use 3 nested loops to efficiently generate the value of the Cartesian product and filter items. A dictionary can be used to track already seen values. The tuple of 3 items can be packed into one integer so to reduce the memory footprint and improve performance (so the dictionary can better fit in CPU caches and avoid the creation of slow tuple objects).
Here is the resulting code:
import numba as nb
# Signature of the function (parameter types)
# Note: `::1` means the array is contiguous
@nb.njit('(uint8[::1], uint8[::1], uint8[::1])')
def with_numba(a, b, c):
seen = dict()
for va in a:
for vb in b:
for vc in c:
# If the 3 values are different
if va != vb and vc != vb and vc != va:
# Sort the 3 values using a fast sorting network
v1, v2, v3 = va, vb, vc
if v1 > v2: v1, v2 = v2, v1
if v2 > v3: v2, v3 = v3, v2
if v1 > v2: v1, v2 = v2, v1
# Compact the 3 values into one 32-bit integer
packedKey = (np.uint32(v1) << 16) | (np.uint32(v2) << 8) | np.uint32(v3)
# Is the sorted tuple (v1,v2,v3) already seen?
if packedKey not in seen:
# Add the value and remember the ordered tuple (va,vb,vc)
packedValue = (np.uint32(va) << 16) | (np.uint32(vb) << 8) | np.uint32(vc)
seen[packedKey] = packedValue
res = np.empty((len(seen), 3), dtype=np.uint8)
for i, packed in enumerate(seen.values()):
res[i, 0] = np.uint8(packed >> 16)
res[i, 1] = np.uint8(packed >> 8)
res[i, 2] = np.uint8(packed)
return res
with_numba(a, b, c)
Benchmark
Here are results on my i5-9600KF processor:
numpy: 122.1 ms (x 1.0)
no_numpy: 49.6 ms (x 2.5)
AlainT's solution: 49.0 ms (x 2.5)
mathfux's solution 34.2 ms (x 3.5)
mathfux's optimized solution 7.5 ms (x16.2)
with_numba: 4.9 ms (x24.9)
The provided solution is about 25 times faster than the slowest implementation and about 1.5 time faster than the fastest provided implementation so far.
The current Numba code is bounded by the speed of the Numba dictionary operations. The code can be optimized using more low-level tricks. On solution is to replace the dictionary by a packed boolean array (1 item = 1 bit) of size 256**3/8
to track the values already seen (by checking the packedKey
th bit). The packed values can be directly added in res
if the fetched bit is not set. This requires res
to be preallocated to the maximum size or to implement an exponentially growing array (like list
in Python or std::vector
in C++). Another optimization is to sort the list and use a tiling strategy so to improve cache locality. Such optimization are far from being easy to implement but I expect them to drastically speed up the execution.
If you plan to use more arrays, then the hash-map can become a bottleneck and a bit-array can be quite big. While using tiling certainly help to reduce the memory footprint, you can speed up the implementation by a large margin using Bloom filters. This probabilist data structure can speed up the execution by skipping many duplicates without causing any cache misses and with a low memory footprint. You can remove most of the duplicates and then sort the array so to then remove the duplicates. Regarding your problem, a radix sort may be faster than usual sorting algorithms.
how can I construct a list of NumPy arrays from these two arrays?
Short and simple:
import numpy as np
x = np.array([111, 303, 405, 513])
y = np.array([523, 546 , 569 , 603 ])
def coupler(lst): return [[c1,c2] for c1,c2 in zip(lst,lst[1:])]
x_couples,y_couples=coupler(x),coupler(y)
boxes=[]
for x1,x2 in x_couples:
for y1,y2 in y_couples:
boxes.append(np.array([[x1,y1],[x2,y1],[x2,y2],[x1,y2]]))
The output:
[array([[111, 523],
[303, 523],
[303, 546],
[111, 546]]),
array([[111, 546],
[303, 546],
[303, 569],
[111, 569]]),...]
Numpy project array of m X coords and array of n Y coords into mxn list of (X,Y) coords
Pythonic comprehension is probably good enough for this one:
coords = [(x,y) for x in xs for y in ys]
Cartesian product but remove duplicates up to cyclic permutations
Here is a naive solution in python:
- Generate all combinations from the Cartesian product of
{1, 2, ...,n}
with itselfr
times; - Only keep one representative combination for each equivalency class; drop all other combinations that are equivalent to this representative combination.
This means we must have some way to compare combinations, and for instance, only keep the smallest combination of every equivalency class.
from itertools import product
def is_representative(comb):
return all(comb[i:] + comb[:i] >= comb
for i in range(1, len(comb)))
def cartesian_product_up_to_cyclic_permutations(n, r):
return filter(is_representative,
product(range(n), repeat=r))
print(list(cartesian_product_up_to_cyclic_permutations(3, 3)))
# [(0, 0, 0), (0, 0, 1), (0, 0, 2), (0, 1, 1), (0, 1, 2), (0, 2, 1), (0, 2, 2), (1, 1, 1), (1, 1, 2), (1, 2, 2), (2, 2, 2)]
print(list(cartesian_product_up_to_cyclic_permutations(2, 4)))
# [(0, 0, 0, 0), (0, 0, 0, 1), (0, 0, 1, 1), (0, 1, 0, 1), (0, 1, 1, 1), (1, 1, 1, 1)]
You mentioned that you wanted to implement the algorithm in C++. The product
function in the python code behaves just like a big for
-loop that generates all the combinations in the Cartesian product. See this related question to implement Cartesian product in C++: Is it possible to execute n number of nested "loops(any)" where n is given?.
Related Topics
How to Concatenate Str and Int Objects
How to Break Out of Multiple Loops
Get a List of Numbers as Input from the User
Getting the Name of a Variable as a String
How to Get a Cron Like Scheduler in Python
How to Watch a File For Changes
Difference Between '/' and '//' When Used For Division
Does Pandas Iterrows Have Performance Issues
How to Parse Data in Json Format
Why Are Multiple Instances of Tk Discouraged
How to Sort a List/Tuple of Lists/Tuples by the Element At a Given Index
Selenium Using Python - Geckodriver Executable Needs to Be in Path
How to Get Monotonic Time Durations in Python
How to Copy a Dictionary and Only Edit the Copy
"Least Astonishment" and the Mutable Default Argument
Is It Pythonic to Use List Comprehensions For Just Side Effects