Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

Thread Subject:
inconsistent results of matrix determinant

Subject: inconsistent results of matrix determinant

From: Deepak

Date: 27 Aug, 2007 21:55:52

Message: 1 of 12

I have a matrix intA that changes extremely slightly on each
iteration and I have to calculate

stabval = det(eye(size(intA)) + expm(intA))

To simulate this problem, I added an extremely small random
matrix to intA and then I calculate this determinant. This
is what I get:

K>> stabval = det(eye(size(intA)) + expm(intA +
1e-12*rand(size(intA))))

stabval =

  4.8858e+064

K>> stabval = det(eye(size(intA)) + expm(intA +
1e-12*rand(size(intA))))

stabval =

 -2.6021e+063

K>> stabval = det(eye(size(intA)) + expm(intA +
1e-12*rand(size(intA))))

stabval =

 -1.7182e+064

K>> stabval = det(eye(size(intA)) + expm(intA +
1e-12*rand(size(intA))))

stabval =

 -2.7059e+064

K>> stabval = det(eye(size(intA)) + expm(intA +
1e-12*rand(size(intA))))

stabval =

 -2.7574e+064

K>> stabval = det(eye(size(intA)) + expm(intA +
1e-12*rand(size(intA))))

stabval =

 -3.4601e+064

K>> stabval = det(eye(size(intA)) + expm(intA +
1e-12*rand(size(intA))))

stabval =

  1.7884e+064

K>> stabval = det(eye(size(intA)) + expm(intA +
1e-12*rand(size(intA))))

stabval =

  2.0292e+064

K>> stabval = det(eye(size(intA)) + expm(intA +
1e-12*rand(size(intA))))

stabval =

 -5.4093e+064

K>> stabval = det(eye(size(intA)) + expm(intA +
1e-12*rand(size(intA))))

stabval =

  4.6570e+063

K>> stabval = det(eye(size(intA)) + expm(intA +
1e-12*rand(size(intA))))

stabval =

 -8.1755e+063


The elements of intA are:

intA =

  1.0e+003 *

         0 0.0021 -0.0015 0 0 -0.0000
    0.0002 0 -0.0000 0 0 -0.0000
   -3.3873 -0.0004 0 0 0 0
   -0.0001 0.0005 0 0 0 -0.0001
    0.0005 0.0001 0 0 0 0.0005
         0 0 0.0010 0 0 0


Why does the value of the determinant vary so wildly for a
negligible change in the matrix?

Subject: inconsistent results of matrix determinant

From: Yang Zhang

Date: 28 Aug, 2007 15:35:37

Message: 2 of 12

"Deepak " <deepak.trivedi@gmail.com> wrote in message
<favh98$8s8$1@fred.mathworks.com>...
> I have a matrix intA that changes extremely slightly on each
> iteration and I have to calculate
>
> stabval = det(eye(size(intA)) + expm(intA))
>
> To simulate this problem, I added an extremely small random
> matrix to intA and then I calculate this determinant. This
> is what I get:
>
> K>> stabval = det(eye(size(intA)) + expm(intA +
> 1e-12*rand(size(intA))))
>
> stabval =
>
> 4.8858e+064
>
> K>> stabval = det(eye(size(intA)) + expm(intA +
> 1e-12*rand(size(intA))))
>
> stabval =
>
> -2.6021e+063
>
> K>> stabval = det(eye(size(intA)) + expm(intA +
> 1e-12*rand(size(intA))))
>
> stabval =
>
> -1.7182e+064
>
> K>> stabval = det(eye(size(intA)) + expm(intA +
> 1e-12*rand(size(intA))))
>
> stabval =
>
> -2.7059e+064
>
> K>> stabval = det(eye(size(intA)) + expm(intA +
> 1e-12*rand(size(intA))))
>
> stabval =
>
> -2.7574e+064
>
> K>> stabval = det(eye(size(intA)) + expm(intA +
> 1e-12*rand(size(intA))))
>
> stabval =
>
> -3.4601e+064
>
> K>> stabval = det(eye(size(intA)) + expm(intA +
> 1e-12*rand(size(intA))))
>
> stabval =
>
> 1.7884e+064
>
> K>> stabval = det(eye(size(intA)) + expm(intA +
> 1e-12*rand(size(intA))))
>
> stabval =
>
> 2.0292e+064
>
> K>> stabval = det(eye(size(intA)) + expm(intA +
> 1e-12*rand(size(intA))))
>
> stabval =
>
> -5.4093e+064
>
> K>> stabval = det(eye(size(intA)) + expm(intA +
> 1e-12*rand(size(intA))))
>
> stabval =
>
> 4.6570e+063
>
> K>> stabval = det(eye(size(intA)) + expm(intA +
> 1e-12*rand(size(intA))))
>
> stabval =
>
> -8.1755e+063
>
>
> The elements of intA are:
>
> intA =
>
> 1.0e+003 *
>
> 0 0.0021 -0.0015 0 0 -0.0000
> 0.0002 0 -0.0000 0 0 -0.0000
> -3.3873 -0.0004 0 0 0 0
> -0.0001 0.0005 0 0 0 -0.0001
> 0.0005 0.0001 0 0 0 0.0005
> 0 0 0.0010 0 0 0
>
>
> Why does the value of the determinant vary so wildly for a
> negligible change in the matrix?
>
>

