Octave criterion for matrix singularity

classic Classic list List threaded Threaded
8 messages Options
Reply | Threaded
Open this post in threaded view
|

Octave criterion for matrix singularity

mabalenk
Dear all,

I'm repeatedly solving a linear system with a Vandermonde matrix. I'm using
a "backslash" operator:

  y = -W'/b'

At each iteration I'm increasing the matrix order by one. At some point
Octave issues a warning:

  warning: matrix singular to machine precision, rcond = 1.56508e-17

If possible, would you please direct me to the Octave source code, where
this decision is made? I would like to understand the criterion, that
triggers the warning. I found a post on StackOverflow that says MATLAB would
issue a similar warning, once rcond < 1e-12. I would like to know, what
threshold does Octave use and how is it defined. Thank you!

--
Best wishes,
Maxim



--
Sent from: http://octave.1599824.n4.nabble.com/Octave-General-f1599825.html


Reply | Threaded
Open this post in threaded view
|

Re: Octave criterion for matrix singularity

tmacchant
-- maxim.abalenkov

> Dear all,
>
> I'm repeatedly solving a linear system with a Vandermonde matrix. I'm using
> a "backslash" operator:
>
>   y = -W'/b'
>
> At each iteration I'm increasing the matrix order by one. At some point
> Octave issues a warning:
>
>   warning: matrix singular to machine precision, rcond = 1.56508e-17
>
> If possible, would you please direct me to the Octave source code, where
> this decision is made? I would like to understand the criterion, that
> triggers the warning. I found a post on StackOverflow that says MATLAB would
> issue a similar warning, once rcond < 1e-12. I would like to know, what
> threshold does Octave use and how is it defined. Thank you!
>
> --
> Best wishes,
> Maxim
>

You wrote
y=W’/b’

/ is slash.

Is the above a type miss?

Tatsuro


Reply | Threaded
Open this post in threaded view
|

Re: Octave criterion for matrix singularity

mabalenk
Yes, you are right. It is a typo. I meant "backslash":

  y = -W'\b';



--
Sent from: http://octave.1599824.n4.nabble.com/Octave-General-f1599825.html


Reply | Threaded
Open this post in threaded view
|

Re: Octave criterion for matrix singularity

Carlo de Falco-2
In reply to this post by mabalenk


> On 25 Feb 2019, at 19:34, mabalenk <[hidden email]> wrote:
>
> Dear all,
>
> I'm repeatedly solving a linear system with a Vandermonde matrix. I'm using
> a "backslash" operator:
>
>  y = -W'/b'
>
> At each iteration I'm increasing the matrix order by one. At some point
> Octave issues a warning:
>
>  warning: matrix singular to machine precision, rcond = 1.56508e-17
>
> If possible, would you please direct me to the Octave source code, where
> this decision is made? I would like to understand the criterion, that
> triggers the warning. I found a post on StackOverflow that says MATLAB would
> issue a similar warning, once rcond < 1e-12. I would like to know, what
> threshold does Octave use and how is it defined. Thank you!
>
> --
> Best wishes,
> Maxim



The code for solving A\b with A real double precision, full matrix is here :

 http://hg.savannah.gnu.org/hgweb/octave/file/14815cb62fbe/liboctave/array/dMatrix.cc#l1445

the condition being checked is

 rcond + 1.0 == 1.0

which means rcond is less or equal than eps(1.0)/2

HTH,
c.





Reply | Threaded
Open this post in threaded view
|

Re: Octave criterion for matrix singularity

Eddy Xiao
In reply to this post by mabalenk
mabalenk wrote
> At each iteration I'm increasing the matrix order by one. At some point
> Octave issues a warning:
>
>   warning: matrix singular to machine precision, rcond = 1.56508e-17
>
> If possible, would you please direct me to the Octave source code, where
> this decision is made? I would like to understand the criterion, that
> ...

You can search the source code by key words
$ grep -IR 'matrix singular to machine precision'
This direct you to file "liboctave/util/lo-array-errwarn.cc" and function
"warn_singular_matrix (double rcond)".

With this next key words, search by
$ grep -IR 'warn_singular_matrix'
get you to e.g. file "liboctave/array/dMatrix.cc". I would guess it is the
function "Matrix::fsolve" that accomplishes the work of "back slash". And
the key lines are

    volatile double rcond_plus_one = rcon + 1.0;
    if (rcond_plus_one == 1.0 || octave::math::isnan (rcon))

Alternatively, if you know that the operator "back slash" is called
"mldivide", then you can do "help mldivide" in Octave, it tells you:
    'mldivide' is a built-in function from the file
libinterp/corefcn/data.cc
From there you may trace down what's going on. Though probably not that easy
due to layers of abstraction.

One lesson I learn is: this "rcond" is actually (reciprocal) condition
number of L1 (or L-infinity) norm, instead of L2 norm. Thus it is easier and
faster to compute (by LAPACK routine dgecon()). You can get this rcond by
"rcond" function in Octave. And to my surprise, it is much more capable to
super-low reciprocal condition number, compared to cond() which is based on
SVD. e.g. rcond(vander((1:50)/50))=1.4262e-31,
1/cond(vander((1:50)/50))=4.8811e-20 (which is well under estimated).




--
Sent from: http://octave.1599824.n4.nabble.com/Octave-General-f1599825.html


Reply | Threaded
Open this post in threaded view
|

Re: Octave criterion for matrix singularity

mabalenk
In reply to this post by Carlo de Falco-2
Yes, it does help. Thank you. But I don't understand, where is the division
by 2 hidden. It seems to me, that rcond is compared to eps, the machine
precision.



--
Sent from: http://octave.1599824.n4.nabble.com/Octave-General-f1599825.html


Reply | Threaded
Open this post in threaded view
|

Re: Octave criterion for matrix singularity

Carlo de Falco-2


> On 26 Feb 2019, at 14:21, mabalenk <[hidden email]> wrote:
>
> Yes, it does help. Thank you. But I don't understand, where is the division
> by 2 hidden. It seems to me, that rcond is compared to eps, the machine
> precision.

Try by yourself :

>> eps
ans =    2.2204e-16
>> 1 + eps(1) == 1
ans = 0
>> 1 + eps(1)/2 == 1
ans = 1

or read the manual :

>> help eps

....

     More precisely, 'eps' is the relative spacing between any two
     adjacent numbers in the machine's floating point system.  This
     number is obviously system dependent.  On machines that support
     IEEE floating point arithmetic, 'eps' is approximately 2.2204e-16
     for double precision and 1.1921e-07 for single precision.

....

HTH,
c.




Reply | Threaded
Open this post in threaded view
|

Re: Octave criterion for matrix singularity

spoorti
In reply to this post by tmacchant
Thanks a lot. It helped a lot



-----
app development companies in Canada
mobile app development company in Russia



--
Sent from: http://octave.1599824.n4.nabble.com/Octave-General-f1599825.html