void mul12_SSE(double x, double y, double* z, double* zz){
        const double factor = 134217729;
        __m128d factors = {factor, factor};
        __m128d xy      = {x, y};
        __m128d hh, tt, pp, qq;

        pp = _mm_mul_pd(xy, factors);
        hh = _mm_sub_pd(xy, pp);
        hh = _mm_add_pd(hh, pp);
        tt = _mm_sub_pd(xy, hh);

        pp = _mm_mul_pd(hh, _mm_shuffle_pd(hh, hh, 1));
        qq = _mm_mul_pd(hh, _mm_shuffle_pd(tt, tt, 1));
        qq = _mm_add_pd(qq, _mm_shuffle_pd(qq, qq, 1));

        double p[2] __attribute__ ((aligned (16)));
        _mm_store_pd(p, pp);
        double q[2] __attribute__ ((aligned (16)));
        _mm_store_pd(q, qq);
        double txty[2] __attribute__ ((aligned (16)));
        _mm_store_pd(txty, _mm_mul_pd(tt, _mm_shuffle_pd(tt, tt, 1)));

        *z = p[1] + q[1];
        *zz = p[1] - *z + q[1] + txty[1];
}