quaternions faster?

//quat rotate

// 16 mul 12 add (done twice)
// total is 32 mul 24 add
quaternion quaternion::operator*(const quaternion &q) const
{
    quaternion result;

    //[s1,v1][s2,v2] = [s1s2 - v1 dot v2, s1v2 + s2v1 + v1 cross v2]
    result.s = (s * q.s) - x * q.x - y * q.y - z * q.z; // 3 add 4 mul
    result.x = s * q.x + q.s * x + y * q.z - z * q.y;   // 3 add 4 mul
    result.y = s * q.y + q.s * y + z * q.x - x * q.z;   // 3 add 4 mul
    result.z = s * q.z + q.s * z + x * q.y - y * q.x;   // 3 add 4 mul
    return result;
}


//Matrix rotate

// sincos 24 mul 10 add
void rotate_vector(float rad, vec3 &vec, const vec3 &axis)
{
    matrix3 m;
    float   sinVal, cosVal, minusVal;

    sinVal = (float)sinf(rad); // sincos
    cosVal = (float)cosf(rad);
    minusVal = 1.0f - cosVal; // 1 add

    m.m[0] = cosVal + minusVal * axis.x * axis.x;       // 2 mul 1 add
    m.m[1] = minusVal * axis.x * axis.y - sinVal * axis.z;  // 3 mul 1 add
    m.m[2] = minusVal * axis.x * axis.z + sinVal * axis.y;  // 3 mul 1 add

    m.m[3] = minusVal * axis.y * axis.x + sinVal * axis.z;  // 3 mul 1 add
    m.m[4] = cosVal + minusVal * axis.y * axis.y;       // 2 mul 1 add
    m.m[5] = minusVal * axis.y * axis.z - sinVal * axis.x;  // 3 mul 1 add

    m.m[6] = minusVal * axis.z * axis.x - sinVal * axis.y;  // 3 mul 1 add
    m.m[7] = minusVal * axis.z * axis.y + sinVal * axis.x;  // 3 mul 1 add
    m.m[8] = cosVal + minusVal * axis.z * axis.z;       // 2 mul 1 add

    vec = m * vec;      // matrix multiply by vec3
}


// 9 mul 6 add total
vec3 matrix3::operator*(const vec3 &vec) const
{
    vec3 result;

    result.x = vec.x * m[0] + vec.y * m[3] + vec.z * m[6]; // 3 mul 2 add
    result.y = vec.x * m[1] + vec.y * m[4] + vec.z * m[7]; // 3 mul 2 add
    result.z = vec.x * m[2] + vec.y * m[5] + vec.z * m[8]; // 3 mul 2 add
    return result;
}



vec3 quaternion::rotate(float delta, const vec3 &axis, const vec3 &vector)
{
    quaternion  v(0.0f, vector);
    quaternion  result_quaternion;
    vec3        result_vector;
    float       cosval = (float)fcos(delta / 2);
    float       sinval = (float)fsin(delta / 2);

    s = cosval;
    x = sinval * axis.x;
    y = sinval * axis.y;
    z = sinval * axis.z;

    quaternion qinv(s, -x, -y, -z);
    result_quaternion = (*this * v * qinv);
    result_vector.x = result_quaternion.x;
    result_vector.y = result_quaternion.y;
    result_vector.z = result_quaternion.z;
    return result_vector;
}