/******************************************************************************
 * This program computes the roots of cubic equations
 * of the form Ax^3 + Bx^2 + Cx + D = 0
 * 
 * Copyright © 2020 Richard Lesh.  All rights reserved.
 *****************************************************************************/

#![allow(dead_code)]
#![allow(non_snake_case)]
#![allow(non_upper_case_globals)]
use std::f64::consts::{PI};

fn main() {
    let A:f64 = args[1].trim().parse::<f64>().unwrap_or(0);
    let B:f64 = args[2].trim().parse::<f64>().unwrap_or(0);
    let C:f64 = args[3].trim().parse::<f64>().unwrap_or(0);
    let D:f64 = args[4].trim().parse::<f64>().unwrap_or(0);

/******************************************************************************
 * General cubic can be changed to a depressed cubic 
 * t^3 + pt + q = 0
 * by substituting
 * x = t - B / 3A
 *****************************************************************************/

    let mut p:f64 = (3.f64 * A * C - B * B) / (3.f64 * A * A);
    let mut q:f64 = (2.f64 * B * B * B - 9.f64 * A * B * C + 27.f64 * A * A * D) / (27.f64 * A * A * A);

/******************************************************************************
 * Roots of the depressed cubic are...
 * tk = a cos(b - 2πk/3) for k = 0, 1, 2
 * with a = 2 sqrt(-p/3)
 * b = (1/3) arccos(3q/2p * sqrt(-3/p))
 *****************************************************************************/
    let mut a:f64 = 2.f64 * -p / 3.f64.sqrt();
    let mut b:f64 = (1.f64 / 3.f64) * acos(3.f64 * q / 2.f64 / p * -3.f64 / p.sqrt());

    let mut t1:f64 = a * cos(b);
    let mut t2:f64 = a * cos(b - 2.f64 * PI / 3.f64);
    let mut t3:f64 = a * cos(b - 4.f64 * PI / 3.f64);
    let mut root1:f64 = t1 - B / (3.f64 * A);
    let mut root2:f64 = t2 - B / (3.f64 * A);
    let mut root3:f64 = t3 - B / (3.f64 * A);
    println!("root #1: {}", root1);
    println!("root #2: {}", root2);
    println!("root #3: {}", root3);
}