1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105

//! # `ideal_gas`
//!
//! Run this example with:
//! ```sh
//! cargo run --example ideal_gas
//! ```
//!
//! ## Setup
//!
//! Find pressure that satisfies the ideal gas equation
//!
//! $PV = nT$
//!
//! where
//!
//! | Notation | Meaning         | value               |
//! |---------:|:----------------|:--------------------|
//! |      $P$ | pressure        | 1.0 (initial guess) |
//! |      $V$ | volume          | 1.618               |
//! |      $n$ | number of moles | 1.0                 |
//! |      $T$ | temperature     | 300.0               |
//!
//! ## Solution
//!
//! In this example, we use Newton method to solve the equation.
//!

use autodj::*;

fn main() {
    let [pressure, volume, temperature, quantity]: [f64; 4] = [1., 1.618, 300., 1.];

    let calc_residual_dual =
        |x| calc_ideal_gas_generic(x, volume.par(), temperature.par(), quantity.par());

    let residual_generic = pressure.derive(&calc_residual_dual);

    print_state_linearization(residual_generic.val, residual_generic.dual, pressure);

    let pressure_newtoned = newton_iterations(calc_residual_dual, pressure, 1e-3, 10);

    match pressure_newtoned {
        Ok(pressure_refined) => {
            println!("{pressure} refined to {pressure_refined} using Newton method")
        }
        Err(ConvergenceError(err)) => {
            println!("Not converged Newton iterations:");
            match err {
                Some(err) => println!("----with an error of {}", err),
                None => println!("----function was not evaluated"),
            }
        }
    }
}

fn print_state_linearization(value: f64, deriv: f64, origin: f64) {
    println!("Linearization: {value} + {deriv} * (pressure - {origin})");
}

fn _calc_ideal_gas(pressure: f64, volume: f64, temperature: f64, quantity: f64) -> f64 {
    pressure * volume - quantity * temperature
}

use std::ops::{Mul, Sub};
fn calc_ideal_gas_generic<T>(pressure: T, volume: T, temperature: T, quantity: T) -> T
where
    T: Mul<Output = T> + Sub<Output = T>,
{
    pressure * volume - quantity * temperature
}

fn _calc_ideal_gas_deriv(volume: f64) -> f64 {
    volume
}

fn newton_iterations<Resid>(
    func: Resid,
    initial: f64,
    tolerance: f64,
    max_iter: u8,
) -> Result<f64, ConvergenceError>
where
    Resid: DualFunction,
{
    let mut result = initial;

    let mut calc = None;

    for _ in 0..=max_iter {
        calc = Some(result.derive(&func));

        let error = (calc.unwrap().val - tolerance).abs();

        if error <= tolerance {
            return Ok(result);
        }

        let delta = -calc.unwrap().val / calc.unwrap().dual;

        result += delta;
    }

    Err(ConvergenceError(calc.map_or(None, |x| Some(x.val))))
}
struct ConvergenceError(Option<f64>);