Center and radius of a circle passing through 3 points in a plane

[PascalVollmer]

Hi,
I'm a Rust beginner and as a first exercise I chose Centre and radius of a circle passing through 3 points in a plane - Rosetta Code. There is no Rust version available and maybe I could contribute to one.
Kind regards, Pascal

// Version
//
// V 0.1.3 2025-08-19 Code coverage check with cargo-llvm-cov (note: unclear how to ignore functions and modules)
// V 0.1.2 2025-08-18 check for float overflow when coordinate values are squared;
//                    visibility of parent module from as seen from module tests
// V 0.1.1 2025-08-17 use of Result in fn circle_through_3_points_2_d(..)
// V 0.1.0 2025-08-17 sketch with two tests that pass
//
// Write a function which returns the centre and radius of a circle passing through three points in a plane.
// Demonstrate the function using the points (22.83,2.07) (14.39,30.24) and (33.65,17.31)
// (from https://rosettacode.org/wiki/Centre_and_radius_of_a_circle_passing_through_3_points_in_a_plane# 2025-08-13)

use core::f64;
use std::f64::MAX;
use num::complex::{Complex, ComplexFloat};

struct Circle {
    // wondering why it is not possible to nest structures. I would have preferred to define a Circle from its
    // center (a Point) and radius.
    center_x: f64,
    center_y: f64,
    radius: f64
}

fn circle_through_3_points_2_d(p1x: f64, p1y: f64, p2x: f64, p2y: f64, p3x: f64, p3y: f64) -> Result<Circle, String> {
    // following https://math.stackexchange.com/questions/213658/get-the-equation-of-a-circle-when-given-3-points
    // following https://web.archive.org/web/20060909065617/http://www.math.okstate.edu/~wrightd/INDRA/MobiusonCircles/node4.html

    // interpret the three points (p1,p2,p3) as being in the complex plane (z1,z2,z3) and 
    // map the two points z1 and z2 to 0 and to 1 in C with z -> (z - z1) / (z2 - z1).
    // This maps p1 to z1=0, p2 to z2=1, p3 to z3.
    // >> would execution time bei smaller when avoiding complex numbers?
    let z1 = num::complex::Complex::new(p1x, p1y);
    let z2 = num::complex::Complex::new(p2x, p2y);
    let z3: Complex<f64> = num::complex::Complex::new(p3x, p3y);

    // check for identity of two points
    if (z1 == z2) || (z2 == z3) || (z3 == z1) {
        // these are not three distinct vertices of a triangle
        // assumptions are not met: refuse computation
        return Err("at least two of three points are identical".to_string())
    }

    // check for overflow: squares will be computed later, make sure this doesn't cause an overflow 
    if 
    z1.re().abs() > MAX.sqrt() || z1.im().abs() > MAX.sqrt() || 
    z2.re().abs() > MAX.sqrt() || z2.im().abs() > MAX.sqrt() || 
    z3.re().abs() > MAX.sqrt()|| z3.im().abs() > MAX.sqrt() {
        return Err("at least one coordinat > MAX.sqrt()".to_string())
    }  

    // now we transform in the complex plane zi -> zti such that
    // zt1 = 0, zt2 = 1, zt3 = to be computed
    let zt3: Complex<f64> = (z3 - z1) / (z2 - z1);

    // check for colinearity
    if zt3.im == 0.0 { 
        return Err("three points are colinear".to_string())
    }

    // the transformed center of the circle is on the line x=0.5, all y except y=0.0. 
    // Determine its specific position
    let ct: Complex<f64> = (zt3 - zt3 * zt3.conj()) / (zt3 - zt3.conj());
    let c: Complex<f64> = (z2-z1) * ct + z1;
    let r: f64 = (z1-c).norm();

    let c = Circle { center_x: c.re, center_y: c.im, radius: r };
    return Ok(c)
}

fn main() {
    let result = circle_through_3_points_2_d(22.83, 2.07, 14.39, 30.24, 33.65, 17.31);
    match result {
        Ok(circle) => {
            println!("circle.center_x={}", circle.center_x); 
            println!("circle.center_y={}", circle.center_y); 
            println!("circle.radius={}", circle.radius);
        }
        Err(e) => println!("{}",e)
    }
}

[CodesInChaos]

I don't like all those separate some_x, some_y values. Those should be combined somehow. For example:

• Use a Vec2 / Point2 type from an appropriate library
• Use a (f32, f32) tuple or [f32; 2] array
• Create your own type:
  • Point(pub f32, pub f32)
  • Point([f32; 2])
  • Point{ pub x:f32, pub y:f32 }

[anon75658756]

Not a mathematics comment, just a few simple modifications (that surely can also be improved.)

I removed the comments, sorry, it was just for reading (i was not planning to reply originally.)

The changes I suggest are:

• Use the imports that you have, not the full path
• Use a nested struct.
• Define a Point used also for the centre.
• Move safety checks to the top
• Use eprintln! for printing errors.

(Note that the safety checks can be done perfectly well before converting to complex numbers)

You can also define a Display method for your struct, so that it prints from there.

I added some of those here

