Firstov.com

Вначале немного теории. Как обычно, будем рассматривать уравнение в котором коэффициент при x³ равен нулю:

x4 + px2 + qx + r = 0    (1)

Представим это уравнение в виде произведения двух квадратных:

x4 + px2 + qx + r = ( x2 + ax + b ) ( x2 - ax + c )    (2)

Тогда для коэффициентов a, b, c получим следующую систему уравнений:

b + c - a2 = p
a ( c - b ) = q    (3)
bc = r

Из первых двух уравнений системы (3) получаем:

2c = p + a2 + q/a    (4)
2b = p + a2 - q/a

Затем подставляем выражения для b и c в третье уравнение системы (3) и получаем уравнение шестой степени относительно a:

a6 + 2pa4 + (p2 - 4r)a2 - q2 = 0    (5)

Сделаем подстановку: y = a².

y3 + 2py2 + (p2 - 4r)y - q2 = 0    (6)

Таким образом получается следующий метод решения уравнения четвёртой степени (1).
Вначале решаем уравнение (6). Берём произвольный корень и получаем a. Затем подставляем a в (4) и находим b и c. Остаётся только решить два квадратных уравнения из (2) и получить четыре корня.
Для нахождения действительных корней уравнения четвёртой степени предлагаются следующие две программы:

// x^4 + p * x^2 + q * x + r = 0

int root4s ( double p, double q, double r, double * x )
{
    if ( q == 0 )
    {
        double t[2];
        int n = root2 ( p, r, t );
        if ( n == 0 ) return 0;
        if ( t[0] >= 0 )
        {
            double s = sqrt ( t[0] );
            x[0] = s;
            x[1] = -s;
            n = 2;
        }
        else
        {
            n = 0;
        }
        if ( t[1] >= 0 )
        {
            double s = sqrt ( t[1] );
            x[n] = s;
            x[n+1] = -s;
            n += 2;
        }
        return n;
    }
    int n = root3 ( p+p, p*p - 4.*r, -q*q, x );
    double a = *x;
    if ( n == 3 )
    {
        if ( a < x[1] ) a = x[1];
        if ( a < x[2] ) a = x[2];
    }
    if ( a < 0 ) return 0;
    p += a;
    a = sqrt(a);
    double b = q/a;
    n = root2 (  a, 0.5*(p-b), x   );
    n+= root2 ( -a, 0.5*(p+b), x+n );
    return n;
}

// x^4 + a * x^3 + b * x^2 + c * x + d = 0

int root4 ( double a, double b, double c, double d, double * x )
{
    if ( a == 0 )
    {
        return root4s ( b, c, d, x );
    }
    if ( d == 0 )
    {
        *x = 0;
        return 1 + root3 ( a, b, c, x+1 );
    }
    double e = a / 4;
    double p = -6*e*e + b;
    double q = 8*e*e*e -2*b*e + c;
    double r = -3*e*e*e*e + b*e*e - c*e + d;
    int n = root4s ( p, q, r, x );
    for ( int i = 0; i < n; ++i )
    {
        x[i] -= e;
    }
    return n;
}

Первая из них решает "неполное" уравнение, т.е. у которого коэффициент при x³ равен нулю, а вторая сводит общий случай к предыдущему.