Hi,

Are you doing stability and controllability analysis? If so,
you can notice that your system is not very robust.
you can check the following value norm(I-(intA+noise)*inv(intA))
A small value indicates a robust system.
HTH
Yang

Subject: inconsistent results of matrix determinant

From: Deepak

Date: 28 Aug, 2007 19:43:41

Message: 3 of 12

Yang,

Thanks for your reply. I am not doing controllability
analysis. Please note that intA is not invertible in this case.

Regards
Deepak
"Yang Zhang" <zhyang99@hotmail.com> wrote in message
<fb1fc9$986$1@fred.mathworks.com>...
> "Deepak " <deepak.trivedi@gmail.com> wrote in message
> <favh98$8s8$1@fred.mathworks.com>...
> > I have a matrix intA that changes extremely slightly on each
> > iteration and I have to calculate
> >
> > stabval = det(eye(size(intA)) + expm(intA))
> >
> > To simulate this problem, I added an extremely small random
> > matrix to intA and then I calculate this determinant. This
> > is what I get:
> >
> > K>> stabval = det(eye(size(intA)) + expm(intA +
> > 1e-12*rand(size(intA))))
> >
> > stabval =
> >
> > 4.8858e+064
> >
> > K>> stabval = det(eye(size(intA)) + expm(intA +
> > 1e-12*rand(size(intA))))
> >
> > stabval =
> >
> > -2.6021e+063
> >
> > K>> stabval = det(eye(size(intA)) + expm(intA +
> > 1e-12*rand(size(intA))))
> >
> > stabval =
> >
> > -1.7182e+064
> >
> > K>> stabval = det(eye(size(intA)) + expm(intA +
> > 1e-12*rand(size(intA))))
> >
> > stabval =
> >
> > -2.7059e+064
> >
> > K>> stabval = det(eye(size(intA)) + expm(intA +
> > 1e-12*rand(size(intA))))
> >
> > stabval =
> >
> > -2.7574e+064
> >
> > K>> stabval = det(eye(size(intA)) + expm(intA +
> > 1e-12*rand(size(intA))))
> >
> > stabval =
> >
> > -3.4601e+064
> >
> > K>> stabval = det(eye(size(intA)) + expm(intA +
> > 1e-12*rand(size(intA))))
> >
> > stabval =
> >
> > 1.7884e+064
> >
> > K>> stabval = det(eye(size(intA)) + expm(intA +
> > 1e-12*rand(size(intA))))
> >
> > stabval =
> >
> > 2.0292e+064
> >
> > K>> stabval = det(eye(size(intA)) + expm(intA +
> > 1e-12*rand(size(intA))))
> >
> > stabval =
> >
> > -5.4093e+064
> >
> > K>> stabval = det(eye(size(intA)) + expm(intA +
> > 1e-12*rand(size(intA))))
> >
> > stabval =
> >
> > 4.6570e+063
> >
> > K>> stabval = det(eye(size(intA)) + expm(intA +
> > 1e-12*rand(size(intA))))
> >
> > stabval =
> >
> > -8.1755e+063
> >
> >
> > The elements of intA are:
> >
> > intA =
> >
> > 1.0e+003 *
> >
> > 0 0.0021 -0.0015 0 0 -0.0000
> > 0.0002 0 -0.0000 0 0 -0.0000
> > -3.3873 -0.0004 0 0 0 0
> > -0.0001 0.0005 0 0 0 -0.0001
> > 0.0005 0.0001 0 0 0 0.0005
> > 0 0 0.0010 0 0 0
> >
> >
> > Why does the value of the determinant vary so wildly for a
> > negligible change in the matrix?
> >
> >
>
> Hi,
>
> Are you doing stability and controllability analysis? If so,
> you can notice that your system is not very robust.
> you can check the following value
norm(I-(intA+noise)*inv(intA))
> A small value indicates a robust system.
> HTH
> Yang

