

In [ ]:

    
%matplotlib inline



In [ ]:

    
import numpy as np
import scipy.sparse as sps
import matplotlib.pyplot as plt



In [ ]:

    
mnorm = lambda x: np.max(np.abs(x))
def kron3(A, B, C):
    return sps.kron(A, sps.kron(B, C))

2D einsum tests



In [ ]:

    
n = 10
A = np.random.rand(n, n)
B = np.random.rand(n, n)
b = np.random.rand(n)
c = np.random.rand(n)



In [ ]:

    
np.einsum('i,i', b, c)-np.dot(b, c)



In [ ]:

    
mnorm(np.einsum('ij,j', A, b)-A.dot(b))



In [ ]:

    
mnorm(np.einsum('ij,jk', A, B)-A.dot(B))



In [ ]:

    
mnorm(np.einsum('jk,ij', A, B)-B.dot(A))



In [ ]:

    
mnorm(np.einsum('kj,ij', A, B)-B.dot(A.T))

3D einsum tests



In [ ]:

    
n = 64
I = sps.eye(n)
D = np.random.rand(n, n)
D1 = kron3(I, I, D)
D2 = kron3(I, D, I)
D3 = kron3(D, I, I)



In [ ]:

    
A  = np.random.rand(n,n,n)
Ar = A.ravel()

$I\otimes I\otimes D$



In [ ]:

    
mnorm(np.einsum('kj,...ij', D, A).ravel()-D1.dot(Ar))



In [ ]:

    
mnorm(np.einsum('kj,aij', D, A).ravel()-D1.dot(Ar))



In [ ]:

    
mnorm(A.reshape((-1,n,n)).dot(D.T).ravel()-D1.dot(Ar))



In [ ]:

    
%timeit np.einsum('kj,...ij', D, A)



In [ ]:

    
%timeit np.einsum('kj,aij', D, A)



In [ ]:

    
%timeit D1.dot(Ar)



In [ ]:

    
%timeit A.reshape((-1,n,n)).dot(D.T)

$I\otimes D\otimes I$



In [ ]:

    
mnorm(np.einsum('ij,...jk', D, A).ravel()-D2.dot(Ar))



In [ ]:

    
mnorm(np.einsum('ij,ajk', D, A).ravel()-D2.dot(Ar))



In [ ]:

    
%timeit np.einsum('ij,...jk', D, A)



In [ ]:

    
%timeit np.einsum('ij,ajk', D, A)



In [ ]:

    
%timeit D2.dot(Ar)

$D\otimes I\otimes I$



In [ ]:

    
mnorm(np.einsum('ij,jkl', D, A).ravel()-D3.dot(Ar))



In [ ]:

    
%timeit np.einsum('ij,jkl', D, A)



In [ ]:

    
%timeit D3.dot(Ar)



In [ ]: