Community
cancel
Showing results for 
Search instead for 
Did you mean: 
swillmon
Beginner
165 Views

Matrix inversion with MKL LAPACK

Jump to solution
Test driving my new Composer XE 2011 with a simple FORTRAN code that computes the inverse of a n x n matrix using the getrf and getri calls to LAPACK. Is there any rational reason why the results of the inverted matrix differ from what MATLAB is giving me?
0 Kudos
1 Solution
mecej4
Black Belt
165 Views
With Matlab 7.5 (R2007b) I got results that agreed with your (and my) results from MKL Lapack95:

[bash]>> A=[   1.0658   -0.6313    0.0000    0.0000    0.0000    0.0000    0.0000
   -0.6313    4.2632   -1.8940    0.0000    0.0000    0.0000    0.0000
    0.0000   -1.8940    8.5265   -3.1567    0.0000    0.0000    0.0000
    0.0000    0.0000   -3.1567   23.2603   -3.6892    0.0000    0.0000
    0.0000    0.0000    0.0000   -3.6892   37.5680   -4.8192    0.0000
    0.0000    0.0000    0.0000    0.0000   -4.8192   47.5467   -5.9492
    0.0000    0.0000    0.0000    0.0000    0.0000   -5.9492   57.5255];
>> inv(A)

ans =

    1.0401    0.1719    0.0402    0.0055    0.0006    0.0001    0.0000
    0.1719    0.2902    0.0679    0.0094    0.0009    0.0001    0.0000
    0.0402    0.0679    0.1395    0.0192    0.0019    0.0002    0.0000
    0.0055    0.0094    0.0192    0.0463    0.0046    0.0005    0.0000
    0.0006    0.0009    0.0019    0.0046    0.0274    0.0028    0.0003
    0.0001    0.0001    0.0002    0.0005    0.0028    0.0216    0.0022
    0.0000    0.0000    0.0000    0.0000    0.0003    0.0022    0.0176

>> cond(A)

ans =

   65.3807

>> A-A'

ans =

     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0[/bash]

View solution in original post

5 Replies
mecej4
Black Belt
165 Views
Are you sure that the original matrix being inverted is the same in your Matlab and Fortran runs? How big is the matrix? Is it possible to give an example where the inverse is noticeably different?
swillmon
Beginner
165 Views
Here are the arrays I've been working with / results:

Array
1.0658 -0.6313 0.0000 0.0000 0.0000 0.0000 0.0000
-0.6313 4.2632 -1.8940 0.0000 0.0000 0.0000 0.0000
0.0000 -1.8940 8.5265 -3.1567 0.0000 0.0000 0.0000
0.0000 0.0000 -3.1567 23.2603 -3.6892 0.0000 0.0000
0.0000 0.0000 0.0000 -3.6892 37.5680 -4.8192 0.0000
0.0000 0.0000 0.0000 0.0000 -4.8192 47.5467 -5.9492
0.0000 0.0000 0.0000 0.0000 0.0000 -5.9492 57.5255

GETRI Array
1.0401 0.1719 0.0402 0.0055 0.0006 0.0001 0.0000
0.1719 0.2902 0.0679 0.0094 0.0009 0.0001 0.0000
0.0402 0.0679 0.1395 0.0192 0.0019 0.0002 0.0000
0.0055 0.0094 0.0192 0.0463 0.0046 0.0005 0.0000
0.0006 0.0009 0.0019 0.0046 0.0274 0.0028 0.0003
0.0001 0.0001 0.0002 0.0005 0.0028 0.0216 0.0022
0.0000 0.0000 0.0000 0.0000 0.0003 0.0022 0.0176

MATLAB INV(ARRAY)
1.0042 0.1114 0.0182 0.0025 0.0003 0.0000 0.0000
0.1114 0.1881 0.0308 0.0042 0.0004 0.0000 0.0000
0.0182 0.0308 0.0915 0.0126 0.0013 0.0001 0.0000
0.0025 0.0042 0.0126 0.0454 0.0045 0.0005 0.0000
0.0003 0.0004 0.0013 0.0045 0.0274 0.0028 0.0003
0.0000 0.0000 0.0001 0.0005 0.0028 0.0216 0.0022
0.0000 0.0000 0.0000 0.0000 0.0003 0.0022 0.0176

