[Numpy-discussion] vector to tensor matrix speed up

[Tim Hochberg]

Ferenc.Pintye at eu.decoma.com wrote:
>
>
> Hi users,
>
>       i have some problem in Numpy with indexing speed for array to tensor
> matrix transport.
> With 200000 cycles it's 9sec ! (without sort(eigvals() - functions)
>
> Many thanks
> f.
>   

The big problem you have here is that you are operating on your matrices 
one element at a time. Overhead is going to kill you doing this. The key 
thing is to operate on your data in chunks -- I'll show one way below 
that speeds thing up by a factor of 40 or so *neglecting 
sort(eigvals())*. Unfortunately, sort(eigvals()) is still likely to kill 
you since there's not really a good way to chunk that up. (I'll have 
more to say about that in a more general
sense in another email). Here's the two versions I compared, perhaps 
this will help you out:


    from numpy import *
    import timeit
    out = zeros((200000,11),float)
   
    def f0(out):
        for j in arange(0,200000):
            pass
   
    def f1(out):
        # stress values/matrix in array form for 200000 points
        #
        #.....out = filling matrix ...etc...
        #
        # stress tensor matrix
        eig = zeros((3,3),float)
        #
        #output for eigvalues
        eigwert = array([0,0,0])
        #
        for j in arange(0,200000):
              eig[0,0] = out[j,1]
              eig[1,1] = out[j,2]
              eig[2,2] = out[j,3]
              #
              eig[0,1] = out[j,4]
              eig[0,2] = out[j,6]
              eig[1,0] = out[j,4]
              eig[1,2] = out[j,5]
              eig[2,0] = out[j,6]
              eig[2,1] = out[j,5]
              #
              #~ eigwert = sort(eigvals(eig))
              out[j,7] = eigwert[2]
              out[j,8] = eigwert[1]
              out[j,9] = eigwert[0]
              out[j,10] = abs(eigwert[0]-eigwert[2])
        #
    def f2(out):
        # stress values/matrix in array form for 200000 points
       
        #
        #.....out = filling matrix ...etc...
        #
        # stress tensor matrix
   
        eig = zeros((100,3,3),float)
        #
        #output for eigvalues
        eigwert = zeros([100,3], dtype=float)
        #
        local_abs = abs
        rangen = range(0,200000,100)
        for j in rangen:
              eig[:,0,0] = out[j:j+100,1]
              eig[:,1,1] = out[j:j+100,2]
              eig[:,2,2] = out[j:j+100,3]
              #
              eig[:,1,0] = eig[:,0,1] = out[j:j+100,4]
              eig[:,2,0] = eig[:,0,2] = out[j:j+100,6]
              eig[:,2,1] = eig[:,1,2] = out[j:j+100,5]
              #
              #~ eigwert = sort(eigvals(eig))
        for j in rangen:
              out[j:j+100,7:10] = eigwert # Changed order of out here
              out[j:j+100,10] = abs(eigwert[:,0]-eigwert[:,2])
       
    if __name__ == '__main__':
        print timeit.Timer("f0(out)", "from scratch import f0, 
out").timeit(1)
        print timeit.Timer("f1(out)", "from scratch import f1, 
out").timeit(1)
        print timeit.Timer("f2(out)", "from scratch import f2, 
out").timeit(1)
>
> # stress values/matrix in array form for 200000 points
> out = zeros((200000,11),Float32)
> #
> #.....out = filling matrix ...etc...
> #
> # stress tensor matrix
> eig = zeros((3,3),Float32)
> #
> #output for eigvalues
> eigwert = array([0,0,0])
> #
> for j in arange(0,200000):
>       eig[0,0] = out[j,1]
>       eig[1,1] = out[j,2]
>       eig[2,2] = out[j,3]
>       #
>       eig[0,1] = out[j,4]
>       eig[0,2] = out[j,6]
>       eig[1,0] = out[j,4]
>       eig[1,2] = out[j,5]
>       eig[2,0] = out[j,6]
>       eig[2,1] = out[j,5]
>       #
>       eigwert = sort(eigvals(eig))
>       out[j,7] = eigwert[2]
>       out[j,8] = eigwert[1]
>       out[j,9] = eigwert[0]
>       out[j,10] = abs(eigwert[0]-eigwert[2])
> #
>
>
> -------------------------------------------------------------------------
> Take Surveys. Earn Cash. Influence the Future of IT
> Join SourceForge.net's Techsay panel and you'll get the chance to share your
> opinions on IT & business topics through brief surveys -- and earn cash
> http://www.techsay.com/default.php?page=join.php&p=sourceforge&CID=DEVDEV
> _______________________________________________
> Numpy-discussion mailing list
> Numpy-discussion at lists.sourceforge.net
> https://lists.sourceforge.net/lists/listinfo/numpy-discussion
>
>
>