# editing entries of a sparse array + sparse Smith normal form

Why does the following code

boundaries = SparseArray /@ {{{1,4}->0}, {{4,6}->0}, {{6,4}->0}, {{4,1}->0}};
MatrixForm /@ boundaries
boundaries[]
boundaries[[3,6,4]]
boundaries[[3,6,4]] += 100;


return an error

"Part 6 of boundaries[[3,6,4]] does not exist."?

The program is able to read the 6,4-th entry in the third matrix, yet it is unable to change it. How so? But it is able to read and change the 4,4th entry in the third matrix, which boggles me...

Also, is there a package for Mathematica, that can compute the Smith normal form (actually just the diagonals, i.e. invariant factors) of a sparse array? If not, what are other (highly effective) software systems that could do that. I don't think GAP can handle SNF of sparse matrices. Chomp and RedHom only work on simplicial/cubical complexes. I'm not sure about LinBox.

A workaround for the first question

boundaries = MapAt[# + 100 &, boundaries, {3, 6, 4}];

MatrixForm /@ boundaries • Hmm, in my program I construct a list boundaries' of matrices with zeros and then have a large for loop in which I edit the entries of each matrix in the list (I'm computing the homology of a chain complex). So editing entries in each sparse array is vital. Does your solution allow that? Also, why does Mathematica return an error to my code above? It is so strange that editing some entries is allowed, but others are untouchable. – Leo Nov 7 '14 at 12:07
• @Leon, array of sparse arrays is unusual for Mathematica. It doesn't allow to apply +=. However, there are two another ways: one 3-dimensional sparse array or several distinct variables instead of boundaries. – ybeltukov Nov 7 '14 at 12:50
• Hmm, how can I create a 3-dimensional sparse array? My current code: boundaries = Prepend[Table[ SparseArray@{{Length@bases[[k-]],Length@bases[[k]]}->0},{k,2,m}], SparseArray@{{1,Length@basis}->0}]; – Leo Nov 7 '14 at 19:58
• @Leon You can use three indexes like {k,i,j}->0 and then SparseArray@Flatten[rules]. – ybeltukov Nov 7 '14 at 20:02
• Thank you! Do you know perhaps something about SNF of sparse matrices? – Leo Nov 12 '14 at 19:17

Perhaps not ideal using Normal converts to list of lists rather than list of sparsearray objects:

norm = Normal@boundaries;
norm[[3, 6, 4]] += 100;
boundaries = SparseArray /@ norm;
MatrixForm /@ boundaries (* to demonstrate result*)
`
• So normal removes the SparseArray structure, so that I obtain ordinary matrices? My matrices are very large and very sparse, so this is not very useful. – Leo Nov 7 '14 at 12:10