[Numpy-discussion] how to vectorize a row of eigenvalue computations ?
Charles R Harris
charlesr.harris@gmail....
Thu Oct 16 21:54:44 CDT 2008
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
More information about the Numpy-discussion
mailing list