thimni 0.2.4

efficient SDF collision without discretizatio, neural nets, or interval arithmetic
Documentation 
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
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168

use crate::vector::Vector;
use kiddo::{KdTree, SquaredEuclidean};

#[derive(Copy, Clone)]
pub struct AABB<const N: usize, T: Vector<N>> {
    pub min: T,
    pub max: T,
}

impl<const N: usize, T: Vector<N>> AABB<N, T> {
    pub fn get_dim(&self) -> T {
        self.max - self.min
    }

    pub fn overlap(&self, other: &Self) -> bool {
        self.min.less_than(&other.max) && self.max.greater_than(&other.min)
    }

    pub fn intersect(&self, other: &Self) -> Option<Self> {
        if self.overlap(other) {
            Some(Self {
                min: self.min.cmax(other.min),
                max: self.max.cmin(other.max),
            })
        } else {
            None
        }
    }

    pub fn get_center(&self) -> T {
        self.max.mid_point(self.min)
    }
}

#[derive(Copy, Clone)]
pub struct Circle<const N: usize, T: Vector<N>> {
    pub c: T,
    pub r: f32,
}

impl<const N: usize, T: Vector<N>> Circle<N, T> {
    pub fn new(c: T, r: f32) -> Self {
        Self { c, r }
    }
}

fn factorial(number: usize) -> usize {
    let mut factorial: usize = 1;

    for i in 1..(number + 1) {
        factorial *= i;
    }

    return factorial;
}

pub fn generate_circles<const N: usize, T: Vector<N>>(
    rect: &AABB<N, T>,
    k: f32,
    max_iter: usize,
    min_radius: f32,
    candidates_per_iter: usize,
) -> Vec<Circle<N, T>> {
    let dim = rect.get_dim();

    let target_area = k * dim.cprod();
    let mut covered_area = 0.0;
    let mut circles: Vec<Circle<N, T>> = vec![];

    let mut tree: KdTree<f32, N> = KdTree::new();

    let mut itera = 0;

    while itera < max_iter && covered_area < target_area {
        let mut best_candidate = None;
        let mut best_radius = 0.0;

        for _ in 0..candidates_per_iter {
            let np = dim.cmap_mut(|x: f32| quad_rand::gen_range(0.0, x));

            let dp = np.cmin(dim - np);

            let edge_dist = dp.min_elem();

            let max_radius = if circles.len() > 0 {
                let nearest = tree.nearest_one::<SquaredEuclidean>(&np.to_arr());
                let nearest_circle = circles[nearest.item as usize];

                let radius_from_circles = nearest.distance - nearest_circle.r;
                edge_dist.min(radius_from_circles)
            } else {
                edge_dist
            };

            if max_radius >= min_radius && max_radius > best_radius {
                best_radius = max_radius;
                best_candidate = Some(np.to_arr());
            }
        }

        if let Some(best) = best_candidate {
            circles.push(Circle::new(T::from_arr(best), best_radius));

            if N == 2 {
                covered_area += std::f32::consts::PI * best_radius * best_radius;
            } else if N == 3 {
                covered_area +=
                    (4.0 / 3.0) * std::f32::consts::PI * best_radius * best_radius * best_radius;
            } else {
                covered_area += (std::f32::consts::PI.powf(N as f32 / 2.0)
                    * best_radius.powi(N as i32))
                    / factorial(N.div_ceil(2) + 1) as f32;
            }

            tree.add(&best, tree.size());
        }

        itera += 1;
    }

    circles
}

#[derive(Clone)]
pub struct CollisionParameters {
    /// step size used when approximating the gradient
    pub normal_epsilon: f32,
    /// minimum distance from surface that is considered a hit
    pub collision_epsilon: f32,
    /// learning rate for gradient descent
    pub learning_rate: f32,
    /// maximum number of iterations of gradeint descent
    pub max_gds_iter: usize,
    /// maximum number of iterations of sphere/circle packing
    pub max_packing_iter: usize,
    /// target area percentage to be covered by sphere/circle packing. range 0 to 1
    pub area_percentage: f32,
    /// minimum radius to be used during sphere/circle packing
    pub minimum_radius: f32,
    /// number of samples used for monte-carlo sampling when approximating
    /// exact collision depth
    pub monte_carlo_samples: usize,
    /// [0,1] how strongly the random values used for monte-carlo are biased
    /// towards the point of collision
    pub monte_carlo_bias_factor: f32,
}

impl Default for CollisionParameters {
    fn default() -> Self {
        Self {
            normal_epsilon: 2.0,
            collision_epsilon: 1.0,
            learning_rate: 0.2,
            max_gds_iter: 100,
            max_packing_iter: 50,
            area_percentage: 0.95,
            minimum_radius: 0.1,
            monte_carlo_samples: 200,
            monte_carlo_bias_factor: 0.6,
        }
    }
}

pub(crate) fn random_biased(min: f32, max: f32, bias: f32, factor: f32) -> f32 {
    let rnd = quad_rand::gen_range(min, max);
    let mix = quad_rand::gen_range(0.0, 1.0) * factor;
    rnd * (1.0 - mix) + bias * mix
}