# Eigenvalue computation
**URL:** https://discourse.paraview.org/t/eigenvalue-computation/9307
**Category:** ParaView Support
**Tags:** python
**Created:** [March 29, 2022, 1:10pm UTC](https://discourse.paraview.org/t/eigenvalue-computation/9307 "2022-03-29T13:10:05Z")
**Posts on this page:** 4
**Page:** 1
### Author:  [@ArroyoAndrea](https://discourse.paraview.org/u/ArroyoAndrea)
#### Post date: [March 29, 2022, 1:10pm UTC](https://discourse.paraview.org/t/eigenvalue-computation/9307/1 "2022-03-29T13:10:06Z")
Hello all,
I am implementing Lambda-2 criterion in Paraview 5.6.0 mounted on Linux. What it does is to provide a view of the presence of a vortex in a 3D flow by knowing when at least two of the eigenvalues of (||S||^2 + ||Omega||^2) are negative. Here, S and Omega are the symmetric and antisymmetric parts of the Velocity gradient tensor.
By following the instructions of some other posts where they dealt with this implementation, I have made the following:
- Gradient of unstructured data set. Inputting Velocity as Scalar Array and Compute Gradients
- Python calculator. Expression “eigenvalue(strain(Velocity)\*\*2 + (Gradients - strain(Velocity))\*\*2)”. Array name “Lambda”. This computes the three eigenvalues of the expression ||S||^2 + ||Omega||^2
- Calculator. “Lambda\_Y”. Result Array Name “Lambda2”. Get the second eigenvalue of the previous result. Since the three eigenvalues are ordered such that Lambda1 \>= Lambda2 \>= Lambda3, if Lambda2 is negative, Lambda3 is negative too, so a vortex is found.
1- I am having some issues since the results that I obtain are not the expected ones. So I took one step back and I have been looking terms by terms what is happening there. I have realized that if I compute the Gradient of Unstructured Data Set to get the velocity gradient and I compute my means of a Python calculator “strain(Velocity)”, the tensor I am obtaining in both cases does not share the trace. Both of them should share the trace since strain(Velocity) = 1/2\*(Gradient(Velocity) + transpose(Gradient(Velocity)), and it is not the case.
2- Also, I am evaluating if the issues comes from the eigenvalue calculation itself. To do so, I tried to simple compute the eigenvalues of the velocity gradient by means of the function eigenvalue(Gradients) as well as by means of numpy.linalg.eigvals(Gradients). Think this second version to compute should work better, however, it draws the following error (it does not recognize the Gradients tensor as a tensor):
Traceback (most recent call last):
File “”, line 4, in
File “/home/softs/softs/ParaView-5.6.0-MPI-Linux-64bit/lib/python2.7/site-packages/paraview/calculator.py”, line 192, in execute
retVal = compute(inputs, expression, ns=variables)
File “/home/softs/softs/ParaView-5.6.0-MPI-Linux-64bit/lib/python2.7/site-packages/paraview/calculator.py”, line 140, in compute
retVal = eval(expression, globals(), mylocals)
File “”, line 1, in
File “/home/softs/softs/ParaView-5.6.0-MPI-Linux-64bit/lib/python2.7/site-packages/numpy/linalg/linalg.py”, line 1021, in eigvals
\_assertRankAtLeast2(a)
File “/home/softs/softs/ParaView-5.6.0-MPI-Linux-64bit/lib/python2.7/site-packages/numpy/linalg/linalg.py”, line 204, in \_assertRankAtLeast2
‘at least two-dimensional’ % a.ndim)
numpy.linalg.linalg.LinAlgError: 0-dimensional array given. Array must be at least two-dimensional
Thank you in advance,
Andrea Arroyo
---
### Author:  [@MiguelAngelAGUIRRE](https://discourse.paraview.org/u/MiguelAngelAGUIRRE)
#### Post date: [April 19, 2022, 9:07am UTC](https://discourse.paraview.org/t/eigenvalue-computation/9307/2 "2022-04-19T09:07:52Z")
Hello,
Here you have a code snippet to compute Lambda 2 in a Python Programmable Filter:
import numpy as np
from vtk.numpy\_interface import algorithms as algs
INPUT=inputs[0]
Vx = INPUT.PointData[“X\_VELOCITY”]
Vy = INPUT.PointData[“Y\_VELOCITY”]
Vz = INPUT.PointData[“Z\_VELOCITY”]
Velocity\_Vector = algs.make\_vector(Vx,Vy,Vz)
STRAIN = strain(Velocity\_Vector)
V\_gradient = gradient(Velocity\_Vector)
AAA=STRAIN\*\*2+(V\_gradient - STRAIN)\*\*2
Lambda\_2 = np.linalg.eigvals(AAA[:])
Lambda\_2 = np.real(Lambda\_2)
output.PointData.append(Lambda\_2, “Lambda\_2”)
Regards
Miguel
---
### Author:  [@Gening\_Dong](https://discourse.paraview.org/u/Gening_Dong)
#### Post date: [February 24, 2023, 10:13pm UTC](https://discourse.paraview.org/t/eigenvalue-computation/9307/3 "2023-02-24T22:13:57Z")
Hello,
I’d like to calculate the principal strains so I modified this code into
```auto
import numpy as np
from vtk.numpy_interface import algorithms as algs
INPUT=inputs[0]
Dx = INPUT.PointData["MS_Displacement_X"]
Dy = INPUT.PointData["MS_Displacement_Y"]
Dz = INPUT.PointData["MS_Displacement_Z"]
Displacement_Vector = algs.make_vector(Dx,Dy,Dz)
STRAIN = strain(Displacement_Vector)
principal_strain = np.linalg.eigvals(STRAIN[:])
outarray0.Arrays[ii] = principal_strain[:,0]
outarray1.Arrays[ii] = principal_strain[:,1]
outarray2.Arrays[ii] = principal_strain[:,2]
output.PointData.append(outarray0, "principal_strain_0")
output.PointData.append(outarray0, "principal_strain_1")
output.PointData.append(outarray0, "principal_strain_2")
```
However, there’s an error shows:
```auto
AttributeError: 'VTKNoneArray' object has no attribute 'shape'
Traceback (most recent call last):
File "", line 22, in
File "", line 10, in RequestData
File "C:\Program Files\ParaView 5.10.1-MPI-Windows-Python3.9-msvc2017-AMD64\bin\Lib\site-packages\vtkmodules\numpy_interface\algorithms.py", line 861, in make_vector
return algs.make_vector(arrayx, arrayy, arrayz)
File "C:\Program Files\ParaView 5.10.1-MPI-Windows-Python3.9-msvc2017-AMD64\bin\Lib\site-packages\vtkmodules\numpy_interface\internal_algorithms.py", line 558, in make_vector
if len(ax.shape) != 1 or len(ay.shape) != 1 or (az is not None and len(az.shape) != 1):
AttributeError: 'VTKNoneArray' object has no attribute 'shape'
```
Any ideas? Thanks in advance!
Best,
Gening
---
### Author:  [@bastian](https://discourse.paraview.org/u/bastian)
#### Post date: [August 18, 2023, 4:58pm UTC](https://discourse.paraview.org/t/eigenvalue-computation/9307/4 "2023-08-18T16:58:25Z")
> [@MiguelAngelAGUIRRE](#):
>
> Hello,
>
> Here you have a code snippet to compute Lambda 2 in a Python Programmable Filter:
>
> ```python
> import numpy as np
> from vtk.numpy_interface import algorithms as algs
>
> INPUT=inputs[0]
> Vx = INPUT.PointData[“X_VELOCITY”]
> Vy = INPUT.PointData[“Y_VELOCITY”]
> Vz = INPUT.PointData[“Z_VELOCITY”]
>
> Velocity_Vector = algs.make_vector(Vx,Vy,Vz)
>
> STRAIN = strain(Velocity_Vector)
> V_gradient = gradient(Velocity_Vector)
> AAA=STRAIN**2+(V_gradient - STRAIN)**2
> Lambda_2 = np.linalg.eigvals(AAA[:])
> Lambda_2 = np.real(Lambda_2)
> output.PointData.append(Lambda_2, “Lambda_2”)
>
> ```
I think that you have a few errors in your code:
1. `STRAIN**2` is probably an elementwise square. What you actually want is the matrix product that matches the signature _(n,k),(k,m)→(n,m)_, which matmul should provide
2. `np.linalg.eigvals` does not order the eigenvalues, but you want lambda\_2 to correspond to the second largest eigenvalue.
This filter should give the proper result, but I haven’t validated it:
```python
import numpy as np
from vtk.numpy_interface import algorithms as algs
# Compute the lambda_2 Vortex criterion, see pp.76-77 in
# J. Jeong and F. Hussain, “On the identification of a vortex,”
# Journal of Fluid Mechanics, vol. 285, pp. 69–94, 1995, doi: 10.1017/S0022112095000462.
vvector = inputs[0].PointData['velocity']
vstrain = strain(vvector)
vskew = gradient(vvector) - vstrain
aaa = matmul(vstrain, vstrain) + matmul(vskew, vskew)
# since aaa is symmetric, it only has real eigenvalues
lambdas = eigenvalue(aaa)
lambdas = real(lambdas)
lambda2 = sort(lambdas)[:,1]
output.PointData.append(lambda2, 'lambda2')
```