A fast solution to cubic equations with three real roots

fastCubic

http://momentsingraphics.de/CubicRoots.html#_Blinn06a

Returns the three real roots of the cubic polynomial (of form 1 + 2x + 3x^2 + 4x^3)

float3 SolveCubic(float4 Coefficient){
    // Normalize the polynomial
    Coefficient.xyz/=Coefficient.w;
    // Divide middle coefficients by three
    Coefficient.yz/=3.0f;
    // Compute the Hessian and the discrimant
    float3 Delta=float3(
        mad(-Coefficient.z,Coefficient.z,Coefficient.y),
        mad(-Coefficient.y,Coefficient.z,Coefficient.x),
        dot(float2(Coefficient.z,-Coefficient.y),Coefficient.xy)
    );
    float Discriminant=dot(float2(4.0f*Delta.x,-Delta.y),Delta.zy);
    // Compute coefficients of the depressed cubic 
    // (third is zero, fourth is one)
    float2 Depressed=float2(
        mad(-2.0f*Coefficient.z,Delta.x,Delta.y),
        Delta.x
    );
    // Take the cubic root of a normalized complex number
    float Theta=atan2(sqrt(Discriminant),-Depressed.x)/3.0f;
    float2 CubicRoot;
    sincos(Theta,CubicRoot.y,CubicRoot.x);
    // Compute the three roots, scale appropriately and 
    // revert the depression transform
    float3 Root=float3(
        CubicRoot.x,
        dot(float2(-0.5f,-0.5f*sqrt(3.0f)),CubicRoot),
        dot(float2(-0.5f, 0.5f*sqrt(3.0f)),CubicRoot)
    );
    Root=mad(2.0f*sqrt(-Depressed.y),Root,-Coefficient.z);
    return Root;
}