Code shared from the Rust Playground

playground.rs

// Filename: archimedes_pi.rs

/// Approximates pi using Archimedes' method by iteratively calculating
/// the semi-perimeter of an inscribed regular polygon with 2^(k+1) sides
/// in a unit circle (r=1).
///
/// The formula used here is an iterative approach for the semi-perimeter (P_n),
/// where P_n is the semi-perimeter of the n-sided polygon.
///
/// # Arguments
/// * `iterations`: The number of times (k) to double the polygon's sides,
///                 starting from a square (k=1).
///
/// # Returns
/// The approximated value of pi as an `f64`.
pub fn archimedes_pi_approximation(iterations: u32) -> f64 {
    // 1. Initial State: Square (n=4), which is k=1 in the limit formula.

    // The side length of a square inscribed in a unit circle (r=1) is sqrt(2).
    // The semi-perimeter (P_4) is 4 * side_length / 2 = 2 * sqrt(2).
    let s_4 = 2.0_f64.sqrt();
    let mut current_semi_perimeter: f64 = 4.0 * s_4 / 2.0; // P_4 = 2 * sqrt(2)
    let mut current_sides: f64 = 4.0; // n=4

    // The approximation is already P_n, which approaches pi.
    let mut pi_approx = current_semi_perimeter;

    // 2. Iteratively double the number of sides (k)
    for k in 1..=iterations {
        // The side length of the 2n-gon (s_2n) can be found from the side length of the n-gon (s_n).
        // The side length s_n relates to the semi-perimeter P_n by: s_n = 2 * P_n / n
        // The current_semi_perimeter (P_n) is the current approximation for pi.

        // This is a simplified recurrence for the side length for unit radius (r=1):
        // s_{2n} = sqrt(2 - 2 * sqrt(1 - (s_n/2)^2))
        // P_{2n} = n * s_{2n}

        // Let x = s_n / 2.
        let x = pi_approx / current_sides;

        // Calculate s_2n.
        let s_2n = (2.0 - 2.0 * (1.0 - x * x).sqrt()).sqrt();

        // Calculate the new semi-perimeter P_2n = (2n) * s_2n / 2 = n * s_2n
        // We update the number of sides and the semi-perimeter.
        current_sides *= 2.0;
        pi_approx = current_sides * s_2n / 2.0;
    }

    pi_approx
}

// Main function to demonstrate the approximation for various iteration values.
fn main() {
    let exact_pi = std::f64::consts::PI;
    println!("Exact π (f64): {}", exact_pi);
    println!("--------------------------------------");

    // Iterations from k=1 (square) to k=10 (2^11 = 2048 sides)
    for k in 1..=30 {
        let pi_approx = archimedes_pi_approximation(k);
        let error = (exact_pi - pi_approx).abs();

        // The number of sides n is 2^(k+1)
        let n = 2_u32.pow(k + 1);

        println!(
            "k = {:2} (n = {:4}): π_k = {:.15}, Error = {:.15e}",
            k, n, pi_approx, error
        );
    }
}