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