Subject: inconsistent results of matrix determinant

From: Tim Davis

Date: 28 Aug, 2007 20:52:43

Message: 4 of 12

"Deepak " <deepak.trivedi@gmail.com> wrote in message ...
> I have a matrix intA that changes extremely slightly on each
> iteration and I have to calculate
>
> stabval = det(eye(size(intA)) + expm(intA))

"det" is fundamentally a very very numerically sensitive
thing to compute. You should "never" use det. Try to find
another method for computing what you want.

Subject: inconsistent results of matrix determinant

From: Tim Davis

Date: 28 Aug, 2007 21:25:33

Message: 5 of 12

"Tim Davis" <davis@cise.ufl.edu> wrote in message
<fb21ur$akl$1@fred.mathworks.com>...
> "Deepak " <deepak.trivedi@gmail.com> wrote in message ...
> > I have a matrix intA that changes extremely slightly on each
> > iteration and I have to calculate
> >
> > stabval = det(eye(size(intA)) + expm(intA))
>
> "det" is fundamentally a very very numerically sensitive
> thing to compute. You should "never" use det. Try to find
> another method for computing what you want.

For example, det(expm(A)) equals exp(trace(A)). That only
gets you part-way, though...

Subject: inconsistent results of matrix determinant

From: Yang Zhang

Date: 28 Aug, 2007 21:50:20

Message: 6 of 12

"Deepak " <deepak.trivedi@gmail.com> wrote in message
<fb1ttd$217$1@fred.mathworks.com>...
> Yang,
>
> Thanks for your reply. I am not doing controllability
> analysis. Please note that intA is not invertible in this
case.

Deepak,

I think the main reason is intA:
expm(A) =
  4.5405e+030 1.3430e+029 -9.5545e+028 0
    0 0
  1.2739e+028 3.7680e+026 -2.6807e+026 0
    0 0
 -2.1576e+032 -6.3817e+030 4.5401e+030 0
    0 0
 -2.0342e+027 -6.0167e+025 4.2805e+025 1.0000e+000
    0 -1.0000e-001
  1.0636e+028 3.1458e+026 -2.2380e+026 0
1.0000e+000 5.0000e-001
 -3.0268e+030 -8.9525e+028 6.3691e+028 0
    0 1.0000e+000

which is already violate the working precision of MATLAB.

rank(expm(A)) = 1 and
eig(expm(A)) =
  1.0000e+000
  1.0000e+000
  1.0000e+000
  9.0810e+030
 -1.3740e+015
 -1.8184e+014

You may notice there are 6 non zero eigenvalues but why
MATLAB say the rank is 1. The working precision of MATLAB is:
tol = max(size(A)) * eps(norm(A)) = 2.1e17
All the values below this value is treated as imprecision by
MATLAB. When expm(A) is so inaccurate, you cannot expect a
robust det calculation result.

Is it make sense?

Yang

