How to make matrix/piddle/array with PDL or perl.

[Feynman]

There are quite a number of algorithms that use banded matricies and I would like to implement a few. I am using PDL because of its speed with large quantities of data.

Is there any particularly efficient way of initializing a banded matrix using PDL? It would be very useful to be able to define the entries of the matrix in terms of off-diagonals.

Eventually, I intend to implement my program with arbitrary dimensions of data, so I do not suppose anyone could help me write a program for creating a banded piddle with arbitrary (but symmetric) dimensionality?

For a banded matrix represented as matrix(i1,i2), it would be useful to specify the matrix elements (i1,i2) such that i1+n=i2. For the N dimensional generalization, one would specify the elements (i1,...iN) such that i1+n1=i2+n2=...=iN.

I am still quite new to PDL so I am not sure how best to go about doing this, so any help would be greatly appreciated.

[BrowserUk]

Neither Perl, nor (as far as I'm aware), PDL provides a mechanism for creating banded (sparse) matrices. You could use normal rectangular matrix, and just ignore the zeros, but for large dimensions, that would waste a lot of space.

However, the "normal" way of dealing with them, is to map the active elements to an appropriately sized (smaller) rectangular matrix. The mapping is very simple & regular. A simple sub that maps the original coordinates to the reduced matrix should be pretty simple to write.

[Anonymous Monk]

Thank you.

Perhaps I should post this as a new question, but I was going to write a spline interpolation code so I could modify it for other forms of interpolation. Hence the reason I wanted to make the band matrix so I could invert it. How might I do that without the matrix itself?

[BrowserUk]

Sorry. Programmer here, not a mathematician.

You're going to have to couch that question in English for me to be able to understand it :)

I'd suggest that you supply a (say) 5x6 3-band matrix as input; and the required output of the inversion.

(Don't forget the <code></code> tags.)

I suspect, that once you have a subroutine that maps the sparse matrix coordinates to the compact rectangular matrix coords, that applying the inversion algorithm in terms of the pre-map coordinates will effect the require transformation upon the compact representation. But I'd need to see it spelt out in english to be sure.

[etj]

For common algorithms, they are likely to be in one of LAPACK or GSL. PDL has bindings to both LAPACK (PDL::LinearAlgebra) and GSL (PDL::GSL).