Thread Subject: create large sparse FEM matrix --> out of memory

Hi,

I'm writing a FEM application and need to create large sparse matrices.

Unfortunately I run out of memory pretty soon, so I need to know how to write a more memory efficient code in MATLAB.

Here's what I do at the moment to create the stiffness matrix:

I know the number of elements and the size of a single element matrix. I use this information to preallocate the memory for three vectors, storing
*row index
*column index
*value

Then I calculate the element matrices and start filling up those three vectors with the appropriate information and finally K = sparse(row,col,val) for assembling the global matrix.

Now I've got a system with 25 000 elements resulting in a length of about 80million entries for each vector. With 1300MB of free memory MATLAB can't even create the vectors and stops with a "out of memory" error.

Is there a better way both spped and memor-wise?
The one I use is as suggested for a fast creation of the sparse matrix and I also don't think it's to bad memory-wise as the siye scales linear with the number of elements.

Any help is appreciated! Thanks a lot!

"A B" <gitsnedbutzi@hotmail.com> wrote in message <h5pbe7$nc2$1@fred.mathworks.com>...
> Hi,
>
> I'm writing a FEM application and need to create large sparse matrices.
>
> Unfortunately I run out of memory pretty soon, so I need to know how to write a more memory efficient code in MATLAB.
>
> Here's what I do at the moment to create the stiffness matrix:
>
> I know the number of elements and the size of a single element matrix. I use this information to preallocate the memory for three vectors, storing
> *row index
> *column index
> *value
>
> Then I calculate the element matrices and start filling up those three vectors with the appropriate information and finally K = sparse(row,col,val) for assembling the global matrix.
>
> Now I've got a system with 25 000 elements resulting in a length of about 80million entries for each vector. With 1300MB of free memory MATLAB can't even create the vectors and stops with a "out of memory" error.

There is something odd about the description. With a system of 25000 elements and 80e6 non-null entries, that means each of the element has in average 80e6/25e3 = 3200. That doesn't sound right. The average ratio should be about 10 depends on the element type you use, but in no way it could reaches 3200.

To store " three double vectors of 80e9, you need 80e6*3*8 bytes = 1.9 Gbytes, and the resulting sparse matrix need about 80e6*2*8 = 1.3 Gbytes (i think). The total required memory is 3.2 Gbytes. I'm not surprised it does not fit into you 1.3 Gbytes memory.

The 80e6 non-zeros entries does not sound right to me.

Bruno

> There is something odd about the description. With a system of 25000 elements and 80e6 non-null entries, that means each of the element has in average 80e6/25e3 = 3200. That doesn't sound right. The average ratio should be about 10 depends on the element type you use, but in no way it could reaches 3200.
>
> To store " three double vectors of 80e9, you need 80e6*3*8 bytes = 1.9 Gbytes, and the resulting sparse matrix need about 80e6*2*8 = 1.3 Gbytes (i think). The total required memory is 3.2 Gbytes. I'm not surprised it does not fit into you 1.3 Gbytes memory.
>
> The 80e6 non-zeros entries does not sound right to me.
>
> Bruno

The element type in the example is a 20node hexahedron, with 3degrees of freedom (x,y,z) and since a FE matrix is always square it has 20*3x20*3 =3600 elements.

as far as i know, some commercial FE packages write to the disk to temporarly store data and also the price of the license depends on how many nodes it can handle. however i dont know how to write such a code that transfers data to the hard disk and i also think it's very hard to make it efficient. however, i'm not sure wether that's really the way it works. just heard it somewhere...

"A B" <gitsnedbutzi@hotmail.com> wrote in message <h5pm76$ou1$1@fred.mathworks.com>...
> > There is something odd about the description. With a system of 25000 elements and 80e6 non-null entries, that means each of the element has in average 80e6/25e3 = 3200. That doesn't sound right. The average ratio should be about 10 depends on the element type you use, but in no way it could reaches 3200.
> >
> > To store " three double vectors of 80e9, you need 80e6*3*8 bytes = 1.9 Gbytes, and the resulting sparse matrix need about 80e6*2*8 = 1.3 Gbytes (i think). The total required memory is 3.2 Gbytes. I'm not surprised it does not fit into you 1.3 Gbytes memory.
> >
> > The 80e6 non-zeros entries does not sound right to me.
> >
> > Bruno
>
> The element type in the example is a 20node hexahedron, with 3degrees of freedom (x,y,z) and since a FE matrix is always square it has 20*3x20*3 =3600 elements.

I see, but many of them must have the same indices. You don't need to list them all and build the stiffness matrix. I bet the final matrix has much less than 80 millions non-zero elements.

You might be able to allocate the matrix and fill it progressively (or by chunk). Might be loop differently when assembling would help.

