Putting together a big stiffness matrix from several small ones

Question

This is a trivial problem, but I'm just starting with matlab, and haven't got used to their way of thinking yet (and syntax).

What I'm asking will be obvious to anyone who's ever done anything with FEM or the like.

How do you put together a big stiffness matrix from several small ones. Say, you got for (element 1) a local stiffness matrix 4x4, the same for (element 2) - only different matrix, of course, but still 4x4.

What is the easiest way to do this:

[|--------| 0  0 ]  
[|        | 0  0 ]  
[|     |--|-----|]  
[|-----|--|     |]  
[0  0  |        |]  
[0  0  |--------|]  

(a33+b11, a34+b12,
(a43+b12, a44+b22, ...)

i.e. make a 'big one' ?

Answer 1

I think your question is this:

A = 4x4 B = 4x4

C = final matrix where A and B overlap and should be summed in overlap.

Do this:

C = zeros(6);
C(1:4,1:4) = A;
C(3:6,3:6) = C(3:6,3:6) + B;

Answer 2

Stiffness matrices WANT to be sparse. Desperately so. The point is that these big stiffness matrices will end up with a huge number of zero elements. So you need to build the matrix with that in mind. The nice thing is that sparse matrices are very nicely supported in matlab. But you need to work with them properly.

The idea is to build your matrix up as a list of row and column indices, plus the value to be inserted into the matrix. Build up this list in advance. Only at the very end do you call sparse to actually build the matrix itself. Sparse will automatically sum up elements where there is overlap.

Once that sparse matrix is built, all operations, like matrix multiplies and backslash are fully supported, and can be very fast compared to the same operations on a full matrix. This is especially important when your global stiffness matrix might be 1e5x1e5 or larger.

In some cases for huge matrices, you may need to migrate to an iterative solver when factorizations of your matrix become impossible.

Answer 3

Just for fun, here's a one-line solution using BLKDIAG:

C = blkdiag(A,zeros(2)) + blkdiag(zeros(2),B);