use core::f64;
use std::f64::MAX;
use num::complex::{Complex, ComplexFloat};

#[derive(PartialEq,Debug)]
struct Point {
    x: f64,
    y: f64,
}

struct Circle {
    center: Point,
    radius: f64,
}

fn circle_through_3_points_2_d(p1: Point, p2: Point, p3: Point) -> Result<Circle, String> {
    fn is_unsafe_square(p: &Point) -> bool {
        p.x.abs() > MAX.sqrt() || p.y.abs() > MAX.sqrt()
    }
    // check for overflow: squares will be computed later, make sure this doesn't cause an overflow
    if is_unsafe_square(&p1) || is_unsafe_square(&p2) || is_unsafe_square(&p3) {
        return Err("at least one coordinate > MAX.sqrt()".to_string());
    }

    // check for identity of two points
    if (p1 == p2) || (p2 == p3) || (p3 == p1) {
        // these are not three distinct vertices of a triangle
        // assumptions are not met: refuse computation
        return Err("at least two of three points are identical".to_string());
    }

    let z1 = Complex::new(p1.x, p1.y);
    let z2 = Complex::new(p2.x, p2.y);
    let z3 = Complex::new(p3.x, p3.y);


    // now we transform in the complex plane zi -> zti such that
    // zt1 = 0, zt2 = 1, zt3 = to be computed
    let zt3: Complex<f64> = (z3 - z1) / (z2 - z1);

    // check for colinearity
    if zt3.im == 0.0 {
        return Err("three points are colinear".to_string());
    }

    // the transformed center of the circle is on the line x=0.5, all y except y=0.0.
    // Determine its specific position
    let ct: Complex<f64> = (zt3 - zt3 * zt3.conj()) / (zt3 - zt3.conj());
    let c: Complex<f64> = (z2 - z1) * ct + z1;
    let r: f64 = (z1 - c).norm();

    let c = Circle {
        center: Point{x:c.re, y:c.im},
        radius: r,
    };
    return Ok(c);
}

fn main() {
    let result = circle_through_3_points_2_d(Point {x:22.83, y:2.07}, Point{x: 14.39, y:30.24}, Point{x:33.65, y:17.31});
    match result {
        Ok(circle) => {
            println!("circle.center_x={}", circle.center.x);
            println!("circle.center_y={}", circle.center.y);
            println!("circle.radius={}", circle.radius);
        }
        Err(e) => eprintln!("{}", e),
    }
}

You can implement a printer as well:

use core::f64;
use num::complex::{Complex, ComplexFloat};
use std::f64::MAX;
use std::fmt;

#[derive(PartialEq)]
struct Point {
    x: f64,
    y: f64,
}

struct Circle {
    center: Point,
    radius: f64,
}
impl fmt::Display for Circle {
    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
        let Circle {
            center:
                Point {
                    x: cx,
                    y: cy,
                },
            radius: r,
        } = &self;
        let result = format!(
            "circle.center_x={cx},\ncircle.center_y={cy},\ncircle.radius={r}
"
        );
        write!(f, "{}", result)
    }
}
fn circle_through_3_points_2_d(p1: Point, p2: Point, p3: Point) -> Result<Circle, String> {
    fn is_unsafe_square(p: &Point) -> bool {
        p.x.abs() > MAX.sqrt() || p.y.abs() > MAX.sqrt()
    }
    // check for overflow: squares will be computed later, make sure this doesn't cause an overflow
    if is_unsafe_square(&p1)
        ||(is_unsafe_square(&p2))
        ||(is_unsafe_square(&p3))
    {
        return Err("at least one coordinate > MAX.sqrt()".to_string());
    }

    // check for identity of two points
    if (p1 == p2) || (p2 == p3) || (p3 == p1) {
        // these are not three distinct vertices of a triangle
        // assumptions are not met: refuse computation
        return Err("at least two of three points are identical".to_string());
    }

    let z1 = Complex::new(p1.x, p1.y);
    let z2 = Complex::new(p2.x, p2.y);
    let z3 = Complex::new(p3.x, p3.y);

    // now we transform in the complex plane zi -> zti such that
    // zt1 = 0, zt2 = 1, zt3 = to be computed
    let zt3: Complex<f64> = (z3 - z1) / (z2 - z1);

    // check for colinearity
    if zt3.im == 0.0 {
        return Err("three points are colinear".to_string());
    }

    // the transformed center of the circle is on the line x=0.5, all y except y=0.0.
    // Determine its specific position
    let ct: Complex<f64> = (zt3 - zt3 * zt3.conj()) / (zt3 - zt3.conj());
    let c: Complex<f64> = (z2 - z1) * ct + z1;
    let r: f64 = (z1 - c).norm();

    let c = Circle {
        center: Point { x: c.re, y: c.im },
        radius: r,
    };
    return Ok(c);
}

fn main() {
    let result = circle_through_3_points_2_d(
        Point { x: 22.83, y: 2.07 },
        Point { x: 14.39, y: 30.24 },
        Point { x: 33.65, y: 17.31 },
    );
    match result {
        Ok(circle) => println!("{circle}"),
        Err(e) => eprintln!("{}", e),
    }
}