Creating overlapped block matrices

I'm doing an FEM assignment using Mathematica.

(EK1 = {{a11, a12}, {a21, a22}}) // MatrixForm
(EK2 = {{b22, b23}, {b32, b33}}) // MatrixForm


I don't know how the best way to create new matrix like this:

(K = {{a11, a12, 0}, {a21, a22 + b22, b23}, {0, b32, b33}}) // MatrixForm


MATLAB does it quite easily but how about Mathematica?

-
Is there a pattern to this? With the b23 where it is and not a21 on the left side of that row I can't think of one. –  Pickett Aug 19 '13 at 8:52
Anon, there are stiffness matrices. There is a coupling between elements and where they join, or overlap, the stiffness is added. But this all depends on the elements and the structure and can get very complicated. –  Nasser Aug 19 '13 at 10:00
Viet Tran, you haven't responded to Anon's comment or any of the answerers. It seems like the 0 in position {2,1} was a typo (as everyone seems to agree) and I've edited the question to change it to a21, so that the answers make sense. If you disagree, please leave a comment, else I'll assume that you are in agreement with this change. –  rm -rf Dec 31 '13 at 1:18

I wrote a small function once to build block-wise diagonal matrices. The zero matrix is just 0 and identity is 1. You could use it in this case like this:

diagonalize[list_] := ArrayFlatten@(DiagonalMatrix[Array[x, Length@list]] /.
Table[x[i] -> list[[i]], {i, Length@list}]);
diagonalize[{EK1, 0}] + diagonalize[{0, EK2}] // MatrixForm


Oh, this is assuming that the second zero in your K is actually a21. Otherwise I don't see a pattern here.

-

You can use ArrayPad:

ArrayPad[EK2, {1, 0}] + ArrayPad[EK1, {0, 1}]


and general approach (Edit: overlap specification added):

n = {3, 3};
arrays = Array[#, n] & /@ {a, b, c, d, e}; (*arrays to work with*)

app[a1_, a2_, overlap_: 1] := With[{dim = Dimensions@a1},
ArrayPad[a1, {0, n[[1]] - overlap}] +
ArrayPad[a2, Transpose@{dim - overlap, {0, 0}}]];

Fold[app, First@arrays, Rest@arrays] // MatrixForm
Fold[app[##, 2] &, First@arrays, Rest@arrays] // MatrixForm
Fold[app[##, 3] &, First@arrays, Rest@arrays] // MatrixForm


It is not so much general, only for constant n$x$n set but I don't know what should be an expected result for set of arrays of different dimensions.

-
This has a21 where there is a zero in the OP's example. I think that might be a typo on the OP's part though. –  Pickett Aug 19 '13 at 8:51
@Anon ups, you are right. But I think it is a mistake in question :p –  Kuba Aug 19 '13 at 8:51

(* I am assuming this is meant to be a global stiffness matrix construction question and OP made type for element (2,1) *)

There are many ways to build global stiffness matrix. For your case, here is a quick hack. But if you simply google the topic of building global stiffness matrix from elements stiffness matrices, you'll find many many methods on it in finite elements books.

Notice, I gave element names inside the local stiffness as b11,b12, etc.... different than what you had to be consistent with the first matrix.

EK1 = {{a11, a12}, {a21, a22}}; (*local stiffness matrix*)
EK2 = {{b11, b12}, {b21, b22}}; (*local stiffness matrix*)
n = 2; (*size of local stiffness matrix *)
locals = {EK1, EK2}; (*local stiffness matrices*)
m = Length[locals]; (*number of elements*)
r = (n m) - (m - 1); (*size of global stiffness matrix*)
K1 = Table[0, {r}, {r}];
z = 1;

Do[ K1[[z ;; z + n - 1, z ;; z + n - 1]] = locals[[i]];
If[i > 1, K1[[z, z]] += locals[[i - 1]][[-1, -1]]];
z = z + n - 1,
{i, 1, m} ];


For example, for 5 elements, you'll get this:

-
I posted similar answer but @Anon has pointed out that it does not fit the question. Take a look at position {2, 1} of OP's result. –  Kuba Aug 19 '13 at 10:02
@Kuba. I see. Oh well. I just guess the OP is either wrong (made typo), or they should then say what does this pattern represents. I assumed it is a global stiffness matrix. If not, then I will also delete this answer. –  Nasser Aug 19 '13 at 10:04
I think you are right, I think so. But we can't be sure :) –  Kuba Aug 19 '13 at 10:08
@Kuba I also read at the top OP saying FEM assignment using Mathematica. and that is why I thought right away it is stiffness matrix. (FEM== finite elements) –  Nasser Aug 19 '13 at 10:10
Good point. I could have checked this. –  Kuba Aug 19 '13 at 10:13