>
> Regards
> Deepak
> "Yang Zhang" <zhyang99@hotmail.com> wrote in message
> <fb1fc9$986$1@fred.mathworks.com>...
> > "Deepak " <deepak.trivedi@gmail.com> wrote in message
> > <favh98$8s8$1@fred.mathworks.com>...
> > > I have a matrix intA that changes extremely slightly
on each
> > > iteration and I have to calculate
> > >
> > > stabval = det(eye(size(intA)) + expm(intA))
> > >
> > > To simulate this problem, I added an extremely small
random
> > > matrix to intA and then I calculate this determinant. This
> > > is what I get:
> > >
> > > K>> stabval = det(eye(size(intA)) + expm(intA +
> > > 1e-12*rand(size(intA))))
> > >
> > > stabval =
> > >
> > > 4.8858e+064
> > >
> > > K>> stabval = det(eye(size(intA)) + expm(intA +
> > > 1e-12*rand(size(intA))))
> > >
> > > stabval =
> > >
> > > -2.6021e+063
> > >
> > > K>> stabval = det(eye(size(intA)) + expm(intA +
> > > 1e-12*rand(size(intA))))
> > >
> > > stabval =
> > >
> > > -1.7182e+064
> > >
> > > K>> stabval = det(eye(size(intA)) + expm(intA +
> > > 1e-12*rand(size(intA))))
> > >
> > > stabval =
> > >
> > > -2.7059e+064
> > >
> > > K>> stabval = det(eye(size(intA)) + expm(intA +
> > > 1e-12*rand(size(intA))))
> > >
> > > stabval =
> > >
> > > -2.7574e+064
> > >
> > > K>> stabval = det(eye(size(intA)) + expm(intA +
> > > 1e-12*rand(size(intA))))
> > >
> > > stabval =
> > >
> > > -3.4601e+064
> > >
> > > K>> stabval = det(eye(size(intA)) + expm(intA +
> > > 1e-12*rand(size(intA))))
> > >
> > > stabval =
> > >
> > > 1.7884e+064
> > >
> > > K>> stabval = det(eye(size(intA)) + expm(intA +
> > > 1e-12*rand(size(intA))))
> > >
> > > stabval =
> > >
> > > 2.0292e+064
> > >
> > > K>> stabval = det(eye(size(intA)) + expm(intA +
> > > 1e-12*rand(size(intA))))
> > >
> > > stabval =
> > >
> > > -5.4093e+064
> > >
> > > K>> stabval = det(eye(size(intA)) + expm(intA +
> > > 1e-12*rand(size(intA))))
> > >
> > > stabval =
> > >
> > > 4.6570e+063
> > >
> > > K>> stabval = det(eye(size(intA)) + expm(intA +
> > > 1e-12*rand(size(intA))))
> > >
> > > stabval =
> > >
> > > -8.1755e+063
> > >
> > >
> > > The elements of intA are:
> > >
> > > intA =
> > >
> > > 1.0e+003 *
> > >
> > > 0 0.0021 -0.0015 0 0
-0.0000
> > > 0.0002 0 -0.0000 0 0
-0.0000
> > > -3.3873 -0.0004 0 0 0
     0
> > > -0.0001 0.0005 0 0 0
-0.0001
> > > 0.0005 0.0001 0 0 0
0.0005
> > > 0 0 0.0010 0 0
     0
> > >
> > >
> > > Why does the value of the determinant vary so wildly for a
> > > negligible change in the matrix?
> > >
> > >
> >
> > Hi,
> >
> > Are you doing stability and controllability analysis? If so,
> > you can notice that your system is not very robust.
> > you can check the following value
> norm(I-(intA+noise)*inv(intA))
> > A small value indicates a robust system.
> > HTH
> > Yang
>

Subject: inconsistent results of matrix determinant

From: Deepak

Date: 28 Aug, 2007 22:36:55

Message: 7 of 12

My main aim is to find out whether this determinant is
positive or negative. I am not so interested in its value.
Is there any way to find just that?

Regards
Deepak

