So you want to compute an outer product? For two 1-dimensional terms (i.e. vectors) you have:
[a1...an]' * [b1...bn] = [a1*b1 ... a1*bn; ...; an*b1 ... an*bn]
which can be written "a_i * b_j" using physicists' tensor notation.
For two 2-dimensional terms (i.e. matrices) you have a tensor product like "a_ij * b_kl", and you can generalize this to any rank (and to other operations besides "*"). I think PDL
can compute outer products for low ranks, but you'll probably need to create an iterator or some such to avoid running out of memory with higher ranks.