SciRuby / nmatrix

Dense and sparse linear algebra library for Ruby via SciRuby
Other
469 stars 133 forks source link

yale matrix's ndnz isn't initialized during the initialization of yale matrix #337

Open spacewander opened 9 years ago

spacewander commented 9 years ago

When I try to implement det_exact, I find a strange behavior:

The ndnz of yale matrix(YALE_STORAGE) will always be 0, regardless of what initial value is given. But, if the yale matrix is cast from other matrix, the ndnz will be a correct value.

I read the code and find out that two special functions will be called to initialize yale matrix:

  1. storage/yale/class.cpp:init
  2. storage/yale/class.cpp:create

None of them initialize the value of ndnz. I am newbie to the source code of nmatrix, have I missed something?

translunar commented 9 years ago

It's been a really long time since I looked at the source code for Yale, unfortunately. You could look at the changelog for that file and try to figure out what's going on with ndnz, or search around in the comments, but off the top of my head I don't know. =(

Is this actually a bug? Is there anything that isn't working as a result? Or is it just a note that ndnz is not used consistently?

spacewander commented 9 years ago

It seems that ndnz is only set to 0 in the initialization of yale matrix from the time the functions were written...

A quick look at: https://github.com/SciRuby/nmatrix/blob/c9ff3d0e68bea9213550bdc984a1f3a368633077/ext/nmatrix/storage/yale.cpp#L1445 And https://github.com/SciRuby/nmatrix/blob/master/ext/nmatrix/storage/yale/class.h#L556

translunar commented 9 years ago

How would you go about fixing it?

spacewander commented 9 years ago

What about calculating the ndnz after a and ija are filled in the init function?

translunar commented 9 years ago

What happens when elements get deleted? What happens when the matrix's non-diagonal default is non-zero?

spacewander commented 9 years ago

Well, I added a function called yale_ndnz just now. It can be checkout from https://github.com/spacewander/nmatrix/tree/ndnz This function works like yale_ija and yale_a, it prints the value of ndnz. With its help, I can test the ndnz of some calculation results. It seems that some calculation don't update ndnz.

For example:

pry(main)> n = NMatrix.new(2, [1,2,3,4], stype: :dense, dtype: :int64).cast(:yale)
pry(main)> n.extend NMatrix::YaleFunctions
pry(main)> n.yale_ndnz
=> 2
pry(main)> a = n.clone
pry(main)> a.yale_ndnz
=> 2
pry(main)> b = n * a
=>
[
  [1,  4]   [9, 16] ]
pry(main)> b.extend NMatrix::YaleFunctions
pry(main)> b.yale_ndnz
=> 0

ndnz should be equal to ija[shape[0]] - ija[0]. If the ndnz is unreliable, we can use ija[shape[0]] - ija[0] instead. In other way, we can also use ija[shape[0]] - ija[0] to update ndnz after a new yale matrix is created or a yale matrix is modified.

translunar commented 9 years ago

It seems like ndnz is defined as being part of the IJA array, isn't it? So have we duplicated this by including it in the Yale struct also?

Argh, frustrating trying to figure this out. It's been so long since I wrote this code.

spacewander commented 9 years ago

Well, ndnz can be calculated out each time we need. Keeping the value of ndnz in Yale struct can save a little calculation, but it will bite us if we forget to update ndnz after changing elements. Considering with the effort to maintain ndnz, keeping ndnz may cost more.

IMHO, ndnz is not so useful, and may be duplicate.

translunar commented 9 years ago

I suspect one reason I didn't remove it before is that it causes the data in yale, list, and dense to be in the same relative memory location.