Thanks for taking a look!
styc
Beginner
165 Views
[plain]$ matlab
Warning: No display specified.  You will not be able to display graphics on the screen.

                            < M A T L A B  >
                  Copyright 1984-2008 The MathWorks, Inc.
                         Version 7.6.0.324 (R2008a)
                             February 10, 2008


  To get started, type one of these: helpwin, helpdesk, or demo.
  For product information, visit www.mathworks.com.

>> A=[    1.0658   -0.6313    0.0000    0.0000    0.0000    0.0000    0.0000
   -0.6313    4.2632   -1.8940    0.0000    0.0000    0.0000    0.0000
    0.0000   -1.8940    8.5265   -3.1567    0.0000    0.0000    0.0000
    0.0000    0.0000   -3.1567   23.2603   -3.6892    0.0000    0.0000
    0.0000    0.0000    0.0000   -3.6892   37.5680   -4.8192    0.0000
    0.0000    0.0000    0.0000    0.0000   -4.8192   47.5467   -5.9492
    0.0000    0.0000    0.0000    0.0000    0.0000   -5.9492   57.5255
]
A =
    1.0658   -0.6313         0         0         0         0         0
   -0.6313    4.2632   -1.8940         0         0         0         0
         0   -1.8940    8.5265   -3.1567         0         0         0
         0         0   -3.1567   23.2603   -3.6892         0         0
         0         0         0   -3.6892   37.5680   -4.8192         0
         0         0         0         0   -4.8192   47.5467   -5.9492
         0         0         0         0         0   -5.9492   57.5255
>> inv(A)
ans =
    1.0401    0.1719    0.0402    0.0055    0.0006    0.0001    0.0000
    0.1719    0.2902    0.0679    0.0094    0.0009    0.0001    0.0000
    0.0402    0.0679    0.1395    0.0192    0.0019    0.0002    0.0000
    0.0055    0.0094    0.0192    0.0463    0.0046    0.0005    0.0000
    0.0006    0.0009    0.0019    0.0046    0.0274    0.0028    0.0003
    0.0001    0.0001    0.0002    0.0005    0.0028    0.0216    0.0022
    0.0000    0.0000    0.0000    0.0000    0.0003    0.0022    0.0176
>> quit
$
[/plain]
mecej4
Black Belt
166 Views
With Matlab 7.5 (R2007b) I got results that agreed with your (and my) results from MKL Lapack95:

[bash]>> A=[   1.0658   -0.6313    0.0000    0.0000    0.0000    0.0000    0.0000
   -0.6313    4.2632   -1.8940    0.0000    0.0000    0.0000    0.0000
    0.0000   -1.8940    8.5265   -3.1567    0.0000    0.0000    0.0000
    0.0000    0.0000   -3.1567   23.2603   -3.6892    0.0000    0.0000
    0.0000    0.0000    0.0000   -3.6892   37.5680   -4.8192    0.0000
    0.0000    0.0000    0.0000    0.0000   -4.8192   47.5467   -5.9492
    0.0000    0.0000    0.0000    0.0000    0.0000   -5.9492   57.5255];
>> inv(A)

ans =

    1.0401    0.1719    0.0402    0.0055    0.0006    0.0001    0.0000
    0.1719    0.2902    0.0679    0.0094    0.0009    0.0001    0.0000
    0.0402    0.0679    0.1395    0.0192    0.0019    0.0002    0.0000
    0.0055    0.0094    0.0192    0.0463    0.0046    0.0005    0.0000
    0.0006    0.0009    0.0019    0.0046    0.0274    0.0028    0.0003
    0.0001    0.0001    0.0002    0.0005    0.0028    0.0216    0.0022
    0.0000    0.0000    0.0000    0.0000    0.0003    0.0022    0.0176

>> cond(A)

ans =

   65.3807

>> A-A'

ans =

     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0[/bash]

View solution in original post

swillmon
Beginner
165 Views
mecej4 - after your first reply, I went back and scrubbed my input data between what I was running w/ FORTRAN and what I put into MATLAB -- apples and oranges! Of course I get the same answer now that I'm using the same input data!

Thanks!
Reply