Difference between revisions of "Numpy linalg"

From relax wiki
Jump to navigation Jump to search
Line 107: Line 107:
 
http://wiki.scipy.org/Cookbook/SegmentAxis
 
http://wiki.scipy.org/Cookbook/SegmentAxis
  
 +
Simple, with window 3
 +
<source lang="python">
 +
row, col = 3, 5
 +
 +
x=np.arange(row*col).reshape((row, col))
 +
print x
 +
 +
window = 3
 +
# Shape should end out in tuple.
 +
shape = tuple([row, col - window + 1,  window])
 +
print "shape is:", shape
 +
 +
# Get strides
 +
x_str = x.strides
 +
print "strides is", x_str
 +
 +
# Get itemsize, Length of one array element in bytes. Depends on dtype. float64=8, complex128=16.
 +
x_itz = x.itemsize
 +
print "itemsize is", x_str
 +
 +
# Again, strides should be in tuple
 +
cut_strides = tuple(list(x_str) + [x_itz])
 +
print "cut strides are", cut_strides
 +
 +
y = as_strided(x, shape=shape, strides=cut_strides)
 +
</source>
 +
 +
And then
  
 
<source lang="python">
 
<source lang="python">

Revision as of 09:48, 21 June 2014

How to transpose higher dimension arrays

http://jameshensman.wordpress.com/2010/06/14/multiple-matrix-multiplication-in-numpy/

Faster dot product using BLAS

http://www.huyng.com/posts/faster-numpy-dot-product/

http://stackoverflow.com/questions/5990577/speeding-up-numpy-dot

http://wiki.scipy.org/PerformanceTips

http://thread.gmane.org/gmane.comp.python.numeric.general/28135/

Multi dot

http://wiki.scipy.org/Cookbook/MultiDot

Einsum

http://chintaksheth.wordpress.com/2013/07/31/numpy-the-tricks-of-the-trade-part-ii/

http://stackoverflow.com/questions/14758283/is-there-a-numpy-scipy-dot-product-calculating-only-the-diagonal-entries-of-the

a = np.arange(4).reshape(2,2)
print a
print "np.einsum('ii', a), row i multiplied downwards"
print np.einsum('ii', a)

print "np.einsum('ij', a), same matrix ?"
print np.einsum('ij', a)

print "np.einsum('ji', a), transpose"
print np.einsum('ji', a)

print "np.einsum('ij,jk', a, a), dot product"
print np.einsum('ij,jk', a, a)
print np.dot(a, a)

Ellipsis broadcasting in numpy.einsum

http://stackoverflow.com/questions/16591696/ellipsis-broadcasting-in-numpy-einsum

http://comments.gmane.org/gmane.comp.python.numeric.general/53705

http://stackoverflow.com/questions/118370/how-do-you-use-the-ellipsis-slicing-syntax-in-python

http://stackoverflow.com/questions/772124/what-does-the-python-ellipsis-object-do

"..." Is designed to mean at this point, insert as many full slices (:) to extend the multi-dimensional slice to all dimensions.

print ""
a = np.arange(4).reshape(2,2)
print "a is"
print a
print "dot a"
print np.dot(a, a)
# Expand one axis in start, and tile up 2 times.
a2 = np.tile(a[None,:], (2, 1, 1))
print "a2 shape", a2.shape
print "einsum dot product over higher dimensions"
a2_e = np.einsum('...ij,...jk', a2, a2)
print a2_e

# Expand one axis in start, and tile up 2 times.
a3 = np.tile(a2[None,:], (2, 1, 1, 1))
print "a3 shape", a3.shape
print "einsum dot product over higher dimensions"
a3_e = np.einsum('...ij,...jk', a3, a3)
print a3_e

# With Ellipsis and axis notation
a = np.arange(4).reshape(2,2)
a = np.arange(4).reshape(2,2)
print "a is"
print a
print "dot a"
print np.dot(a, a)
# Expand one axis in start, and tile up 2 times.
a2 = np.tile(a[None,:], (2, 1, 1))
print "a2 shape", a2.shape
print "einsum dot product over higher dimensions"
a2_e = np.einsum(a2, [Ellipsis, 0, 1], a2, [Ellipsis, 1, 2])
print a2_e
 
# Expand one axis in start, and tile up 2 times.
a3 = np.tile(a2[None,:], (2, 1, 1, 1))
print "a3 shape", a3.shape
print "einsum dot product over higher dimensions"
a3_e =  np.einsum(a3, [Ellipsis, 0, 1], a3, [Ellipsis, 1, 2])
print a3_e


Stride tricks

http://chintaksheth.wordpress.com/2013/07/31/numpy-the-tricks-of-the-trade-part-ii/

http://stackoverflow.com/questions/8070349/using-numpy-stride-tricks-to-get-non-overlapping-array-blocks

http://stackoverflow.com/questions/4936620/using-strides-for-an-efficient-moving-average-filter

http://www.rigtorp.se/2011/01/01/rolling-statistics-numpy.html

http://stackoverflow.com/questions/8070349/using-numpy-stride-tricks-to-get-non-overlapping-array-blocks

http://wiki.scipy.org/Cookbook/GameOfLifeStrides

http://wiki.scipy.org/Cookbook/SegmentAxis

Simple, with window 3

row, col = 3, 5

x=np.arange(row*col).reshape((row, col))
print x

window = 3
# Shape should end out in tuple.
shape = tuple([row, col - window + 1,  window])
print "shape is:", shape

# Get strides
x_str = x.strides
print "strides is", x_str

# Get itemsize, Length of one array element in bytes. Depends on dtype. float64=8, complex128=16.
x_itz = x.itemsize
print "itemsize is", x_str

# Again, strides should be in tuple
cut_strides = tuple(list(x_str) + [x_itz])
print "cut strides are", cut_strides

y = as_strided(x, shape=shape, strides=cut_strides)

And then

from numpy.lib.stride_tricks import as_strided
import numpy as np

NE, NS, NM, NO, ND, Row, Col = 1, 2, 2, 1, 2, 2, 2

mat = np.arange(1,NE*NS*NM*NO*ND*Row*Col+1).reshape(NE, NS, NM, NO, ND, Row, Col)
print "mat is:"
print mat

sz = mat.itemsize
print "itemsize is:"
print sz

print "strides is"
print mat.strides

print "height and width are"
print h, w

bh,bw = Row,Col

shape = (h/bh, w/bw, bh, bw)
print shape