"Tim Davis" <davis@cise.ufl.edu> wrote in message
<fb23sc$cka$1@fred.mathworks.com>...
> "Tim Davis" <davis@cise.ufl.edu> wrote in message
> <fb21ur$akl$1@fred.mathworks.com>...
> > "Deepak " <deepak.trivedi@gmail.com> wrote in message ...
> > > I have a matrix intA that changes extremely slightly
on each
> > > iteration and I have to calculate
> > >
> > > stabval = det(eye(size(intA)) + expm(intA))
> >
> > "det" is fundamentally a very very numerically sensitive
> > thing to compute. You should "never" use det. Try to find
> > another method for computing what you want.
>
> For example, det(expm(A)) equals exp(trace(A)). That only
> gets you part-way, though...

Subject: inconsistent results of matrix determinant

From: Yang Zhang

Date: 29 Aug, 2007 22:09:16

Message: 8 of 12

"Deepak " <deepak.trivedi@gmail.com> wrote in message
<fb2826$ppi$1@fred.mathworks.com>...
> My main aim is to find out whether this determinant is
> positive or negative. I am not so interested in its value.
> Is there any way to find just that?
>
> Regards
> Deepak

Deepak,

Are there any possibility to use analytical method instead
of numerical method in your problem? Or can you rescale your
inputs/outputs such that the A matrix will behave better?

Yang


>
> "Tim Davis" <davis@cise.ufl.edu> wrote in message
> <fb23sc$cka$1@fred.mathworks.com>...
> > "Tim Davis" <davis@cise.ufl.edu> wrote in message
> > <fb21ur$akl$1@fred.mathworks.com>...
> > > "Deepak " <deepak.trivedi@gmail.com> wrote in message ...
> > > > I have a matrix intA that changes extremely slightly
> on each
> > > > iteration and I have to calculate
> > > >
> > > > stabval = det(eye(size(intA)) + expm(intA))
> > >
> > > "det" is fundamentally a very very numerically sensitive
> > > thing to compute. You should "never" use det. Try to
find
> > > another method for computing what you want.
> >
> > For example, det(expm(A)) equals exp(trace(A)). That only
> > gets you part-way, though...
>

Subject: inconsistent results of matrix determinant

From: Deepak

Date: 31 Aug, 2007 00:59:40

Message: 9 of 12

Apparently, the problem is with expm. Consider this:



intA =

         0 0.6620 -2.6724 0 0 0
    0.0736 0 0 0 0 0
 -190.3624 -0.0002 0 0 0 0
   -0.0110 0.0992 0 0 0 -0.0111
    0.0992 0.0110 0 0 0 0.1000
         0 0 1.0000 0 0 0

>> expm(intA)*expm(-intA)

ans =

  1.0e+005 *

    0.1931 -0.0054 0.0237 0 0 0
    0.0006 -0.0000 0.0001 0 0 0
   -1.5066 0.0417 -0.1852 0 0 0
   -0.0001 0.0000 -0.0000 0.0000 0 0
    0.0005 -0.0000 0.0001 0 0.0000 0
   -0.0701 0.0019 -0.0086 0 0 0.0000

expm(intA)*expm(-intA) should have evaluated to I.

"Yang Zhang" <zhyang99@hotmail.com> wrote in message
<fb4qqc$cis$1@fred.mathworks.com>...
> "Deepak " <deepak.trivedi@gmail.com> wrote in message
> <fb2826$ppi$1@fred.mathworks.com>...
> > My main aim is to find out whether this determinant is
> > positive or negative. I am not so interested in its value.
> > Is there any way to find just that?
> >
> > Regards
> > Deepak
>
> Deepak,
>
> Are there any possibility to use analytical method instead
> of numerical method in your problem? Or can you rescale your
> inputs/outputs such that the A matrix will behave better?
>
> Yang
>
>
> >
> > "Tim Davis" <davis@cise.ufl.edu> wrote in message
> > <fb23sc$cka$1@fred.mathworks.com>...
> > > "Tim Davis" <davis@cise.ufl.edu> wrote in message
> > > <fb21ur$akl$1@fred.mathworks.com>...
> > > > "Deepak " <deepak.trivedi@gmail.com> wrote in
message ...
> > > > > I have a matrix intA that changes extremely slightly
> > on each
> > > > > iteration and I have to calculate
> > > > >
> > > > > stabval = det(eye(size(intA)) + expm(intA))
> > > >
> > > > "det" is fundamentally a very very numerically sensitive
> > > > thing to compute. You should "never" use det. Try to
> find
> > > > another method for computing what you want.
> > >
> > > For example, det(expm(A)) equals exp(trace(A)). That only
> > > gets you part-way, though...
> >
>

