Thread Subject: Any way to store and operate on large matrices built from smaller

Hi, I've used matlab for a long time and this is a question that has
me stumped.

I'm doing some Cramer-Rao bounds on jointly Gaussian random vectors
and I can represent the covariance matrixes of the random vectors
quite easily. I can also represent their derivatives easily
numerically, so I have a matrix-arithmetic formula that chugs through
(it includes some inverses, some matrix-multiplies, some transposes,
etc.).

My problem is, beyond toy examples, the different matrixes I need
become extremely large. E.g., when my data vector is 575x1, its
covariance matrix needs to be 575x575. To represent just one 575x575
complex matrix needs 5 GB of memory (575^2 * 16 / 1024 / 1024), or
half as much if I use complex singles, but I don't have that much RAM.

However, that 575x575 matrix is generated by operating on a few 575x1
vectors. That is, I have some ten vectors 575-elements long and they
combine to produce this huge 575x575 matrix.

So my question is: is there any way in matlab to build matrices out of
(simple) operations on vectors?

I think my problem can be solved by answering this meta-question: is
there a way to represent repmat(vec, 1, N) using just "vec" and "N"
and not copy "vec" so many times? And then is there a way to propagate
that through matrix equations that at the end yield scalars, so that
when matlab goes to do matrix multiplication, e.g., it just pulls the
appropriate number out of the original vectors?

I'm going to try and see if Maple does symbolic matrix algebra and see
if I can simplify the Cramer-Rao bound equation, which is in the form
of 1/2 * trace(A / B * C / B) + x' / B * x, where A, B, C are large
square matrices, and x is a column vector. (A, B, and C are results of
matrix-arithmetic with other vectors.)

Thanks in advance!!!

wuzzyview@gmail.com wrote:
...
> ...just one 575x575 complex matrix needs 5 GB of memory ...

I get ~40 _M_B, quite a bit less???

--

Doh, you're right, dpb, I made a mistake in my notes and another in my
arithmetic.

My data vectors can get as long as 38'928 elements, which would make a
single covariance matrix >22 GB (need to divide by 1024 three times,
38928^2 * 16 / 1024 / 1024 / 1024).

I know this is a problem because my workstation with 4GB of RAM dies
very quickly!

On Jun 21, 7:19 pm, dpb <n...@non.net> wrote:
> wuzzyv...@gmail.com wrote:
>
> ...
>
> > ...just one 575x575 complex matrix needs 5 GB of memory ...
>
> I get ~40 _M_B, quite a bit less???
>
> --

I guess this question basically boils down to: how do I solve a system
of equations, Ax=y so large that A can't fit in memory, but where A is
constructed from fairly simple arithmetic expressions on a handful of
vectors? How can I sacrifice CPU speed to conserve memory?

I don't remember now, but when I was working on Finite Element Methods, I
had a similar problem and I remember having found some techniques from a
Numerical Analysis book.

Another alternative would be to expand the equation Ax = y by expanding A
since it is derived from simple operations on vectors as you say. Can you
write out A in terms of those vectors? It may be helpful.

<wuzzyview@gmail.com> wrote in message
news:41d82244-c1f3-4b4f-ab85-6b4f09ee6608@c58g2000hsc.googlegroups.com...
>I guess this question basically boils down to: how do I solve a system
> of equations, Ax=y so large that A can't fit in memory, but where A is
> constructed from fairly simple arithmetic expressions on a handful of
> vectors? How can I sacrifice CPU speed to conserve memory?

<wuzzyview@gmail.com> wrote in message
news:41d82244-c1f3-4b4f-ab85-6b4f09ee6608@c58g2000hsc.googlegroups.com...
>I guess this question basically boils down to: how do I solve a system
> of equations, Ax=y so large that A can't fit in memory, but where A is
> constructed from fairly simple arithmetic expressions on a handful of
> vectors? How can I sacrifice CPU speed to conserve memory?

1) Is A sparse? If so, then try creating it as a sparse matrix.

2) If you know how to compute the matrix-vector product A*x without actually
constructing A, try using one of the iterative solvers like GMRES or PCG:

http://www.mathworks.com/access/helpdesk/help/techdoc/ref/gmres.html

http://www.mathworks.com/access/helpdesk/help/techdoc/ref/pcg.html

3) Can you break apart your large system of equations into smaller
subsystems, where no two of the smaller subsystems have a variable in
common?

4) If your matrix is REALLY big, you could try using Parallel Computing
Toolbox to create your A matrix as a distributed array:

http://www.mathworks.com/access/helpdesk/help/toolbox/distcomp/bqi9fln-1.html

--
Steve Lord
slord@mathworks.com

Thank you Mr Lord! Excellent suggestions!

1_ A is not sparse unfortunately. A is in the form of [A1 A2; conj(A2)
conj(A1)], and where A1=(m * m') .* toeplitz(w) +
eye(length(m))*alpha, A2=(m * transpose(m)) .* toeplitz(w), where m
and w are Nx1 vectors (N very large, 30-60k+), and where w is
sufficiently complicated that the Toeplitz part nor the outer product
give much sparseness.

2_ I will definitely will look into the iterative solvers. I'm fairly
sure it's more feasible to compute A*x without having to store A.

3_ You may be able to tell from 1_ but A is pretty dense, though it
has a lot of structure. Looking at imagesc(abs(A)) for example shows a
fair number of very small values, I may be able to analytically
triangle-block the larger matrix.

In my introductory DSP class, my professor talked about the very
elegant results of Butterworth et al. filter design from the pre-
computer era and hypothesized that if they had Matlab back then, they
probably wouldn't have discovered any of the beautiful math. That
lecture is resonating in me---there's probably a very smart
mathematical way to solve this particular problem and I'm chastising
myself... because just this morning I got an account on machines with
64 GB RAM (http://www.osc.edu/supercomputing/hardware/). I think this
might be the most worthwhile way forward :P thanks!

On Jun 22, 10:00 pm, "Steven Lord" <sl...@mathworks.com> wrote:
> <wuzzyv...@gmail.com> wrote in message
>
> news:41d82244-c1f3-4b4f-ab85-6b4f09ee6608@c58g2000hsc.googlegroups.com...
>
> >I guess this question basically boils down to: how do I solve a system
> > of equations, Ax=y so large that A can't fit in memory, but where A is
> > constructed from fairly simple arithmetic expressions on a handful of
> > vectors? How can I sacrifice CPU speed to conserve memory?
>
> 1) Is A sparse? If so, then try creating it as a sparse matrix.
>
> 2) If you know how to compute the matrix-vector product A*x without actually
> constructing A, try using one of the iterative solvers like GMRES or PCG:
>
> http://www.mathworks.com/access/helpdesk/help/techdoc/ref/gmres.html
>
> http://www.mathworks.com/access/helpdesk/help/techdoc/ref/pcg.html
>
> 3) Can you break apart your large system of equations into smaller
> subsystems, where no two of the smaller subsystems have a variable in
> common?
>
> 4) If your matrix is REALLY big, you could try using Parallel Computing
> Toolbox to create your A matrix as a distributed array:
>
> http://www.mathworks.com/access/helpdesk/help/toolbox/distcomp/bqi9fl...
>
> --
> Steve Lord
> sl...@mathworks.com