>
> as far as i know, some commercial FE packages write to the disk to temporarly store data and also the price of the license depends on how many nodes it can handle. however i dont know how to write such a code that transfers data to the hard disk and i also think it's very hard to make it efficient. however, i'm not sure wether that's really the way it works. just heard it somewhere...

Yes, this is absolutely needed 10 years ago where RAM is limited. For modern computer running on 64-bits OS, I believe writing on disk is no-longer required.

I bet your matrix can be stored entirely in the memory.

Take a look at this FEX, it might help to assemble matrix by chunk. http://www.mathworks.com/matlabcentral/fileexchange/23488

Bruno

> I see, but many of them must have the same indices. You don't need to list them all and build the stiffness matrix. I bet the final matrix has much less than 80 millions non-zero elements.
>
> You might be able to allocate the matrix and fill it progressively (or by chunk). Might be loop differently when assembling would help.

it's true that many nodes will share the same index and therefore the final matrix has much less than 80million non-zero elements. however then the assembly and preallocation gets more difficult...

i guess, i'll have to figure out another way to know the final number of non-zero elements and then write the assemble the matrix.

> > as far as i know, some commercial FE packages write to the disk to temporarly store data and also the price of the license depends on how many nodes it can handle. however i dont know how to write such a code that transfers data to the hard disk and i also think it's very hard to make it efficient. however, i'm not sure wether that's really the way it works. just heard it somewhere...
>
> Yes, this is absolutely needed 10 years ago where RAM is limited. For modern computer running on 64-bits OS, I believe writing on disk is no-longer required.
>
> I bet your matrix can be stored entirely in the memory.

i want it to also work on "older" 32bit systems. also, i'm a little confused why commercial packages can handle 1 million nodes, while i can only handle ~100 000nodes. i think they still must have some tricks.

> Take a look at this FEX, it might help to assemble matrix by chunk. http://www.mathworks.com/matlabcentral/fileexchange/23488

thanks for the link, i'll have a look at it and hope that it can solve my problems

btw. i'm also thinking of switching to a c++ spare matrix library (preferably mtl4). have you got any idea on how fast matlab sparse is on operations as additions and multiplication compared to one of the common c++ libraries?

A B

I have two hints.

The first is to search for an author on the file exchange named Tim Davis. He does a lot of work with sparse matrices, and has published mex functions in C (and I believe c++ in places?) which make MATLAB about as fast as any other system for handling sparse matrices. His Sparse2 command might be very useful is speed of assembly is your issue.

The second is to have a look at the way in which MATLAB stores sparse matrices (comb through the sparse documentation for 'ir', 'jc', 'pr' vectors). You don't need to store both row and column indices. For most sparsity patterns (i.e. anything where you have on average more than one element per column), it is more memory-efficient to store just the row index. Then, store a vector of length ncolumns+1, containing the cumulative number of elements in each column.

Storing like that could save you almost 1/3 of the memory space you're currently allocating.

Kind regards, and good luck!

Tom Clark




"A B" <gitsnedbutzi@hotmail.com> wrote in message <h5pp33$67h$1@fred.mathworks.com>...
> btw. i'm also thinking of switching to a c++ spare matrix library (preferably mtl4). have you got any idea on how fast matlab sparse is on operations as additions and multiplication compared to one of the common c++ libraries?

Hi,

currently speed of the assembly isn't my primary issue, but I'll look into it to increase the speed.

I know that I can store the sparse matrix i another format like compressed column. However interfacing it with MATLAB is the problem, because as far as I have seen it in the documentation sparse only works with row and column indices. Internally it then is treated in a different storage format to use less memory.

One of the most simple solution to save some memory is to make use of the symmetry of the matrices. That can save me almost one half of the required memory, but still I think that more can be achieved.

"Thomas Clark" <t.clark@remove.spamcantab.net> wrote in message <h5q7u0$60t$1@fred.mathworks.com>...
> A B
>
> I have two hints.
>
> The first is to search for an author on the file exchange named Tim Davis. He does a lot of work with sparse matrices, and has published mex functions in C (and I believe c++ in places?) which make MATLAB about as fast as any other system for handling sparse matrices. His Sparse2 command might be very useful is speed of assembly is your issue.
>
> The second is to have a look at the way in which MATLAB stores sparse matrices (comb through the sparse documentation for 'ir', 'jc', 'pr' vectors). You don't need to store both row and column indices. For most sparsity patterns (i.e. anything where you have on average more than one element per column), it is more memory-efficient to store just the row index. Then, store a vector of length ncolumns+1, containing the cumulative number of elements in each column.
>
> Storing like that could save you almost 1/3 of the memory space you're currently allocating.
>
> Kind regards, and good luck!
>
> Tom Clark