blackholesfactory 0.0.2

Black hole factory — create and simulate stellar, intermediate-mass, and supermassive black holes with full Kerr physics
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

use crate::config::parameters::{C, G, gravitational_radius, schwarzschild_radius};

pub struct Singularity {
    pub mass: f64,
    pub spin: f64,
    pub charge: f64,
}

impl Singularity {
    pub fn new(mass: f64) -> Self {
        Self {
            mass,
            spin: 0.0,
            charge: 0.0,
        }
    }

    pub fn with_spin(mut self, spin: f64) -> Self {
        self.spin = spin.clamp(-1.0, 1.0);
        self
    }

    pub fn with_charge(mut self, charge: f64) -> Self {
        self.charge = charge;
        self
    }

    pub fn gravitational_radius(&self) -> f64 {
        gravitational_radius(self.mass)
    }

    pub fn schwarzschild_radius(&self) -> f64 {
        schwarzschild_radius(self.mass)
    }

    pub fn event_horizon(&self) -> f64 {
        let rg = self.gravitational_radius();
        let a = self.spin * rg;
        let rq2 = G * self.charge * self.charge / (C * C * C * C);
        rg + (rg * rg - a * a - rq2).max(0.0).sqrt()
    }

    pub fn cauchy_horizon(&self) -> f64 {
        let rg = self.gravitational_radius();
        let a = self.spin * rg;
        let rq2 = G * self.charge * self.charge / (C * C * C * C);
        rg - (rg * rg - a * a - rq2).max(0.0).sqrt()
    }

    pub fn ergosphere_radius(&self, theta: f64) -> f64 {
        let rg = self.gravitational_radius();
        let a = self.spin * rg;
        rg + (rg * rg - a * a * theta.cos().powi(2)).max(0.0).sqrt()
    }

    pub fn is_naked(&self) -> bool {
        let rg = self.gravitational_radius();
        let a = self.spin * rg;
        let rq2 = G * self.charge * self.charge / (C * C * C * C);
        a * a + rq2 > rg * rg
    }

    pub fn surface_gravity(&self) -> f64 {
        let rh = self.event_horizon();
        let rg = self.gravitational_radius();
        let a = self.spin * rg;
        let r_inner = (rh * rh - a * a).max(0.0).sqrt();
        if rh + r_inner < 1e-30 {
            return 0.0;
        }
        (rh - rg) / (2.0 * rh * (rh + r_inner))
    }

    pub fn angular_momentum(&self) -> f64 {
        self.spin * G * self.mass * self.mass / C
    }

    pub fn irreducible_mass(&self) -> f64 {
        let a2 = self.spin * self.spin;
        self.mass / 2.0_f64.sqrt() * (1.0 + (1.0 - a2).max(0.0).sqrt()).sqrt()
    }

    pub fn penrose_energy_extraction(&self) -> f64 {
        (self.mass - self.irreducible_mass()) * C * C
    }

    pub fn entropy(&self) -> f64 {
        use crate::config::parameters::{HBAR, K_B};
        let rh = self.event_horizon();
        let rg = self.gravitational_radius();
        let a = self.spin * rg;
        let area = 4.0 * std::f64::consts::PI * (rh * rh + a * a);
        K_B * C * C * C * area / (4.0 * HBAR * G)
    }
}