Interpolation between spaces exploiting tensor product structure

The number of floating point operations per cell is currently `O((P0+1)^d *(P1+1)^d) =~ O(P^2*d)` where P0 is the degree of interpolation from and P1 the degree for interpolation to.

If the pair of elements have tensor product factorization, this number can be reduced to ~ `P^d+1` , if only factor matrices are computed.
I'm not sure if basix should return this 1d matrices, or the user should request the interpolation matrix for a pair of factors instead.