[Insight-users] [Fwd: SVD compiled with optimziation]
Karthik Krishnan
Karthik.Krishnan at kitware.com
Tue Dec 20 11:08:50 EST 2005
Did you ever get a chance to see if this happens when you replace the
-O3 flag in Utilities/vxl/netlib/Makefile to -O2 ?
[ To be more specific, I ended up manually editing the netlib/Makefile
to add
-fno-inline-functions to the already existing flags
-w -O3 -fPIC
and that took care of the issue.
And these flags were added only to the makefile rules of dsvdc.o,
csvdc.o, zsvdc.o, ssvdc.o ]
They are to prevent inlining of that function below, which I believe is
the miscreant somehow.. See below..
See the pasted email before.. I ran into a similar problem.... I did not
look into it further, but now that you mention it, I'll add it to the
bug tracker.
Thanks
karthik
/------------Mon Nov 28 20:30:25 EST
2005/--------------------------------------------------------------------------------------
I ran into the same problem today, using itk::ImageMomemtsCalculator.
(I'm using gcc40 on debian)
Turns out, that uses vnl_svd which uses netlib/dsvdc.c and all these are
related.
Funnily enough, turns out that is some numerical issue in netlib/dsvdc
when optimization is turned on.
Did you compile with the -O3 flag (usually the Release mode on CMake) ?
For now, I changed the Utilities/vxl/netlib/Makefile to -O2 on the build
line that compiles dsvdc.c
and everything was working great.
I haven't snooped around to see what exactly is causing the problem, but
there seem to be others in the VXL community who've run into similar
issues.
The following code snippet is from dsvdc.c and its cvs logs from the vxl
sources.
<snip>
----------------------------
revision 1.6
date: 2002/03/13 15:51:28; author: fsm; state: Exp; lines: +12 -2
hack to avoid excess precision
----------------------------
/*
* Calling this ensures that the operands are spilled to
* memory and thus avoids excessive precision when compiling
* for x86 with heavy optimization (gcc). It is better to do
* this than to turn on -ffloat-store.
*/
static int fsm_ieee_doubles_differ(double *x, double *y)
{
return *x != *y;
}
</snip>
In any case comparing two doubles for equality sounds questionable.
What's bothersome is that the inverse with -O3 builds is not marginally
wrong, but quite wrong for a 3x3 rotation matrix ?
Has anybody run into this issue. Is it limited only to gcc ?
I've changed my nightly builds to reflect this. All the associated (23
tests or so fail now on the release builds)
http://www.itk.org/Testing/Sites/Sabbath.kitware/zRel22-Linux-gcc40/20051220-0100-Nightly/Test.html
Thanks
karthik
-------------------------------------------------------------------------------------------------------
Suyash P. Awate wrote:
> Dear ITK Users,
>
> I found the cause of the error. It deals with my build of ITK in the
> SUSE 10 environment that uses gcc 4.0 by default. When I build itk and
> the applications on SUSE 9.3 using gcc 3.4 (definitely not gcc 4.0) or
> Fedora Core 2 (again not using gcc 4.0), I did NOT get any vnl error.
>
> I confirmed that precision is NOT an issue because the condition
> number for the matrix with truncated values is still very close to 1.
> I found this via matlab. So numerical "instabilities" in
> inverse-finding should not occur.
>
> So it seems that there is something strange about the way the vnl
> code/libraries interact with gcc 4.0.
>
> Thanks again for all your help.
>
> Best regards,
> Suyash.
>
> ----- Original Message -----
> *From:* Miller, James V (Research) <mailto:millerjv at crd.ge.com>
> *To:* Suyash P. Awate <mailto:suyash at cs.utah.edu> ; jiafucang
> <mailto:jiafucang at asisz.com> ; jjomier at cs.unc.edu
> <mailto:jjomier at cs.unc.edu> ; insight-users at itk.org
> <mailto:insight-users at itk.org>
> *Sent:* Friday, November 11, 2005 11:24 AM
> *Subject:* RE: [Insight-users] Re: On incorrect matrix inverse
>
> Suyash,
>
> You might not have enough precision in the matrix values that you
> are reading from the command line. Try generating the values with
> more
> precision from whatever tool is generating them. Try going out to
> about 12 decimal places.
>
> I have seen this before where one program computes a matrix that
> is supposed to be invertible, writes the matrix in ascii, another
> program reads that matrix and it is no longer invertible.
>
> Jim
>
>
>
> -----Original Message-----
> *From:* insight-users-bounces+millerjv=crd.ge.com at itk.org
> <mailto:insight-users-bounces+millerjv=crd.ge.com at itk.org>
> [mailto:insight-users-bounces+millerjv=crd.ge.com at itk.org]*On
> Behalf Of *Suyash P. Awate
> *Sent:* Friday, November 11, 2005 1:06 PM
> *To:* jiafucang; jjomier at cs.unc.edu; insight-users at itk.org
> *Subject:* [Insight-users] Re: On incorrect matrix inverse
>
>
> Hi everybody,
>
> Thanks for the generous help.
>
> I was also trying to get the inverse concerning with the
> registration application, like Fucang (see email attached).
>
> Here is my piece of code that gave me the wrong result:
>
> ------------
> const unsigned int Dimension = 3;
> typedef double PixelType;
> typedef itk::Matrix <PixelType, Dimension, Dimension> MatrixType;
> MatrixType matrix;
> matrix.Fill (0);
> matrix (0,0)= atof (argv[4]);
> matrix (0,1)= atof (argv[5]);
> matrix (0,2)= atof (argv[6]);
> matrix (1,0)= atof (argv[7]);
> matrix (1,1)= atof (argv[8]);
> matrix (1,2)= atof (argv[9]);
> matrix (2,0)= atof (argv[10]);
> matrix (2,1)= atof (argv[11]);
> matrix (2,2)= atof (argv[12]);
> vnl_matrix_fixed <PixelType, Dimension, Dimension> inverse;
> inverse= matrix.GetInverse();
> -----------
> matrix
> 0.99925 0.0371382 -0.0109982
> -0.0373772 0.99905 -0.0223916
> 0.0101562 0.0227859 0.999689
>
> warnings:
> InsightToolkit-2.2.0/Utilities/vxl/core/vnl/algo/vnl_svd.txx:
> suspicious return value (2) from SVDC
> InsightToolkit-2.2.0/Utilities/vxl/core/vnl/algo/vnl_svd.txx:
> M is 3x3
> M = [ ...
> 0.9992500000000 0.0371382000000 -0.0109982000000
> -0.0373772000000 0.9990500000000 -0.0223916000000
> 0.0101562000000 0.0227859000000 0.9996890000000 ]
>
> matrix inverse using vnl in itk (incorrect)
> -0.982557 -0.0644525 -0.174434
> -0.0677511 0.997618 0.0130162
> -0.17318 -0.0246073 0.984583
>
> matrix inverse in matlab (correct)
> 0.9992 -0.0374 0.0102
> 0.0371 0.9991 0.0228
> -0.0110 -0.0224 0.9997
> ------------
>
> I had also tried using the AffineTransform but for some reason
> could not get that to work and so I wrote this simple matrix
> code. Perhaps that was due to the same problem that Fucang
> reported.
> Currently, however, my problem is different. As we can see I
> am not using any Transform class at all. Moreover I am getting
> warnings from the vnl_svd class, so vnl is definitely involved
> here.
>
> Julien: I am giving 9 numbers (truncated: 0.99925 0.0371382
> -0.0109982 -0.0373772 0.99905 -0.0223916 0.0101562 0.0227859
> 0.999689) to the matrix in itk and also the same numbers to
> the matrix in matlab. I got these numbers from the output of a
> separate program.
> Just FYI, I have compiled ITK on Suse 10.0 using GCC 4.0. I
> wonder if that has anything to do with this !
>
> Best regards,
> Suyash.
>
>
>
> ----- Original Message -----
> *From:* jiafucang <mailto:jiafucang at asisz.com>
> *To:* suyash at cs.utah.edu <mailto:suyash at cs.utah.edu> ;
> jjomier at cs.unc.edu <mailto:jjomier at cs.unc.edu>
> *Sent:* Friday, November 11, 2005 7:10 AM
> *Subject:* On incorrect matrix inverse
>
> Hello, Suyash, Jomier,
>
> I have reported similar problem on
> http://public.kitware.com/pipermail/insight-users/2005-October/015441.html
>
>
> and added a bug 2450
>
> http://www.itk.org/Bug/bug.php?op=show&bugid=2450&pos=24
> <http://www.itk.org/Bug/bug.php?op=show&bugid=2450&pos=24>
>
> I think it is because ITK does not computer inverse when
> setparameters used, and not VNL problem.
>
> HTH
>
> Fucang
>
>------------------------------------------------------------------------
>
>_______________________________________________
>Insight-users mailing list
>Insight-users at itk.org
>http://www.itk.org/mailman/listinfo/insight-users
>
More information about the Insight-users
mailing list