Octave hangup for expm with non-finite arguments

Ben Abbott bpabbott at mac.com
Thu Jan 17 07:38:19 CST 2008


On Jan 17, 2008, at 7:28 AM, Rolf Fabian wrote:
>
> Rolf Fabian wrote:
>>
>> octave-3.0.0.exe:40> expm(y=[-Inf,0;0,-Inf])      # FAIL
>> HANGUP : Had to press Ctrl-C to get back Octave prompt
>>
>> In order to check what we can expect,
>> y might be approximated by using realmax.
>>
>> octave-3.0.0.exe:41> z = -realmax*eye(2)
>> x =
>>  -1.7977e+308  -0.0000e+000
>>  -0.0000e+000  -1.7977e+308
>>
>> octave-3.0.0.exe:42> expm(z)
>> ans =
>>   1   0
>>   0   1
>>
>> Hmm ...
>> expm ( - Inf ) = 0   and
>> expm ( [-Inf,0;0,-Inf] ) = eye(2) ?
>> Is this consistent ?
>> Have to think about this point.
>>

Matlab gives

 >> expm([-Inf,0;0,-Inf])
Warning: Matrix is singular, close to singular or badly scaled.
          Results may be inaccurate. RCOND = NaN.
 > In expm>PadeApproximantOfDegree at 123
   In expm at 39

ans =

    NaN   NaN
    NaN   NaN

 >>

> I've analyzed this issue in the meantime.
>
> expm ( [-Inf,0;0,-Inf] ) should return zeros(2), and not ones(2) !
>
> More generally for arguments  X of expm with repeated -Inf
> on the main diagonal ( and zeros elsewhere )
> X = diag( -Inf*ones(1,n) )
> ===========================================
> expm( X ) should return zeros(n)  for all dimensions n>0 .
> ===========================================
>
> Thus X can be considered to be the matrix logarithm of
> zeros(n) matrices.
> ==============================
> logm( zeros(n) ) =  diag( -Inf*ones(1,n) )
> ==============================
>
> Proof:
> This definition is reasonable, because it also
> results from the limiting process outlined
> below.
>
> I figured out the following general argument scaling
> formula for the matrix logarithm function
> ( non-zero real or complex scaling factor s required ):
>
> ==================================
> logm ( s*x ) = logm ( x ) +  log (s) .* eye (n)
> ==================================
> ( it may be checked with any arbitrary non-zero x
>  but keep in mind that Octave's logm(x) may fail
>  if your selection of x is 'defective'. E.g. for
>  random matrices this is extremely unlikely. )
>
> If we apply the general scaling formula to the particular
> matrix x = eye(n) and look at the limiting process
> s -> 0, we'll get
>
> lim s-> 0   logm ( eye(n)*s ) =  logm ( eye(n) ) + log (s) .* eye(n)
>            =  logm ( zeros(n) ) =   zeros(n) + log(0).*eye(n) )  =
>            =  -Inf .* diag( ones(1,n) )             q.e.d
>
> Alignment with ML ?
> Can anybody with access to MatLab please check what
> it returns for  the following code :
>
> -Inf * eye ( 2 )
> expm(  [ -Inf, 0; 0, -Inf ] )
> logm ( zeros (2 ) )

 >> logm(zeros(2))
Warning: Principal matrix logarithm is not defined for A with
          nonpositive real eigenvalues. A non-principal matrix
          logarithm is returned.
 > In funm at 156
   In logm at 27

ans =

    NaN   NaN
    NaN   NaN

> Thanks in advance.
>

Regarding logm(0) being -Inf, as that unique?

What I'm asking is, although exp(-Inf) = 0, is it proper to conclude  
that log(0) must be -Inf?

It would appear to me that there is no manner to determine a unique  
solution. All we can be sure of is that the real part must be negative  
and infinite, but what of the imaginary part? ... it can take any  
finite or infinite value correct?

I'm not familiar with how the IEEE defines such things, but it would  
appear to me to be a NaN, but what does the IEEE say?

Ben




More information about the Bug-octave mailing list