Groups | Search | Server Info | Keyboard shortcuts | Login | Register [http] [https] [nntp] [nntps]

Groups > comp.lang.python > #71252

Re: How can this assert() ever trigger?

Show all headers | View raw

On 05/10/2014 08:24 AM, Albert van der Horst wrote:
> I have the following code for calculating the determinant of
> a matrix. It works inasfar that it gives the same result as an
> octave program on a same matrix.
>
> / ----------------------------------------------------------------
>
> def determinant( mat ):
>      ''' Return the determinant of the n by n matrix mat
>          i row j column
>          Destroys mat ! '''
>      #print "getting determinat of", mat
>      n=len(mat)
>      nom = 1.
>      if n == 1: return mat[0][0]
>      lastr = mat.pop()
>      jx=-1
>      for j in xrange(n):
>         if lastr[j]:
>             jx=j
>             break
>      if jx==-1: return 0.
>      result = lastr[jx]
>      assert(result<>0.)
>      # Make column jx zero by subtracting a multiple of the last row.
>      for i in xrange(n-1):
>          pivot = mat[i][jx]
>          if 0. == pivot: continue
>          assert(result<>0.)
>          nom *= result   # Compenstate for multiplying a row.
>          for j in xrange(n):
>              mat[i][j] *= result
>          for j in xrange(n):
>              mat[i][j] -= pivot*lastr[j]
>      # Remove colunm jx
>      for i in xrange(n-1):
>         x= mat[i].pop(jx)
>         assert( x==0 )
>
>      if (n-1+jx)%2<>0: result = -result
>      det = determinant( mat )
>      assert(nom<>0.)
>      return result*det/nom
>
> /-----------------------------------------
>
> Now on some matrices the assert triggers, meaning that nom is zero.
> How can that ever happen? mon start out as 1. and gets multiplied
> with a number that is asserted to be not zero.

Easily due to *underflow* precision trouble.  Your "result" may never be 
zero, but it can be very small.  Take the product of many of such tiny 
values, and the result can be less then the smallest value representable 
by a float, at which point it becomes zero.

To see this clearly, try this Python code:
 >>> a = 1.0
 >>> while a > 0:
...   a = a*1.0e-50
...   print(a)
...
1e-50
1e-100
1e-150
1e-200
1e-250
1e-300
0.0

Gary Herron





>
> Any hints appreciated.
>
> Groetjes Albert

Back to comp.lang.python | Previous | Next — Previous in thread | Next in thread | Find similar | Unroll thread

Thread

How can this assert() ever trigger? albert@spenarnc.xs4all.nl (Albert van der Horst) - 2014-05-10 15:24 +0000
  Re: How can this assert() ever trigger? Peter Otten <__peter__@web.de> - 2014-05-10 17:50 +0200
  Re: How can this assert() ever trigger? Ethan Furman <ethan@stoneleaf.us> - 2014-05-10 08:35 -0700
  Re: How can this assert() ever trigger? Gary Herron <gary.herron@islandtraining.com> - 2014-05-10 08:48 -0700
  Re: How can this assert() ever trigger? Alain Ketterlin <alain@dpt-info.u-strasbg.fr> - 2014-05-10 17:56 +0200
    Re: How can this assert() ever trigger? albert@spenarnc.xs4all.nl (Albert van der Horst) - 2014-05-10 16:39 +0000
  Re: How can this assert() ever trigger? Joseph Martinot-Lagarde <joseph.martinot-lagarde@m4x.org> - 2014-05-12 19:05 +0200
    Re: How can this assert() ever trigger? albert@spenarnc.xs4all.nl (Albert van der Horst) - 2014-05-13 09:56 +0000
      Re: How can this assert() ever trigger? Joseph Martinot-Lagarde <joseph.martinot-lagarde@m4x.org> - 2014-05-14 01:14 +0200

csiph-web