rint 0.1.2

#![allow(dead_code)]

use std::f64::consts::PI;

use crate::Integrand;

/* These are the test functions from table 4.1 of the QUADPACK book */

/* f1(x) = x^alpha * log(1/x) */
/* integ(f1,x,0,1) = 1/(alpha + 1)^2 */
pub(crate) struct Function1 {
    pub(crate) alpha: f64,
}

impl Function1 {
    pub(crate) fn new(alpha: f64) -> Self {
        Self { alpha }
    }
}

impl Integrand for Function1 {
    type Point = f64;
    type Scalar = f64;
    fn evaluate(&self, x: &Self::Point) -> Self::Scalar {
        let alpha = self.alpha;
        x.powf(alpha) * (1.0 / x).ln()
    }
}

/* f2(x) = 4^-alpha / ((x-pi/4)^2 + 16^-alpha) */
/* integ(f2,x,0,1) = arctan((4-pi)4^(alpha-1)) + arctan(pi 4^(alpha-1)) */
pub(crate) struct Function2 {
    pub(crate) alpha: f64,
}

impl Function2 {
    pub(crate) fn new(alpha: f64) -> Self {
        Self { alpha }
    }
}

impl Integrand for Function2 {
    type Point = f64;
    type Scalar = f64;
    fn evaluate(&self, x: &Self::Point) -> Self::Scalar {
        let alpha = self.alpha;
        4.0_f64.powf(-alpha) / ((x - PI / 4.0).powi(2) + 16.0_f64.powf(-alpha))
    }
}

/* f3(x) = cos(2^alpha * sin(x)) */
/* integ(f3,x,0,pi) = pi J_0(2^alpha) */
pub(crate) struct Function3 {
    pub(crate) alpha: f64,
}

impl Function3 {
    pub(crate) fn new(alpha: f64) -> Self {
        Self { alpha }
    }
}

impl Integrand for Function3 {
    type Point = f64;
    type Scalar = f64;
    fn evaluate(&self, x: &Self::Point) -> Self::Scalar {
        let alpha = self.alpha;
        (2.0_f64.powf(alpha) * x.sin()).cos()
    }
}

/* f11(x) = log(1/x)^(alpha - 1) */
/* integ(f11,x,0,1) = Gamma(alpha) */
pub(crate) struct Function11 {
    pub(crate) alpha: f64,
}

impl Function11 {
    pub(crate) fn new(alpha: f64) -> Self {
        Self { alpha }
    }
}

impl Integrand for Function11 {
    type Point = f64;
    type Scalar = f64;
    fn evaluate(&self, x: &Self::Point) -> Self::Scalar {
        let alpha = self.alpha;
        f64::ln(1.0 / x).powf(alpha - 1.0)
    }
}

pub(crate) struct Function15 {
    pub(crate) alpha: i32,
}

impl Function15 {
    pub(crate) fn new(alpha: i32) -> Self {
        Self { alpha }
    }
}

impl Integrand for Function15 {
    type Point = f64;
    type Scalar = f64;
    fn evaluate(&self, x: &Self::Point) -> Self::Scalar {
        let alpha = self.alpha;
        x.powi(2) * f64::exp(-2.0f64.powi(-alpha) * x)
    }
}

/* f16(x) = x^(alpha - 1) / (1 + 10 x)^2*/
pub(crate) struct Function16 {
    pub(crate) alpha: i32,
}

impl Function16 {
    pub(crate) fn new(alpha: i32) -> Self {
        Self { alpha }
    }
}

impl Integrand for Function16 {
    type Point = f64;
    type Scalar = f64;
    fn evaluate(&self, x: &Self::Point) -> Self::Scalar {
        let alpha = self.alpha;
        if x.to_bits() == 0f64.to_bits() && alpha == 1 {
            1.0
        } else if x.to_bits() == 0f64.to_bits() && alpha > 1 {
            0.0
        } else {
            x.powi(alpha - 1) / (1.0 + 10.0 * x).powi(2)
        }
    }
}

/* myfn1(x) = exp(-x - x^2) */
/* integ(myfn1,x,-inf,inf) = sqrt(pi) exp(-1/4) */

