[vtkusers] ICP Error: vtkMath::Jacobi: Error extracting eigenfunctions

[Dženan Zukić]

I think it is useful to put this finding on the mailing list for future
reference.

On Mon, Feb 6, 2012 at 18:28, Maarten Beek <beekmaarten at yahoo.com> wrote:

> Hate to admit, but I found a bug in my code as well. Same bug as Dzenan.
> Both vtkPoints objects were created in similar for-loops. So one of the
> arrays was not totally filled. Gives very unpredictable behavior (sometimes
> not even crashing).
> Switching the point sets however always caused crashes and should have
> been the clue I shouldn't have missed.
>
> Remain the facts that vtkMath::vtkJacobiN doesn't really handle this error
> nicely and that the algorithm used is potentially unstable.
>
> Maarten
>
>   ------------------------------
> *From:* Dženan Zukić <dzenanz at gmail.com>
> *To:* David Gobbi <david.gobbi at gmail.com>
> *Cc:* Maarten Beek <beekmaarten at yahoo.com>
> *Sent:* Friday, February 3, 2012 10:55:58 AM
>
> *Subject:* Re: [vtkusers] ICP Error: vtkMath::Jacobi: Error extracting
> eigenfunctions
>
> In trying to create a self-contained test program, I discovered a bug in
> my own code. I was supplying some junk points to the vtkLandmarkTransform
> (one points array not completely filled with values).
>
> Now the execution in release mode doesn't give any warnings either. Maybe
> something similar happened to you Maarten?
>
> 2012/2/1 David Gobbi <david.gobbi at gmail.com>
>
> Would it be possible for you to send me a short, simple program that
> includes vtkLandmarkTransform that triggers the error?  Or can you try
> to modify vtkMath::JacobiN in such a way that the error goes away,
> using the Numerical Recipes discussion on the following link as a
> guideline? http://www.nr.com/forum/showthread.php?p=4913
>
> There isn't much that I can do on my side unless I have some way of
> reproducing the error...
>
> 2012/2/1 Dženan Zukić <dzenanz at gmail.com>:
> > Same: problem is not encountered.
> >
> > 2012/2/1 David Gobbi <david.gobbi at gmail.com>
> >>
> >> Can you modify my test program so that main() uses vtkMath::JacobiN
> >> instead of using its own JacobiN, and then try it again?
> >>
> >> 2012/2/1 Dženan Zukić <dzenanz at gmail.com>:
> >> > This program does not crash in either debug or release mode. The
> >> > optimization options are the same (as CMake sets them).
> >> >
> >> > I have modified it to read input for the files I sent you, but it goes
> >> > through them without errors. I guess that whatever is happening, is
> not
> >> > triggered in this simplified example.
> >> >
> >> > 2012/1/31 David Gobbi <david.gobbi at gmail.com>
> >> >>
> >> >> I haven't been able to reproduce the problem, but at work all I have
> is
> >> >> a Mac and the optimizations it does are probably different from MSVC.
> >> >>
> >> >> Can you try different matrices with the attached program and see if
> >> >> you can get it to fail?
> >> >>
> >> >>  - David
> >> >>
> >> >>
> >> >> 2012/1/31 Dženan Zukić <dzenanz at gmail.com>:
> >> >> > Here is an output from non-release version, same input as before.
> >> >> >
> >> >> > 2012/1/31 Dženan Zukić <dzenanz at gmail.com>
> >> >> >>
> >> >> >> Here are 3 different runs of my program, identical input.
> >> >> >>
> >> >> >> 2012/1/31 David Gobbi <david.gobbi at gmail.com>
> >> >> >>>
> >> >> >>> Thanks for the info.  Can you post your compiler optimization
> flags
> >> >> >>> to
> >> >> >>> the list?  There is definitely a problem with JacobiN that needs
> to
> >> >> >>> be
> >> >> >>> fixed.  The fix should be very easy to do if we can get a matrix
> >> >> >>> that
> >> >> >>> reliably causes JacobiN to fail when certain optimization flags
> are
> >> >> >>> used, so if either of you can print out such a matrix (with all
> >> >> >>> numbers printed with a precision of 20 digits), it would be a big
> >> >> >>> help.
> >> >> >>>
> >> >> >>>  - David
> >> >> >>>
> >> >> >>>
> >> >> >>> 2012/1/31 Dženan Zukić <dzenanz at gmail.com>:
> >> >> >>> > I also saw this warning recently, from the block that "never
> gets
> >> >> >>> > called".
> >> >> >>> > And it only happened in release version of my program.
> >> >> >>> >
> >> >> >>> > On Tue, Jan 31, 2012 at 00:48, Maarten Beek
> >> >> >>> > <beekmaarten at yahoo.com>
> >> >> >>> > wrote:
> >> >> >>> >>
> >> >> >>> >> Hi David,
> >> >> >>> >>
> >> >> >>> >> Ai, seems these crashes are hard to get rid off....
> >> >> >>> >>
> >> >> >>> >> I currently use VTK within Amira modules.
> >> >> >>> >>
> >> >> >>> >> I changed the type of vtkPoints from the default VTK_FLOAT to
> >> >> >>> >> VTK_DOUBLE
> >> >> >>> >> when I convert an Amira object to a VTK object. This reduces
> the
> >> >> >>> >> occurrences
> >> >> >>> >> of the errors, but doesn't solve the issue.
> >> >> >>> >> (Ps: why is the default in vtkPoints still float, while
> >> >> >>> >> everywhere
> >> >> >>> >> else
> >> >> >>> >> floats were changed into doubles a long while ago?)
> >> >> >>> >>
> >> >> >>> >> I don't really know much about different optimizations, so I
> >> >> >>> >> just
> >> >> >>> >> copy-paste what I found in MSVS:
> >> >> >>> >>
> >> >> >>> >> I build my VTK with:
> >> >> >>> >> /O2 /Ob2 /FD /EHsc /MD /Fo"vtkCommon.dir\Release\\"
> >> >> >>> >> /Fd"G:\vtk-5.6.0\build_vs90_noQt\bin\Release/vtkCommon.pdb"
> /W4
> >> >> >>> >> /nologo /c
> >> >> >>> >> /TP /errorReport:prompt/Zm 1000 and a couple of '/I'-s and
> >> >> >>> >> '/D'-s
> >> >> >>> >>
> >> >> >>> >> I build my Amira modules with:
> >> >> >>> >> /O2 /Ob1 /Oi /Oy /GF /FD /EHsc /MD /GS- /Gy /fp:fast /openmp
> >> >> >>> >> /Fo"./../../obj/arch-Win64VC9-Optimize/GraphPackage/"
> >> >> >>> >> /Fd"./../.././bin/arch-Win64VC9-Optimize/GraphPackage.pdb" /W3
> >> >> >>> >> /nologo
> >> >> >>> >> /c
> >> >> >>> >> /Zi /wd4068 /wd4305 /wd4018 /wd4244 /wd4800 /wd4275 /wd4251
> >> >> >>> >> /errorReport:prompt and a couple of '/I'-s and '/D'-s
> >> >> >>> >>
> >> >> >>> >> I am working on getting vtkMath::JacobiN to print the matrix
> in
> >> >> >>> >> the
> >> >> >>> >> if-loop that 'gets never called' using vtkGenericWarningMacro
> to
> >> >> >>> >> avoid
> >> >> >>> >> having to include Amira stuff. Having a hard time getting that
> >> >> >>> >> to
> >> >> >>> >> work. (
> >> >> >>> >> char* test = "hello";vtkGenericWarningMacro( test ); won't do
> >> >> >>> >> the
> >> >> >>> >> trick).
> >> >> >>> >> Any suggestions?
> >> >> >>> >>
> >> >> >>> >> Maarten
> >> >> >>> >>
> >> >> >>> >> ________________________________
> >> >> >>> >> From: David Gobbi <david.gobbi at gmail.com>
> >> >> >>> >> To: Maarten Beek <beekmaarten at yahoo.com>
> >> >> >>> >> Cc: VTK list <vtkusers at vtk.org>
> >> >> >>> >> Sent: Monday, January 30, 2012 2:31:38 PM
> >> >> >>> >> Subject: Re: [vtkusers] ICP Error: vtkMath::Jacobi: Error
> >> >> >>> >> extracting
> >> >> >>> >> eigenfunctions
> >> >> >>> >>
> >> >> >>> >> Hi Maarten,
> >> >> >>> >>
> >> >> >>> >> I think there are two issues here, one is that
> >> >> >>> >> vtkLandmarkTransform
> >> >> >>> >> should have a check as you have described, but another
> >> >> >>> >> potentially
> >> >> >>> >> serious issue is that, as far as I understand it, JacobiN
> should
> >> >> >>> >> never
> >> >> >>> >> fail.  I did a google search and found a message from 2002 by
> >> >> >>> >> the
> >> >> >>> >> Numerical Recipes authors stating that the NR Jacobi code can
> >> >> >>> >> become
> >> >> >>> >> numerically unstable due to certain compiler optimizations.
> I'm
> >> >> >>> >> guessing that certain compiler optimizations could similarly
> >> >> >>> >> make
> >> >> >>> >> the
> >> >> >>> >> VTK Jacobi code become unstable.
> >> >> >>> >> http://www.nr.com/forum/showthread.php?p=4913
> >> >> >>> >>
> >> >> >>> >> Can you provide info about your compiler/flags etc?  Also, is
> it
> >> >> >>> >> possible to get your code to print a matrix that causes
> JacobiN
> >> >> >>> >> to
> >> >> >>> >> fail?
> >> >> >>> >>
> >> >> >>> >> - David
> >> >> >>> >>
> >> >> >>> >>
> >> >> >>> >> On Mon, Jan 30, 2012 at 9:44 AM, Maarten Beek
> >> >> >>> >> <beekmaarten at yahoo.com>
> >> >> >>> >> wrote:
> >> >> >>> >> > Hi all,
> >> >> >>> >> >
> >> >> >>> >> > I get the following warning when using ICP:
> >> >> >>> >> >
> >> >> >>> >> > vtkMath::Jacobi: Error extracting eigenfunctions
> >> >> >>> >> >
> >> >> >>> >> > It comes from vtkMath::vtkJacobiN that calculates the
> >> >> >>> >> > eigenvalues,
> >> >> >>> >> > -vectors
> >> >> >>> >> > of a matrix.
> >> >> >>> >> > The function returns 0, when this error occurs (otherwise
> 1),
> >> >> >>> >> > but
> >> >> >>> >> > vtkLandmarkTransform doesn't check for the return value and
> >> >> >>> >> > thus
> >> >> >>> >> > uses
> >> >> >>> >> > undefined values for the eigenvalues, -vectors.
> >> >> >>> >> > This means my point set becomes invalid.
> >> >> >>> >> > I believe the way to avoid this error is to 'giggle' the
> >> >> >>> >> > coordinates
> >> >> >>> >> > a
> >> >> >>> >> > little when it happens, but how can I intervene when
> >> >> >>> >> > vtkLandmarkTransform
> >> >> >>> >> > stubbornly uses whatever values it gets from
> >> >> >>> >> > vtkMath::vtkJacobiN?
> >> >> >>> >> > Is there a way to intercept the error message?
> >> >> >>> >> >
> >> >> >>> >> > I also noticed that if I apply ICP on A and B1 resulting in
> >> >> >>> >> > the
> >> >> >>> >> > error
> >> >> >>> >> > message, applying ICP on A and B2 (loading the same data
> file
> >> >> >>> >> > that
> >> >> >>> >> > gives
> >> >> >>> >> > me
> >> >> >>> >> > B1) works fine. I don't see vtkMath::vtkJacobiN using a
> random
> >> >> >>> >> > number or
> >> >> >>> >> > doing something similar, which makes me believe the error
> >> >> >>> >> > message
> >> >> >>> >> > I
> >> >> >>> >> > get
> >> >> >>> >> > is
> >> >> >>> >> > due to floating point accuracy? Even more a reason to think
> >> >> >>> >> > that
> >> >> >>> >> > 'giggling'
> >> >> >>> >> > might work (if I were able to intervene before the invalid
> >> >> >>> >> > transform
> >> >> >>> >> > is
> >> >> >>> >> > applied to my data...)
> >> >> >>> >> >
> >> >> >>> >> > Thanks - Maarten
> >> >> >>
> >> >> >>
> >> >> >
> >> >
> >> >
> >
> >
>
>
>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.vtk.org/pipermail/vtkusers/attachments/20120206/17cd2b26/attachment.htm>