Skip to content

Commit

Permalink
Fix circumsphere function -- it works!
Browse files Browse the repository at this point in the history
- Updated the version of the "syn" package from 2.0.42 to 2.0.43 in Cargo.lock.
- Modified the Cell struct in cell.rs to implement the Div trait for T and Sub trait for &'a T, where T is a generic type parameter.
- Added a conversion implementation From<[T; D]> for Point<f64, D> in point.rs to convert an array of type [T; D] to [f64; D].
- Fixed the calculation of b vector in cell.rs by using na::distance_squared instead of manually calculating the sum of squared coordinates.
- Adjusted the solution array by dividing each element by 2.0 before assigning it to solution_point in cell.rs.
- Updated the expected value for circumcenter test case in tests module of cell.rs.

Note: Removed commented out code and debug print statements.
  • Loading branch information
acgetchell committed Dec 28, 2023
1 parent e7204f0 commit 9fac2c1
Show file tree
Hide file tree
Showing 3 changed files with 46 additions and 40 deletions.
6 changes: 3 additions & 3 deletions Cargo.lock

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

68 changes: 31 additions & 37 deletions src/delaunay_core/cell.rs
Original file line number Diff line number Diff line change
@@ -1,6 +1,9 @@
//! Data and operations on d-dimensional cells or [simplices](https://en.wikipedia.org/wiki/Simplex).
use std::{fmt::Debug, ops::Sub};
use std::{
fmt::Debug,
ops::{Div, Sub},
};

use uuid::Uuid;

Expand Down Expand Up @@ -44,6 +47,9 @@ impl<
V,
const D: usize,
> Cell<T, U, V, D>
where
for<'a> &'a T: Div<f64>,
f64: From<T>,
{
/// The function `new` creates a new `Cell`` object with the given vertices.
/// A D-dimensional cell has D + 1 vertices, so the number of vertices must be less than or equal to D + 1.
Expand Down Expand Up @@ -237,49 +243,43 @@ impl<
}
}

let _a_inv = matrix.try_inverse().ok_or("Singular matrix!")?;
let a_inv = matrix.try_inverse().ok_or("Singular matrix!")?;

println!("Matrix A: {:?}", matrix);
for i in 0..dim {
println!("Row sum of {}-th row is {}", i, matrix.row(i).row_sum());
}

// let mut b: na::DVector<<T as ComplexField>::RealField> = na::DVector::zeros(dim);
// for i in 0..dim {
// b[i] = na::distance_squared(
// &na::Point::from(self.vertices[i + 1].point.coords),
// &na::Point::from(self.vertices[0].point.coords),
// );
// }

let mut b: na::SMatrix<T, D, 1> = na::SMatrix::zeros();

for i in 0..dim {
b[i] = self.vertices[i + 1]
.point
.coords
.iter()
.map(|x| x.powi(2))
.sum::<T>()
- self.vertices[0]
.point
.coords
.iter()
.map(|x| x.powi(2))
.sum::<T>();
b[i] = na::distance_squared(
&na::Point::from(self.vertices[i + 1].point.coords),
&na::Point::from(self.vertices[0].point.coords),
);
}

println!("Vector b: {:?}", b);

// let _solution = a_inv * b;
let solution = a_inv * b;

let mut solution_array: [T; D] = solution
.data
.as_slice()
.try_into()
.expect("Failed to convert solution to array");

// let solution_array: [f64; D] = _solution
// .data
// .as_slice()
// .try_into()
// .expect("Failed to convert solution to array");
// Ok(Point::new(solution_array))
for i in 0..dim {
solution_array[i] /= na::convert::<f64, T>(2.0);
}

Ok(Point::<f64, D>::origin())
println!("Solution: {:?}", solution_array);

let solution_point: Point<f64, D> = Point::<f64, D>::try_from(solution_array)
.expect("Failed to convert solution array to Point");

// let solution_point = Point::<f64, D>::try_from(solution_array);
Ok(Point::<f64, D>::new(solution_point.coords))
}

/// The function `circumsphere_contains` checks if a given vertex is contained in the circumsphere of the Cell.
Expand Down Expand Up @@ -437,15 +437,9 @@ mod tests {
let cell: Cell<f64, i32, Option<()>, 3> =
Cell::new(vec![vertex1, vertex2, vertex3, vertex4]).unwrap();

// let circumcenter = cell
// .circumcenter_wrong(&[point1, point2, point3, point4])
// .unwrap();

// assert_eq!(circumcenter, Point::new([0.5, 0.5, 0.5]));

let circumcenter = cell.circumcenter().unwrap();

assert_eq!(circumcenter, Point::origin());
assert_eq!(circumcenter, Point::new([0.5, 0.5, 0.5]));

// Human readable output for cargo test -- --nocapture
println!("Circumcenter: {:?}", circumcenter);
Expand Down
12 changes: 12 additions & 0 deletions src/delaunay_core/point.rs
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,18 @@ pub struct Point<T, const D: usize> {
pub coords: [T; D],
}

impl<T, const D: usize> From<[T; D]> for Point<f64, D>
where
[T; D]: Sized,
[f64; D]: Sized,
T: Into<f64>,
{
fn from(coords: [T; D]) -> Self {
let coords_f64: [f64; D] = coords.map(|coord| coord.into()).into(); // Convert the `coords` array to `[f64; D]`
Self { coords: coords_f64 }
}
}

impl<T: Clone, const D: usize> Point<T, D> {
/// The function `new` creates a new instance of a `Point` with the given coordinates.
///
Expand Down

0 comments on commit 9fac2c1

Please sign in to comment.