parry2d/mass_properties/
mass_properties_convex_polygon.rs

1use crate::mass_properties::MassProperties;
2use crate::math::{Point, Real};
3use crate::shape::Triangle;
4
5impl MassProperties {
6    /// Computes the mass properties of a convex polygon.
7    pub fn from_convex_polygon(density: Real, vertices: &[Point<Real>]) -> MassProperties {
8        let (area, com) = convex_polygon_area_and_center_of_mass(vertices);
9
10        if area == 0.0 {
11            return MassProperties::new(com, 0.0, 0.0);
12        }
13
14        let mut itot = 0.0;
15
16        let mut iterpeek = vertices.iter().peekable();
17        let first_element = *iterpeek.peek().unwrap(); // store first element to close the cycle in the end with unwrap_or
18        while let Some(elem) = iterpeek.next() {
19            let triangle = Triangle::new(com, *elem, **iterpeek.peek().unwrap_or(&first_element));
20            let area = triangle.area();
21            let ipart = triangle.unit_angular_inertia();
22            itot += ipart * area;
23        }
24
25        Self::new(com, area * density, itot * density)
26    }
27}
28
29/// Computes the area and center-of-mass of a convex polygon.
30pub fn convex_polygon_area_and_center_of_mass(
31    convex_polygon: &[Point<Real>],
32) -> (Real, Point<Real>) {
33    let geometric_center = convex_polygon
34        .iter()
35        .fold(Point::origin(), |e1, e2| e1 + e2.coords)
36        / convex_polygon.len() as Real;
37    let mut res = Point::origin();
38    let mut areasum = 0.0;
39
40    let mut iterpeek = convex_polygon.iter().peekable();
41    let first_element = *iterpeek.peek().unwrap(); // Stores first element to close the cycle in the end with unwrap_or.
42    while let Some(elem) = iterpeek.next() {
43        let (a, b, c) = (
44            elem,
45            iterpeek.peek().unwrap_or(&first_element),
46            &geometric_center,
47        );
48        let area = Triangle::new(*a, **b, *c).area();
49        let center = (a.coords + b.coords + c.coords) / 3.0;
50
51        res += center * area;
52        areasum += area;
53    }
54
55    if areasum == 0.0 {
56        (areasum, geometric_center)
57    } else {
58        (areasum, res / areasum)
59    }
60}