[Paraview] MPI-Error during using "Programmable Filter"...
Stefan Melber
Stefan.Melber at DLR.de
Wed Jun 21 00:57:46 EDT 2017
Hi Utkarsh,
the state file is the one with errors. I got a long list of differences
between v_bad and v in the window, where i have started the pvservers
with mpi like
-0.149121856689 0.485999450684
-0.158555145264 0.47960144043
-0.17734588623 0.457775421143
-0.149189071655 0.476917114258
-0.241045379639 0.370332946777
-0.177021789551 0.414346313477
-0.333812713623 0.232836074829
-0.237988204956 0.299920196533
-0.446063957214 0.0598551177979
-0.324894866943 0.146606750488
-0.566599311829 -0.131107559204
-0.42886264801 -0.0301560974121
Running it in serial no differences between v_bad and v occur.
I tried the same binaries of 5.4 from the webpage ... same error.
Best regards,
Stefan
> Stefan,
>
> I am a little confused on how test your state file. Is the attached
> state bad or good? I don't get any errors even when I run in parallel
> with 2 pvservers.
>
> Utkarsh
>
> On Mon, Jun 19, 2017 at 8:34 AM, Stefan Melber <Stefan.Melber at dlr.de> wrote:
>> Hi,
>>
>>
>> Update: found the reason of the mpi-error - maybe a PV-bug:
>>
>> If i use this formula in a programmable filter
>>
>> v_bad = min((1., max((-1., tmp2))))
>>
>>
>> then paraview crashes in parallel. If i use the same but written out in
>>
>> if (tmp2 > -1.):
>>
>> if (tmp2 > 1.):
>>
>> v = 1.
>>
>> else:
>>
>> v = tmp2
>>
>> else:
>>
>> v = -1.
>>
>> the result for v is the same, however paraview does not crash.
>>
>> Further more, the results of v_bad and v in parallel differ, in serial mode
>> they are the same. Find attached a short restart file, which demonstrates
>> the effect.
>>
>> Stefan
>>
>> Hi,
>>
>>
>> i have used a programmable filter (find it below) to calculate the
>> lambda2-criterion from gradients of velocity. In serial mode this works fine
>> - running it in parallel an error message comes up. Any idea what does it
>> mean and whats wrong with the filter?
>>
>> Best regards,
>>
>> Stefan
>>
>>
>>
>>
>> Executing with: 0
>> Executing with: 0
>> Executing with: 0
>> Executing with: 0
>> Executing with: 0
>> Executing with: 0
>> Executing with: 0
>> Executing with: 0
>>
>> 2385958
>> 2406377
>> 2220610
>> 2292001
>> 2144431
>> 2331490
>> 2061524
>> 2344595
>> Traceback (most recent call last):
>> File "<string>", line 22, in <module>
>> File "<string>", line 76, in RequestData
>> File
>> "/opt/PARAVIEW_5_4_0_OpenGL2/ParaView-v5.4.0.bin/lib/site-packages/vtk/numpy_interface/algorithms.py",
>> line 358, in max
>> return _global_func(MaxImpl(), array, axis, controller)
>> File
>> "/opt/PARAVIEW_5_4_0_OpenGL2/ParaView-v5.4.0.bin/lib/site-packages/vtk/numpy_interface/algorithms.py",
>> line 199, in _global_func
>> max_dims, size = _reduce_dims(res, comm)
>> File
>> "/opt/PARAVIEW_5_4_0_OpenGL2/ParaView-v5.4.0.bin/lib/site-packages/vtk/numpy_interface/algorithms.py",
>> line 168, in _reduce_dims
>> comm.Allreduce([dims, mpitype], [max_dims, mpitype], MPI.MAX)
>> File "MPI/Comm.pyx", line 715, in mpi4py.MPI.Comm.Allreduce
>> (/opt/PARAVIEW_5_4_0_OpenGL2/ParaView-v5.4.0/VTK/ThirdParty/mpi4py/vtkmpi4py/src/mpi4py.MPI.c:99224)
>> mpi4py.MPI.Exception: MPI_ERR_OTHER: known error not in list
>> Traceback (most recent call last):
>> File "<string>", line 22, in <module>
>> Traceback (most recent call last):
>> File "<string>", line 76, in RequestData
>> File "<string>", line 22, in <module>
>> File "<string>", line 76, in RequestData
>> File
>> "/opt/PARAVIEW_5_4_0_OpenGL2/ParaView-v5.4.0.bin/lib/site-packages/vtk/numpy_interface/algorithms.py",
>> line 358, in max
>> File
>> "/opt/PARAVIEW_5_4_0_OpenGL2/ParaView-v5.4.0.bin/lib/site-packages/vtk/numpy_interface/algorithms.py",
>> line 358, in max
>> return _global_func(MaxImpl(), array, axis, controller)
>> File
>> "/opt/PARAVIEW_5_4_0_OpenGL2/ParaView-v5.4.0.bin/lib/site-packages/vtk/numpy_interface/algorithms.py",
>> line 199, in _global_func
>> max_dims, size = _reduce_dims(res, comm)
>> return _global_func(MaxImpl(), array, axis, controller)
>> File
>> "/opt/PARAVIEW_5_4_0_OpenGL2/ParaView-v5.4.0.bin/lib/site-packages/vtk/numpy_interface/algorithms.py",
>> line 199, in _global_func
>> File
>> "/opt/PARAVIEW_5_4_0_OpenGL2/ParaView-v5.4.0.bin/lib/site-packages/vtk/numpy_interface/algorithms.py",
>> line 168, in _reduce_dims
>> max_dims, size = _reduce_dims(res, comm)
>> comm.Allreduce([dims, mpitype], [max_dims, mpitype], MPI.MAX)
>> File "MPI/Comm.pyx", line 715, in mpi4py.MPI.Comm.Allreduce
>> (/opt/PARAVIEW_5_4_0_OpenGL2/ParaView-v5.4.0/VTK/ThirdParty/mpi4py/vtkmpi4py/src/mpi4py.MPI.c:99224)
>> File
>> "/opt/PARAVIEW_5_4_0_OpenGL2/ParaView-v5.4.0.bin/lib/site-packages/vtk/numpy_interface/algorithms.py",
>> line 168, in _reduce_dims
>> mpi4py.MPI.Exception: MPI_ERR_OTHER: known error not in list
>> comm.Allreduce([dims, mpitype], [max_dims, mpitype], MPI.MAX)
>> File "MPI/Comm.pyx", line 715, in mpi4py.MPI.Comm.Allreduce
>> (/opt/PARAVIEW_5_4_0_OpenGL2/ParaView-v5.4.0/VTK/ThirdParty/mpi4py/vtkmpi4py/src/mpi4py.MPI.c:99224)
>> mpi4py.MPI.Exception: MPI_ERR_OTHER: known error not in list
>>
>>
>>
>>
>>
>>
>>
>>
>> grad = inputs[0].PointData["Gradients"]
>>
>> npoints = len(grad)
>>
>> print npoints
>>
>> data = []
>>
>> for i in range(npoints):
>>
>> # Gradients
>>
>> dvx = grad[i,:][0]
>>
>> dvy = grad[i,:][1]
>>
>> dvz = grad[i,:][2]
>>
>> # Symmetrical part of flow tensor -> S
>>
>> s0 = dvx[0]
>>
>> s1 = 0.5 * (dvx[1] + dvy[0])
>>
>> s2 = 0.5 * (dvx[2] + dvz[0])
>>
>> s3 = 0.5 * (dvy[0] + dvx[1])
>>
>> s4 = dvy[1];
>>
>> s5 = 0.5 * (dvy[2] + dvz[1])
>>
>> s6 = 0.5 * (dvz[0] + dvx[2])
>>
>> s7 = 0.5 * (dvz[1] + dvy[2])
>>
>> s8 = dvz[2]
>>
>> # Antisymmetrical part of flow tensor -> Omega
>>
>> Omega0 = 0.0
>>
>> Omega1 = 0.5 * (dvx[1] - dvy[0])
>>
>> Omega2 = 0.5 * (dvx[2] - dvz[0])
>>
>> Omega3 = 0.5 * (dvy[0] - dvx[1])
>>
>> Omega4 = 0.0
>>
>> Omega5 = 0.5 * (dvy[2] - dvz[1])
>>
>> Omega6 = 0.5 * (dvz[0] - dvx[2])
>>
>> Omega7 = 0.5 * (dvz[1] - dvy[2])
>>
>> Omega8 = 0.0
>>
>> # Matrix M = (S^2 + Omega^2)
>>
>> # M is symmetric
>>
>> # | m0 m1 m2 |
>>
>> # M = | m1 m3 m4 |
>>
>> # | m2 m4 m5 |
>>
>> m0 = s0*s0 + s1*s3 + s2*s6 + Omega1*Omega3 + Omega2*Omega6
>>
>> m1 = s0*s1 + s1*s4 + s2*s7 + Omega2*Omega7
>>
>> m2 = s0*s2 + s1*s5 + s2*s8 + Omega1*Omega5
>>
>> m3 = s3*s1 + s4*s4 + s5*s7 + Omega3*Omega1 + Omega5*Omega7
>>
>> m4 = s3*s2 + s4*s5 + s5*s8 + Omega3*Omega2
>>
>> m5 = s6*s2 + s7*s5 + s8*s8 + Omega6*Omega2 + Omega7*Omega5
>>
>> # computing now the eigenvalues of M
>>
>> # | M - lambda I | = 0
>>
>> # returns the characteristic equation:
>>
>> # lambda^3 + p*lambda^2 + q*lambda + r = 0
>>
>> # due to symmetric Matrix the following assumption can be made:
>>
>> # all three eigenvalues will be real root values ( no imaginary parts )
>>
>> p = -1*( m0 + m3 + m5)
>>
>> q = -1*( m1*m1 + m2*m2 + m4*m4 - m0*m3 - m0*m5 - m3*m5)
>>
>> r = -1*( m0*m3*m5 + 2*m1*m2*m4 - m2*m2*m3 - m4*m4*m0 - m1*m1*m5)
>>
>> # computing now the reduced equation
>>
>> # lambda^3 + s*lambda + t = 0
>>
>> sx = (3.0 * q - p*p) / 3.0
>>
>> t = (2.0 * math.pow(p, 3.0) / 27.0) - (p * q / 3.0) + r
>>
>> # calculate unordered eigenvalues (Cardano's method)
>>
>> # -> at first check if ( 1.0 / sqrtt(-s^3) ) can be computed
>>
>> l = [0,0,0]
>>
>> if ((sx * sx * sx) < -1e-16):
>>
>> tmp1 = math.sqrt( -1* sx / 3.0 )
>>
>> tmp2 = -1* t / (2.0 * tmp1*tmp1*tmp1)
>>
>> c1 = 2.0 * tmp1
>>
>> c2 = math.acos( min((1., max((-1., tmp2))) ) ) / 3.0
>>
>> l[0] = c1 * math.cos(c2) - p / 3.0
>>
>> l[1] = c1 * math.cos(c2 + 2.*math.pi/3.) - p / 3.0
>>
>> l[2] = c1 * math.cos(c2 + 4.*math.pi/3.) - p / 3.0
>>
>> # sort eigenvalues
>>
>> # lambda1 <= lambda2 <= lambda3
>>
>> l.sort()
>>
>> # save second eigenvalue
>>
>> data.append(l[1])
>>
>> output.PointData.append(data, "lambda2")
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>> _______________________________________________
>> Powered by www.kitware.com
>>
>> Visit other Kitware open-source projects at
>> http://www.kitware.com/opensource/opensource.html
>>
>> Please keep messages on-topic and check the ParaView Wiki at:
>> http://paraview.org/Wiki/ParaView
>>
>> Search the list archives at: http://markmail.org/search/?q=ParaView
>>
>> Follow this link to subscribe/unsubscribe:
>> http://public.kitware.com/mailman/listinfo/paraview
>>
>>
>> _______________________________________________
>> Powered by www.kitware.com
>>
>> Visit other Kitware open-source projects at
>> http://www.kitware.com/opensource/opensource.html
>>
>> Please keep messages on-topic and check the ParaView Wiki at:
>> http://paraview.org/Wiki/ParaView
>>
>> Search the list archives at: http://markmail.org/search/?q=ParaView
>>
>> Follow this link to subscribe/unsubscribe:
>> http://public.kitware.com/mailman/listinfo/paraview
>>
> .
>
More information about the ParaView
mailing list