\\% Kronecker product of matrices.
\\ Author: Joerg Arndt
\\ online at http://www.jjj.de/pari/
\\ version: 2008-July-25 (02:18)

matkronecker(A, B)=
{ /* Kronecker product  C = A (x) B */
    local(nra, nca, nrb, ncb, nr, nc);
    local(t, C);
    local(ro, co, r, c);
    local(ea, eb, p);
    t = matsize(A); nra=t[1]; nca=t[2];
    t = matsize(B); nrb=t[1]; ncb=t[2];
    nr = nra*nrb;  nc = nca*ncb;
    C=matrix(nr,nc);
    for ( ra=0, nra-1,
        ro = ra * nrb;
        for ( ca=0, nca-1,
            co = ca * ncb;
            ea = A[ra+1,ca+1];
            for ( rb=0, nrb-1,
                r = ro + rb;
                for ( cb=0, ncb-1,
                    eb = B[rb+1,cb+1];
                    c = co + cb;
                    p = ea * eb;
                    C[r+1, c+1] = p;  
                );
            );
        );
    );
    return( C );
} /* ----- */

matkronpow(A, n)=
{ /* Kronecker product  C = A (x) A (x) ... (x) A */
\\ Example: matkronpow([1,-1;1,1],n) gives the (2^n) X (2^n) Hadamard matrix
    local(P);
    P = A;
    for (k=2, n, P=matkronecker(A, P) );
    return( P );
} /* ----- */

matkrondiag(A, n)=
{\\ Example: matkrondiag(A,3) gives [A,0,0; 0,A,0; 0,0,A]
    return( matkronecker(matid(n), A) );
} /* ----- */


\\ ==== end of file ====