tara 1.0.0

Tara — stellar astrophysics engine for star classification, evolution, nucleosynthesis, and spectral analysis
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
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297

use tara::atmosphere;
use tara::bridge;
use tara::classification::{self, LuminosityClass};
use tara::evolution;
use tara::luminosity;
use tara::nucleosynthesis;
use tara::spectral;
use tara::sse;
use tara::star::{SpectralClass, Star};

#[test]
fn crate_loads() {
    let _ = std::any::type_name::<tara::error::TaraError>();
}

/// End-to-end: build a Sun, verify classification, evolution, and luminosity agree.
#[test]
fn sun_end_to_end() {
    let sun = Star::sun().unwrap();

    // Classification matches
    assert_eq!(sun.spectral_class, SpectralClass::G);
    assert_eq!(sun.spectral_subclass, 2);
    assert_eq!(sun.luminosity_class, LuminosityClass::V);

    // MK string
    let mk = classification::format_classification(sun.temperature_k, sun.luminosity_class);
    assert_eq!(mk, "G2V");

    // Luminosity from radius + temp should be ~1 L_sun
    let l = luminosity::luminosity_solar(sun.radius_solar, sun.temperature_k);
    assert!((l - 1.0).abs() < 0.02, "Solar luminosity: {l}");

    // Absolute bolometric magnitude
    let m_bol = luminosity::absolute_bolometric_magnitude(l);
    assert!((m_bol - 4.74).abs() < 0.1, "Solar M_bol: {m_bol}");

    // Evolutionary phase at solar age
    let phase = evolution::evolutionary_phase(sun.mass_solar, sun.age_years);
    assert_eq!(phase, evolution::EvolutionaryPhase::MainSequence);

    // Remnant type
    assert_eq!(
        evolution::remnant_type(sun.mass_solar),
        evolution::RemnantType::WhiteDwarf
    );
}

/// Build a massive O-star and verify cross-module consistency.
#[test]
fn massive_star_cross_module() {
    let star = Star::builder(30.0, 10.0, 35_000.0, 200_000.0, 3e6)
        .build()
        .unwrap();

    assert_eq!(star.spectral_class, SpectralClass::O);

    // Short-lived
    let lifetime = evolution::main_sequence_lifetime(star.mass_solar);
    assert!(
        lifetime < 1e8,
        "O-star lifetime should be < 100 Myr: {lifetime:.2e}"
    );

    // Should become a black hole
    assert_eq!(
        evolution::remnant_type(star.mass_solar),
        evolution::RemnantType::BlackHole
    );

    // HR region
    let region = classification::hr_region(star.temperature_k, star.luminosity_solar);
    assert_eq!(region, classification::HrRegion::Supergiant);

    // Nucleosynthesis: CNO dominant at ~30 MK core
    assert_eq!(
        nucleosynthesis::dominant_process(30e6),
        nucleosynthesis::FusionProcess::CnoCycle
    );
}

/// Bridge module: round-trip through falak/prakash/tanmatra bridges.
#[test]
fn bridge_round_trips() {
    // Falak: mass → gravitational parameter
    let mu = bridge::stellar_mass_solar_to_mu(1.0);
    assert!(mu > 1e20, "Solar mu: {mu}");

    // Prakash: temperature → peak wavelength
    let peak = bridge::temperature_to_peak_wavelength_nm(5772.0);
    assert!((peak - 502.0).abs() < 5.0, "Solar peak: {peak} nm");

    // Limb darkening from log(g)
    let log_g = atmosphere::log_surface_gravity(1.0, 1.0);
    let u = bridge::surface_gravity_to_limb_darkening(log_g);
    assert!(u > 0.4 && u < 0.7, "Solar limb darkening: {u}");

    // Tanmatra: composition → lifetime scale
    let scale = bridge::composition_to_lifetime_scale(0.74);
    assert!(
        (scale - 1.0).abs() < 0.01,
        "Solar composition scale: {scale}"
    );

    // Tanmatra: core temp → CNO fraction
    let cno = bridge::core_temp_to_cno_fraction(15.7e6);
    assert!(cno < 0.5, "Solar core should be pp-dominated: {cno}");

    // Bridge: luminosity class → abs magnitude
    let m = bridge::luminosity_class_to_abs_magnitude(LuminosityClass::V);
    assert!((m - 5.0).abs() < 0.1);
}

