It's a long time since I worked on this kind of problem, so please bear with me.

I have an approximate inverse matrix that I'm using as a preconditioner to solve the conjugate gradient method. Unfortunately the matrix is much more dense than the storage I have available and reducing it to a much sparser matrix results in a far less optimal solution.

Here is a subset of one of the matrices being produced to get an idea of the arrangement and layout of the data:

10.880     1.917     1.979     0.669     0.311     0.000     0.000     0.000
 1.917    10.871     0.651     1.924     0.000     0.000     0.000     0.000     
 1.979     0.651    11.596     2.159     2.038     0.330     0.000     0.000    
 0.669     1.924     2.159    11.253     0.350     0.000     0.000     0.000   
 0.311     0.000     2.038     0.350    11.221     1.972     0.320     0.311   
 0.000     0.000     0.330     0.000     1.972    11.171     1.862     1.810   
 0.000     0.000     0.000     0.000     0.320     1.862    10.550     0.310   
 0.000     0.000     0.000     0.000     0.311     1.810     0.310    10.550  
 0.322     0.059     2.062     0.694     0.673     0.057     0.000     0.000   
 0.112     0.331     0.736     1.990     0.123     0.011     0.000     0.000   
 0.000     0.000     0.696     0.063     1.934     0.321     0.000     0.000  
 0.000     0.000     0.239     0.342     0.343     0.000     0.000     0.000  
 0.000     0.000     0.000     0.000     0.000     0.000     0.000     0.000   

Now it looks to me that this kind of matrix could be fairly well approximated in fourier space, and by using a wavelet transform I could store a lot more of the data without resorting to a sparse matrix format where I cull the data that doesn't exceed a certain value.

Does that sound reasonable looking at the example dataset? How would I go about doing this, what existing techniques are involved and what resources should I start looking at to help me become familiar with this area?

I'm interested in using the least amount of storage space to store the best approximation to these matrices.

  • $\begingroup$ The (square) part of the matrix you have shown appears to be perfectly symmetric. Is the entire matrix symmetric, or only approximately so? $\endgroup$ – Steve Huntsman Feb 17 '10 at 17:49
  • $\begingroup$ If it was symmetric and small enough to handle, I guess you'd want to diagonalize and throw away the small eigenvalues, then revert back to the original basis and sparsify. But if you can't do something like this b/c of storage, how will FFTing be possible? $\endgroup$ – Steve Huntsman Feb 17 '10 at 17:54
  • $\begingroup$ diagonalizing the matrix and removing small eigenvalues looks feasible, I don't mind using maybe 4 - 5 times the storage space while forming the preconditioner, but need to compress it down when it comes to applying it in the conjugate gradient (and yes it's symmetric) $\endgroup$ – user2731 Feb 19 '10 at 10:38

(Sorry, I can't post this as a comment as my reputation isn't high enough)

@Steve, he's using conjugate gradient, so the matrix had better be symmetric :)

Can you explain more about the structure of the matrix? It's hard to get a feel for the sparsity from that small slice. From the looks of it there, you might try a tridiagonal or pentadiagonal solver as the preconditioner, but there might be some more interesting bits far off the diagonal that we can't see.

Someone else pointed out (in a deleted comment?) that your matrix looks almost-Toeplitz; you might be able to get away with treating it as a circulant matrix, which you can then solve with an FFT. There's some literature on circulant preconditioners out there, but it's been a while since I read anything on them.

Also, it sounds like you're explicitly forming your approximate inverse matrix; you should really only care about solving M*y = z (where M is the preconditioner), right?

Hope some of that helps...

  • $\begingroup$ i don't think a circulant matrix will work here, it really doesn't match that pattern. I am explicitly forming the approximate inverse, so it's really just the preconditioner M that I need to store and ideally compress to get a little more detail. $\endgroup$ – user2731 Feb 19 '10 at 10:45

Wavelets could work, but they're not the only option. Have you had a go with bandwidth reduction reordering? (Or, more subtle - reordering to reduce "fill in" of the Cholesky factor?)

Have you looked at the inverse? Have you looked for low grade (aka semi-separable) structure?

Both of these could lead a bit far afield, but there are very powerful techniques of model reduction which can be applied; the way the inverse elements decay away from the main diagonal can be diagnostic for some of this.

You might consider some sort of multipole acceleration if your matrix is big enough.

  • $\begingroup$ Not familiar with any of these techniques, but they all look like they're worth investigating, thanks. $\endgroup$ – user2731 Feb 19 '10 at 10:48
  • $\begingroup$ It's OK to show us a bit of the inverse. $\endgroup$ – Andrew Mullhaupt Feb 19 '10 at 20:04

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.