pub(crate) struct MyFunciton1;

impl MyFunciton1 {
    pub(crate) fn new() -> Self {
        Self
    }
}

impl Integrand for MyFunciton1 {
    type Point = f64;
    type Scalar = f64;
    fn evaluate(&self, x: &Self::Point) -> Self::Scalar {
        f64::exp(-x - x.powi(2))
    }
}

pub(crate) struct MyFunciton2 {
    pub(crate) alpha: f64,
}

impl MyFunciton2 {
    pub(crate) fn new(alpha: f64) -> Self {
        Self { alpha }
    }
}

impl Integrand for MyFunciton2 {
    type Point = f64;
    type Scalar = f64;
    fn evaluate(&self, x: &Self::Point) -> Self::Scalar {
        let alpha = self.alpha;
        f64::exp(alpha * x)
    }
}

/* f455(x) = log(x)/(1+100*x^2) */
/* integ(f455,x,0,inf) = -log(10)/20 */
pub(crate) struct Function455;

impl Integrand for Function455 {
    type Point = f64;
    type Scalar = f64;
    fn evaluate(&self, x: &Self::Point) -> Self::Scalar {
        x.ln() / (1.0 + 100.0 * x.powi(2))
    }
}

impl Function455 {
    pub(crate) fn new() -> Self {
        Self
    }
}

pub(crate) fn test_relative_error(
    calculated: f64,
    target: f64,
    relative_error: f64,
    description: &str,
) -> Result<(), String> {
    if calculated.is_nan() || target.is_nan() {
        if calculated.is_nan() != calculated.is_nan() {
            return Err(format!(
                "Failed test {description}: calculated.is_nan() != target.is_nan()"
            ));
        }
    }

    if calculated.is_infinite() || target.is_infinite() {
        if calculated.is_infinite() != calculated.is_infinite() {
            return Err(format!(
                "Failed test {description}: calculated.is_infinite() != target.is_infinite()"
            ));
        }
    }

    if (target > 0.0 && target < f64::MIN_POSITIVE) || (target < 0.0 && target > -f64::MIN_POSITIVE)
    {
        return Err(format!(
            "Failed test {description}: target value smaller than f64::MIN_POSITIVE"
        ));
    }

    if target != 0.0 {
        if (calculated - target).abs() / target.abs() > relative_error {
            return Err(format!(
                "Failed test {description}: calculated relative error is larger than target"
            ));
        }
    } else {
        if calculated.abs() > relative_error {
            return Err(format!("Failed test {description}: calculated integral value is smaller than calculated relative error"));
        }
    }

    Ok(())
}

//pub(crate) fn test_absolute_error(
//    calculated: f64,
//    target: f64,
//    absolute_error: f64,
//    description: &str,
//) -> Result<(), String> {
//    if calculated.is_nan() || target.is_nan() {
//        if calculated.is_nan() != calculated.is_nan() {
//            return Err(format!(
//                "Failed test {description}: calculated.is_nan() != target.is_nan()"
//            ));
//        }
//    }
//
//    if calculated.is_infinite() || target.is_infinite() {
//        if calculated.is_infinite() != calculated.is_infinite() {
//            return Err(format!(
//                "Failed test {description}: calculated.is_infinite() != target.is_infinite()"
//            ));
//        }
//    }
//
//    if (target > 0.0 && target < f64::MIN_POSITIVE) || (target < 0.0 && target > -f64::MIN_POSITIVE)
//    {
//        return Err(format!(
//            "Failed test {description}: target value smaller than f64::MIN_POSITIVE"
//        ));
//    } else {
//        if (calculated - target).abs() > absolute_error {
//            return Err(format!("Failed test {description}: calculated absolute error is smaller than target relative error"));
//        }
//    }
//    Ok(())
//}

pub(crate) fn test_int(calculated: usize, target: usize, description: &str) -> Result<(), String> {
    if calculated == target {
        Ok(())
    } else {
        Err(format!(
            "Failed test {description}: calculated value: {calculated} target value: {target}"
        ))
    }
}