/// Spectral: Planck → Doppler → broadening pipeline.
#[test]
fn spectral_pipeline() {
    // Planck radiance at H-alpha for solar temperature
    let h_alpha_nm = 656.3;
    let radiance = spectral::planck_radiance_nm(h_alpha_nm, 5772.0);
    assert!(radiance > 0.0, "H-alpha radiance: {radiance}");

    // Doppler shift at 100 km/s recession
    let shifted = spectral::doppler_shift(h_alpha_nm, 1e5);
    assert!(shifted > h_alpha_nm, "Redshifted: {shifted}");

    // Thermal broadening of H-alpha
    let sigma = spectral::thermal_broadening(656.3e-9, 5772.0, tara::constants::AMU);
    assert!(sigma > 0.0 && sigma < 1e-9, "Broadening sigma: {sigma}");

    // Voigt profile should produce non-zero value at line center
    let v = spectral::pseudo_voigt_profile(h_alpha_nm, h_alpha_nm, 0.01, 0.005);
    assert!(v > 0.0, "Voigt at center: {v}");
}

/// Atmosphere: T_eff → log(g) → grey atmosphere → opacity chain.
#[test]
fn atmosphere_chain() {
    let t_eff = atmosphere::effective_temperature(tara::constants::L_SUN, tara::constants::R_SUN);
    assert!((t_eff - 5772.0).abs() < 10.0, "Effective temp: {t_eff}");

    let log_g = atmosphere::log_surface_gravity(1.0, 1.0);
    assert!((log_g - 4.44).abs() < 0.02, "log(g): {log_g}");

    // Grey atmosphere: at τ=2/3 should recover ~T_eff
    let t_23 = atmosphere::grey_atmosphere_temperature(t_eff, 2.0 / 3.0);
    // T(τ=2/3) = (3/4 × T_eff⁴ × 4/3)^(1/4) = T_eff
    assert!((t_23 - t_eff).abs() < 1.0, "Grey T at τ=2/3: {t_23}");

    // Electron scattering for solar composition
    let kappa_es = atmosphere::electron_scattering_opacity(0.74);
    assert!((kappa_es - 0.348).abs() < 0.01, "kappa_es: {kappa_es}");

    // Scale height should be order ~100 km
    let g_cgs = 10.0_f64.powf(log_g); // cm/s²
    let g_si = g_cgs / 100.0; // m/s²
    let h = atmosphere::scale_height(t_eff, 1.3, g_si);
    assert!(h > 1e5 && h < 3e5, "Scale height: {h} m");
}

/// Luminosity: distance modulus ↔ apparent/absolute magnitude consistency.
#[test]
fn magnitude_distance_consistency() {
    // Sun at 10 pc: m = M (distance modulus = 0)
    let m_abs = luminosity::absolute_bolometric_magnitude(1.0);
    let m_app = luminosity::apparent_magnitude(m_abs, 10.0);
    assert!(
        (m_app - m_abs).abs() < f64::EPSILON,
        "At 10 pc, m should equal M"
    );

    // At 100 pc, apparent should be 5 mag dimmer
    let m_100 = luminosity::apparent_magnitude(m_abs, 100.0);
    assert!(
        (m_100 - m_abs - 5.0).abs() < 1e-10,
        "100 pc distance modulus: {}",
        m_100 - m_abs
    );

    // Round-trip: magnitude → luminosity → magnitude
    let l = luminosity::luminosity_from_abs_bol_mag(m_abs);
    let m_back = luminosity::absolute_bolometric_magnitude(l);
    assert!(
        (m_back - m_abs).abs() < 1e-10,
        "Magnitude roundtrip: {m_abs} → L={l} → {m_back}"
    );
}

/// SSE: ZAMS → MS evolution → TMS consistency across mass range.
#[test]
fn sse_solar_evolution() {
    let z = sse::Z_SUN;

    // ZAMS properties should be physically reasonable
    let zams = sse::zams_properties(1.0, z);
    assert!(zams.luminosity_solar > 0.5 && zams.luminosity_solar < 1.0);
    assert!(zams.radius_solar > 0.7 && zams.radius_solar < 1.2);
    assert!(zams.temperature_k > 5000.0 && zams.temperature_k < 6500.0);

    // MS lifetime ~10 Gyr
    let t_ms = sse::ms_lifetime(1.0, z);
    assert!(t_ms > 8e9 && t_ms < 13e9, "Solar t_MS: {t_ms:.2e}");

    // Mid-life properties should be between ZAMS and TMS
    let mid = sse::ms_properties(1.0, z, t_ms * 0.5);
    assert!(mid.luminosity_solar >= zams.luminosity_solar);
    assert!(mid.luminosity_solar <= sse::tms_luminosity(1.0, z));

    // TMS should be brighter and larger than ZAMS
    let l_tms = sse::tms_luminosity(1.0, z);
    let r_tms = sse::tms_radius(1.0, z);
    assert!(l_tms > zams.luminosity_solar);
    assert!(r_tms >= zams.radius_solar);
}

