sagittariusas 0.0.1

Simulation engine for Sagittarius A* — Kerr spacetime, accretion, jets, lensing, and shadow observables
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

use crate::config::parameters::{C, G, HBAR, K_B, SIGMA_SB};
use sciforge::hub::prelude::astronomy::stellar;
use std::f64::consts::PI;

#[derive(Debug, Clone, Copy)]
pub struct Singularity {
    pub mass: f64,
    pub spin: f64,
    pub charge: f64,
}

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

    pub fn kerr(mass: f64, spin: f64) -> Self {
        Self {
            mass,
            spin: spin.clamp(-1.0, 1.0),
            charge: 0.0,
        }
    }

    pub fn kerr_newman(mass: f64, spin: f64, charge: f64) -> Self {
        Self {
            mass,
            spin: spin.clamp(-1.0, 1.0),
            charge,
        }
    }

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

    pub fn gravitational_radius(&self) -> f64 {
        self.schwarzschild_radius() / 2.0
    }

    pub fn hawking_temperature(&self) -> f64 {
        HBAR * C * C * C / (8.0 * PI * G * self.mass * K_B)
    }

    pub fn hawking_luminosity(&self) -> f64 {
        let t = self.hawking_temperature();
        let rs = self.schwarzschild_radius();
        let area = 4.0 * PI * rs * rs;
        SIGMA_SB * t.powi(4) * area
    }

    pub fn evaporation_time(&self) -> f64 {
        5120.0 * PI * G * G * self.mass.powi(3) / (HBAR * C.powi(4))
    }

    pub fn bekenstein_hawking_entropy(&self) -> f64 {
        let rs = self.schwarzschild_radius();
        let area = 4.0 * PI * rs * rs;
        K_B * C * C * C * area / (4.0 * G * HBAR)
    }

    pub fn penrose_energy_extraction_max(&self) -> f64 {
        let a = self.spin.abs();
        let efficiency = 1.0 - (0.5 * (1.0 + (1.0 - a * a).sqrt())).sqrt();
        efficiency * self.mass * C * C
    }

    pub fn surface_gravity(&self) -> f64 {
        let rg = self.gravitational_radius();
        let r_plus = rg * (1.0 + (1.0 - self.spin * self.spin).sqrt());
        (r_plus - rg) / (2.0 * r_plus * r_plus + 2.0 * rg * self.spin.powi(2) * rg) * C * C
    }

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

    pub fn is_extremal(&self) -> bool {
        (self.spin.abs() - 1.0).abs() < 1e-12
    }

    pub fn is_naked(&self) -> bool {
        let rg = self.gravitational_radius();
        let a = self.spin * rg;
        let epsilon_0 = 8.854_187_817e-12;
        let r_q2 = G * self.charge.powi(2) / (4.0 * PI * epsilon_0 * C.powi(4));
        a * a + r_q2 > rg * rg
    }
}