For an online quantum mechanics course, I had to construct a matrix representation of the tensor product of two Pauli Spin Matrices. I did it as follows:

In[1]:= p1 = PauliMatrix[1]; p2 = PauliMatrix[2];

In[2]:= Join[Join[p1[[1, 1]]*p2, p1[[1, 2]]*p2, 2], 
 Join[p1[[2, 1]]*p2, p1[[2, 2]]*p2, 2]]

Out[2]= {{0, 0, 0, -I}, {0, 0, I, 0}, {0, -I, 0, 0}, {I, 0, 0, 0}}

I got the right answer but don't like the method I used, as it would be difficult to generalize to large matrices. Any suggestions?