/// SSE: metallicity effects across the mass range.
#[test]
fn sse_metallicity_consistency() {
    // Lower Z → brighter ZAMS (less opacity)
    for mass in [0.5, 1.0, 5.0, 10.0] {
        let l_low = sse::zams_luminosity(mass, 0.001);
        let l_high = sse::zams_luminosity(mass, 0.02);
        assert!(
            l_low > l_high,
            "M={mass}: low-Z L={l_low} should exceed high-Z L={l_high}"
        );
    }

    // Lower Z → shorter MS lifetime (higher L means faster fuel burn)
    let t_low = sse::ms_lifetime(1.0, 0.001);
    let t_high = sse::ms_lifetime(1.0, 0.02);
    assert!(
        t_low < t_high,
        "Low-Z lifetime {t_low:.2e} < high-Z {t_high:.2e}"
    );
}

/// SSE: mass-luminosity ordering holds for ZAMS, TMS, and mid-MS.
#[test]
fn sse_mass_ordering() {
    let masses = [0.5, 1.0, 2.0, 5.0, 10.0, 20.0];
    let z = sse::Z_SUN;

    for w in masses.windows(2) {
        let (m_lo, m_hi) = (w[0], w[1]);
        assert!(sse::zams_luminosity(m_hi, z) > sse::zams_luminosity(m_lo, z));
        assert!(sse::tms_luminosity(m_hi, z) > sse::tms_luminosity(m_lo, z));
        assert!(sse::ms_lifetime(m_hi, z) < sse::ms_lifetime(m_lo, z));
    }
}

/// SSE: evolve() lifecycle — phase transitions from MS through remnant.
#[test]
fn sse_evolve_lifecycle() {
    let z = sse::Z_SUN;
    let t_ms = sse::ms_lifetime(1.0, z);
    let t_bgb = sse::t_bgb(1.0, z) * 1e6;

    // Early MS
    let s = sse::evolve(1.0, z, t_ms * 0.1);
    assert_eq!(s.phase, sse::SsePhase::MainSequence);
    assert!(s.luminosity_solar > 0.0);

    // Mid MS — luminosity should be between ZAMS and TMS
    let s = sse::evolve(1.0, z, t_ms * 0.5);
    assert_eq!(s.phase, sse::SsePhase::MainSequence);
    assert!(s.luminosity_solar >= sse::zams_luminosity(1.0, z));
    assert!(s.luminosity_solar <= sse::tms_luminosity(1.0, z));

    // Just past TMS — should be in HG
    let s = sse::evolve(1.0, z, t_ms + (t_bgb - t_ms) * 0.5);
    assert_eq!(s.phase, sse::SsePhase::HertzsprungGap);
    assert!(s.luminosity_solar > sse::tms_luminosity(1.0, z));

    // Past BGB — should be on RGB
    let s = sse::evolve(1.0, z, t_bgb * 1.01);
    assert_eq!(s.phase, sse::SsePhase::RedGiantBranch);
    assert!(s.radius_solar > 1.5, "RGB radius: {}", s.radius_solar);

    // Very old — should be a remnant
    let s = sse::evolve(1.0, z, 1e11);
    assert_eq!(s.phase, sse::SsePhase::WhiteDwarf);

    // Massive star remnant
    let s = sse::evolve(30.0, z, 1e10);
    assert_eq!(s.phase, sse::SsePhase::BlackHole);

    // Intermediate mass remnant
    let s = sse::evolve(12.0, z, 1e10);
    assert_eq!(s.phase, sse::SsePhase::NeutronStar);

    // Luminosity should increase monotonically through MS → HG → RGB
    let l_ms = sse::evolve(1.0, z, t_ms * 0.99).luminosity_solar;
    let l_hg = sse::evolve(1.0, z, t_ms + (t_bgb - t_ms) * 0.99).luminosity_solar;
    let l_rgb = sse::evolve(1.0, z, t_bgb * 1.5).luminosity_solar;
    assert!(l_hg > l_ms, "L_HG > L_MS: {l_hg} vs {l_ms}");
    assert!(l_rgb > l_hg, "L_RGB > L_HG: {l_rgb} vs {l_hg}");
}