[Numpy-discussion] how to vectorize a row of eigenvalue computations ?

[Charles R Harris]

On Thu, Oct 16, 2008 at 1:22 PM, LB <berthe.loic@gmail.com> wrote:

> I've got an array S of shape (N, 6) with N >> 100000 containing the
> six components of a stress field given in N points.
>
> I need to make a lot of computation of principal stresses, which are
> in fact the eigenvalues of the stress tensors.
>
> I'm using the basic code described below :
>
> import numpy as np
>
> def calc_principal_stresses(S):
>    """ Return the principal stress corresponding to the tensor S.
>    S is interpreted as an array containing [Sx, Sy, Sz, Sxy, Syz,
> Szx]
>
>    Return array([S3, S2, S1])
>    """
>    p_stresses = np.linalg.eigvalsh( np.array(
>        [  [ S[0], S[3], S[5]],
>           [ S[3], S[1], S[4]],
>           [ S[5], S[4], S[2]],
>        ]))
>
>    return p_stresses.sort()
>
> p_stresses = array([ calc_principal_stresses(s) for s in S])
>

Something like this?

In [1]: m = np.zeros((2,6)) + np.arange(6)

In [2]: m = m[:,[0,3,5,3,1,4,5,4,2]].reshape(-1,3,3)

In [3]: p_stresses = np.sort([np.linalg.eigvalsh(i) for i in m], axis = -1)

In [4]: p_stresses
Out[4]:
array([[-4.11645769, -2.04234215,  9.15879984],
       [-4.11645769, -2.04234215,  9.15879984]])


You might want to break this into more statements for clarity.

Chuck
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://projects.scipy.org/pipermail/numpy-discussion/attachments/20081016/22879d26/attachment-0001.html