Polynomial multiplication in M2(R)?

73 views Asked by At

I was trying to implement a FFT-based multiplication algorithm in M2(R). Basically an algorithm that gets as an input two polynoms with elements given as matrices, and builds the product polynom. However, even though the algorithm should work, as it looks exactly identical to a version I wrote earlier on regular number, it doesn't. The coefficients are always off by a little.

I have not found any articles on the roots of unity in M2(C), but I have found(on paper) that choosing eps = ((cos(2PI/n) , i sin(2PI/n)) , ( i sin(2PI/n) , cos(2PI/n))), I get a nice cycle.

Is there something wrong in my approach?

Here is the code:

struct FFT {

    PolyC To, Aux[17][2], Res[17][2], AC, BC, ResC, ResD, ArgA, ArgB;

    void fft(PolyC V, var depth, var n, PolyC To, MatC step) {
        if(n == 1) {
            To[0] = V[0];
        } else {

            MatC eps = matCHelper.I2;

            //We "split" the poly in 2
            for(var i=0; i<n; i++)
                Aux[depth+1][i&1][i>>1] = V[i];

            //We recursively apply FFT to the components
            fft(Aux[depth+1][0], depth+1, n/2, Res[depth+1][0], step*step);
            fft(Aux[depth+1][1], depth+1, n/2, Res[depth+1][1], step*step);

            //We compute the result for the n roots
            for(var i=0; i<n/2; i++) {
                To[i] = Res[depth+1][0][i] + eps * Res[depth+1][1][i];
                To[n/2+i] = Res[depth+1][0][i] - eps * Res[depth+1][1][i];
                eps = eps * step;
            }
        }
    }

    void FFTMultiply(Poly Res, Poly A, Poly B, var n1, var n2) {

        var M;
        for(M = 1; M <= 2*n1 || M <= 2*n2; M <<= 1);

        for(var i=0; i<n1; i++) ArgA[i] = A[i];
        for(var i=n1; i<M; i++) ArgA[i] = matCHelper.O2;

        for(var i=0; i<n2; i++) ArgB[i] = B[i];
        for(var i=n2; i<M; i++) ArgB[i] = matCHelper.O2;

        MatC step( Complex(cos(2*PI/M), 0) , Complex(0, sin(2*PI/M)),
                   Complex(0, sin(2*PI/M)) , Complex(cos(2*PI/M), 0) );

        fft(ArgA, 0, M, AC, step);
        fft(ArgB, 0, M, BC, step);

        for(var i=0; i<M; i++) {
            RezC[i] = AC[i] * BC[i];
        }

        step.b = -step.b;
        step.c = -step.c;

        fft(RezC, 0, M, RezD, step);

        for(var i=0; i<M; i++) {
            // Now I divided everything by M and copied every element of ResD to Res modulo some number
        }
    }
};
1

There are 1 answers

0
Lutz Lehmann On BEST ANSWER

You can not expect this method to work if your coefficient matrices do not commute with your matrix step. To get it working correctly, use the diagonal matrix corresponding to multiplication with the scalar exp(i*2*PI/M).