#include <iostream>
#include <cmath>
#define PI 3.1415926535897932
using namespace std;

double arsinh(double x)
{
    return log(x+sqrt(x*x+1));
}

double arcosh(double x)
{
    return log(x+sqrt(x*x-1));
}

double cubrt(double x)
{
    if (x >= 0)
        return pow(x, 1.0/3);
    else
        return pow(-x, 1.0/3);
}

int main(){
    system("chcp 1250");
    cout << "Rozwiązywanie równania trzeciego stopnia ax^3+bx^2+cx+d = 0\n";
    cout << "Podaj współczynniki równania oddzielone odstępami: ";
    double a, b, c, d;
    cin >> a >> b >> c >> d;
    if (a != 0) {
        double p = (3*a*c-b*b)/(9*a*a);
        double q = (b*b*b)/(27*a*a*a)-(b*c)/(6*a*a)+d/(2*a);
        double D = q*q+p*p*p;
        double r, fi, x1, x2, x3, y1 ,y2, y3;
        if (q >= 0)
            r = sqrt(fabs(p));
        else
            r = -sqrt(fabs(p));
        if (fabs(D) < 1e-15) {
            if (fabs(p) < 1e-15) { /* q = 0 */
                /* pierwiastek potrójny, x = 0 */
                y1 = 0;
                x1 = y1-b/(3*a);
                cout << "Równanie ma pierwiastek potrójny, x = " << x1 << endl;
            } else { /* p < 0 – inaczej być nie może */
                /* fi = arccos(q/...); itd. */
                fi = acos(q/(r*r*r));
                y1 = -2*r*cos(fi/3);
                x1 = y1-b/(3*a);
                cout << "x = " << x1 << endl;
                y2 = 2*r*cos(PI/3-fi/3);
                x2 = y2-b/(3*a);
                cout << "x = " << x2 << endl;
                y3 = 2*r*cos(PI/3+fi/3);
                x3 = y3-b/(3*a);
                cout << "x = " << x3 << endl;
            }
        } else if (D > 0)
            if (p > 0) {
                /* fi = arsinh(q/...); itd. */
                fi = arsinh(q/(r*r*r));
                y1 = -2*r*sinh(fi/3);
                x1 = y1-b/(3*a);
                cout << "x = " << x1 << endl;
                cout << "x = " << r*sinh(fi/3)-b/(3*a) << " + " << r*cosh(fi/3)*sqrt(3)<< "i\n";
                cout << "x = " << r*sinh(fi/3)-b/(3*a) << " - " << r*cosh(fi/3)*sqrt(3)<< "i\n";
            } else if (p < 0) {
                /* fi = arcosh(q/…); itd.*/
                fi = arcosh(q/(r*r*r));
                y1 = -2*r*cosh(fi/3);
                x1 = y1-b/(3*a);
                cout << "x = " << x1 << endl;
                cout << "x = " << r*cosh(fi/3)-b/(3*a) << " + " << r*sinh(fi/3)*sqrt(3)<< "i\n";
                cout << "x = " << r*cosh(fi/3)-b/(3*a) << " + " << r*sinh(fi/3)*sqrt(3)<< "i\n";
            } else { /* pozostał wariant: p = 0 */
                /* pierwiastek rzeczywisty i dwa zespolone*/
                if (q > 0)
                    y1 = -cubrt(2*q);
                else
                    y1 = cubrt(-2*q);
                x1 = y1-b/(3*a);
                cout << "x = " << x1 << endl;
                cout << "x = " << -y1/2-b/(3*a) << " + " << y1*sqrt(3)/2 << "i\n";
                cout << "x = " << -y1/2-b/(3*a) << " + " << y1*sqrt(3)/2 << "i\n";
            } else { /* przypadek D <= 0 */
                /* fi = arccos(q/...); itd. */
                fi = acos(q/(r*r*r));
                y1 = -2*r*cos(fi/3);
                x1 = y1-b/(3*a);
                cout << "x = " << x1 << endl;
                y2 = 2*r*cos(PI/3-fi/3);
                x2 = y2-b/(3*a);
                cout << "x = " << x2 << endl;
                y3 = 2*r*cos(PI/3+fi/3);
                x3 = y3-b/(3*a);
                cout << "x = " << x3 << endl;
            }
    } else
        cout << "a = 0, równanie nie jest równaniem trzeciego stopnia.\n";
    system("pause");
    return 0;
}