Subject: inconsistent results of matrix determinant

From: John D'Errico

Date: 31 Aug, 2007 08:14:50

Message: 10 of 12

"Deepak " <deepak.trivedi@gmail.com> wrote in message <fb7p5s$jh4
$1@fred.mathworks.com>...
> Apparently, the problem is with expm. Consider this:

Look at expm(A). If we can write A = P*D*P', where
D is diagonal, then expm(A) = P*expm(D)*P'. And
the matrix exponential of a diagonal matrix is given
by diag(exp(diag(A))). This is one of the methods
for expm that you would find in "19 Dubious Ways...".
(There are reasons why this is not the method used
inside expm.)

Now, look at your matrix. What are the eigenvalues?

eig(intA)
ans =
            0
            0
            0
       22.556
      -22.556
  -7.7319e-08

What is exp(22.556)? exp(-22.556)? Can you add
these two numbers together and get a result that
is different from exp(22.556)?

(exp(22.556) + exp(-22.556)) == exp(22.556)
ans =
     1

In fact, in the context of a matrix exponential and
double precision floating point arithmetic as
implemented in Matlab, this is a moderately nasty
problem. Just because a problem is mathematically
simple to write down on paper, this does not mean
it is numerically simple to solve, at least not
accurately.

The world is full of people who wail that "the system"
is inadequate to solve their problem, when formulated
as they choose. So rather than find a way to work
within the system to solve a reformulated version
of their problem, they let themselves be stopped in
their tracks, insisting that the problem is the system.
(Yes, you might argue that the system should be
improved, but until that happens, this gets you
nowhere.)

Suppose that, given double precision, these same
people demand quad precision. When given quad
precision, these same individuals then just formulate
larger, nastier problems. And when exponentials are
involved, those few extra digits get eaten up very
quickly. Do I hear a call for oct precision?

John

Subject: inconsistent results of matrix determinant

From: Tim Davis

Date: 31 Aug, 2007 10:19:42

Message: 11 of 12

"John D'Errico" <woodchips@rochester.rr.com> wrote in ...
> This is one of the methods
> for expm that you would find in "19 Dubious Ways...".
> (There are reasons why this is not the method used
> inside expm.)

John's right ... expm is not trivial to compute. See:

Nineteen Dubious Ways to Compute the Exponential of a
Matrix, Cleve Moler, Charles Van Loan, SIAM Review, Vol. 20,
No. 4 (Oct., 1978), pp. 801-836. The abstract states:

"In principle, the exponential of a matrix could be computed
in many ways. Methods involving approximation theory,
differential equations, the matrix eigenvalues, and the
matrix characteristic polynomial have been proposed. In
practice, consideration of computational stability and
efficiency indicates that some of the methods are preferable
to others, but that none are completely satisfactory."

> Do I hear a call for oct precision?
>
> John

Since the matrix is 6-by-6, you might need 6*8 byte
precision, or maybe more ... would that be octoquadragenuple
precision :-) ?

If you're desperate, your might try expm(sym(A)). I gave it
a try just for fun ... but it's still working ...

Subject: inconsistent results of matrix determinant

From: carlos lopez

Date: 5 Sep, 2007 15:03:24

Message: 12 of 12

"Tim Davis" <davis@cise.ufl.edu> wrote in message
<fb8pvt$skv$1@fred.mathworks.com>...
> Since the matrix is 6-by-6, you might need 6*8 byte
> precision, or maybe more ... would that be octoquadragenuple
> precision :-) ?
As a last resort, you might consider the multiple precision
toolbox by Ben Barrowes (available at the FEX) to use
arbitrary precision.
However, as John pointed out, using such a tool does not
preclude to think about the reformulation of the problem.
Regards
Carlos

Tags for this Thread

What are tags?

A tag is like a keyword or category label associated with each thread. Tags make it easier for you to find threads of interest.

Anyone can tag a thread. Tags are public and visible to everyone.

Contact us