use suns::physics::energy_transport::{
ConvectiveTransport, RadiativeTransfer, adiabatic_gradient, core_radiative_transfer,
deep_convection, eddington_flux, electron_scattering_opacity, granulation_power_spectrum,
kramers_opacity, limb_darkening, photon_random_walk_time, radiative_equilibrium_temperature,
radiative_zone_transfer, schwarzschild_criterion, surface_convection,
};
use suns::physics::magnetic_fields::{
MagneticFieldRegion, active_region, coronal_field, cyclotron_frequency,
differential_rotation_period_days, differential_rotation_rate, dynamo_number,
electron_cyclotron_frequency, electron_larmor_radius, hale_cycle_period, larmor_radius,
magnetic_energy, magnetic_flux, magnetic_helicity, omega_effect_shear,
proton_cyclotron_frequency, proton_larmor_radius, quiet_sun, solar_dynamo_number_estimate,
sunspot_umbra, sweet_parker_reconnection_rate,
};
use suns::physics::nuclear_fusion::{
FusionReaction, cno_cycle, cno_rate, core_luminosity, coulomb_barrier_mev, gamow_peak_energy,
gamow_window_width, hydrogen_consumption_rate, main_sequence_lifetime_years,
mass_to_energy_rate, neutrino_luminosity, pp_chain_rate, pp1_chain, pp2_chain, pp3_chain,
total_energy_generation_rate, tunneling_probability,
};
use suns::physics::stellar_structure::{
ADIABATIC_INDEX, StellarLayer, central_pressure_estimate, convective_zone, core,
dynamical_timescale, eddington_luminosity, gravitational_binding_energy, hydrostatic_pressure,
jeans_mass, kelvin_helmholtz_timescale, luminosity_from_temperature, nuclear_timescale,
radiative_zone, roche_limit, schwarzschild_radius, surface_gravity, thermal_timescale,
virial_temperature,
};
#[test]
fn core_temperature_around_15_million_k() {
let c = core();
let t = c.mean_temperature();
assert!(t > 5e6 && t < 2e7, "Core mean T={t} not in expected range");
}
#[test]
fn core_inner_radius_is_zero() {
let c = core();
assert!((c.inner_radius_frac - 0.0).abs() < 1e-10);
}
#[test]
fn radiative_zone_outside_core() {
let c = core();
let rz = radiative_zone();
assert!(rz.inner_radius_frac >= c.outer_radius_frac);
}
#[test]
fn convective_zone_extends_to_surface() {
let cz = convective_zone();
assert!((cz.outer_radius_frac - 1.0).abs() < 0.01);
}
#[test]
fn layer_mass_positive() {
for layer in [core(), radiative_zone(), convective_zone()] {
assert!(
layer.mass_estimate_kg() > 0.0,
"{} mass not positive",
layer.name
);
}
}
#[test]
fn layer_volume_increases_outward() {
let c = core();
let rz = radiative_zone();
assert!(rz.volume_m3() > c.volume_m3());
}
#[test]
fn sound_speed_positive() {
for layer in [core(), radiative_zone(), convective_zone()] {
assert!(layer.sound_speed() > 0.0, "{} cs <= 0", layer.name);
}
}
#[test]
fn hydrostatic_pressure_positive() {
let p1 = hydrostatic_pressure(0.1);
let p2 = hydrostatic_pressure(0.5);
assert!(p1 > 0.0);
assert!(p2 > 0.0);
}
#[test]
fn timescale_ordering() {
let t_dyn = dynamical_timescale();
let t_kh = kelvin_helmholtz_timescale();
let t_nuc = nuclear_timescale();
assert!(t_dyn < t_kh, "dynamical should be << KH");
assert!(t_kh < t_nuc, "KH should be << nuclear");
}
#[test]
fn dynamical_timescale_about_30_min() {
let t = dynamical_timescale();
assert!(t > 1000.0 && t < 5000.0, "t_dyn={t}s not ~1800s");
}
#[test]
fn kh_timescale_about_15_myr() {
let t_yr = kelvin_helmholtz_timescale() / 3.15576e7;
assert!(t_yr > 1e6 && t_yr < 1e8, "t_KH={t_yr}yr not ~1.5e7");
}
#[test]
fn nuclear_timescale_about_10_gyr() {
let t_yr = nuclear_timescale() / 3.15576e7;
assert!(t_yr > 1e9 && t_yr < 1e11, "t_nuc={t_yr}yr not ~1e10");
}
#[test]
fn eddington_luminosity_above_solar() {
assert!(eddington_luminosity() > suns::SOLAR_LUMINOSITY);
}
#[test]
fn virial_temperature_positive() {
assert!(virial_temperature() > 0.0);
}
#[test]
fn central_pressure_very_high() {
let p = central_pressure_estimate();
assert!(p > 1e13, "P_central={p} too low");
}
#[test]
fn luminosity_from_temperature_near_solar() {
let l = luminosity_from_temperature();
let ratio = l / suns::SOLAR_LUMINOSITY;
assert!((ratio - 1.0).abs() < 0.1, "L ratio={ratio} too far from 1");
}
#[test]
fn gravitational_binding_energy_order_of_magnitude() {
let e = gravitational_binding_energy();
assert!(e > 1e40 && e < 1e42, "E_bind={e}");
}
#[test]
fn surface_gravity_about_274() {
let g = surface_gravity();
assert!((g - 274.0).abs() < 10.0, "g={g} not ~274");
}
#[test]
fn schwarzschild_radius_about_3km() {
let r = schwarzschild_radius();
assert!(r > 2000.0 && r < 5000.0, "R_s={r}m");
}
#[test]
fn roche_limit_reasonable() {
let r = roche_limit(5000.0);
assert!(r > 0.0);
assert!(r > 0.1 * suns::SOLAR_RADIUS);
}
#[test]
fn jeans_mass_positive() {
let m = jeans_mass(1e4, 1e-20);
assert!(m > 0.0);
}
#[test]
fn pp_chains_branch_fractions_sum_to_roughly_one() {
let sum =
pp1_chain().branch_fraction + pp2_chain().branch_fraction + pp3_chain().branch_fraction;
assert!((sum - 1.0).abs() < 0.05, "pp branches sum={sum}");
}
#[test]
fn pp1_dominant_branch() {
assert!(pp1_chain().branch_fraction > pp2_chain().branch_fraction);
assert!(pp1_chain().branch_fraction > pp3_chain().branch_fraction);
}
#[test]
fn cno_temperature_sensitivity_higher_than_pp() {
assert!(cno_cycle().temperature_sensitivity > pp1_chain().temperature_sensitivity);
}
#[test]
fn net_energy_positive_all_reactions() {
for rxn in [pp1_chain(), pp2_chain(), pp3_chain(), cno_cycle()] {
assert!(rxn.net_energy_mev() > 0.0, "{} net E <= 0", rxn.name);
assert!(rxn.net_energy_joules() > 0.0);
}
}
#[test]
fn pp_rate_increases_with_temperature() {
let r1 = pp_chain_rate(1e7, 1e5, 0.34);
let r2 = pp_chain_rate(1.5e7, 1e5, 0.34);
assert!(r2 > r1, "pp rate should increase with T");
}
#[test]
fn cno_rate_increases_with_temperature() {
let r1 = cno_rate(1e7, 1e5, 0.34, 0.01);
let r2 = cno_rate(1.5e7, 1e5, 0.34, 0.01);
assert!(r2 > r1);
}
#[test]
fn core_luminosity_positive() {
let l = core_luminosity();
assert!(l > 0.0, "core luminosity should be positive");
}
#[test]
fn mass_loss_rate_about_4e9_kgs() {
let dm = mass_to_energy_rate();
assert!(dm > 1e9 && dm < 1e10, "dm/dt={dm}");
}
#[test]
fn main_sequence_lifetime_gyr_scale() {
let t = main_sequence_lifetime_years();
assert!(t > 1e9 && t < 2e10, "MS lifetime={t}yr");
}
#[test]
fn gamow_peak_positive() {
assert!(gamow_peak_energy(1.5e7) > 0.0);
}
#[test]
fn gamow_window_positive_and_comparable_to_peak() {
let peak = gamow_peak_energy(1.5e7);
let width = gamow_window_width(1.5e7);
assert!(width > 0.0);
assert!(width < 10.0 * peak);
}
#[test]
fn coulomb_barrier_positive() {
assert!(coulomb_barrier_mev(1.0, 1.0, 1.2) > 0.0);
}
#[test]
fn tunneling_probability_less_than_one() {
let p = tunneling_probability(0.01, 1.0);
assert!(p > 0.0 && p < 1.0);
}
#[test]
fn neutrino_luminosity_fraction_of_total() {
let l_nu = neutrino_luminosity();
assert!(l_nu > 0.0 && l_nu < suns::SOLAR_LUMINOSITY);
}
#[test]
fn sunspot_field_stronger_than_quiet_sun() {
assert!(sunspot_umbra().field_strength_t > quiet_sun().field_strength_t);
}
#[test]
fn active_region_between_quiet_and_sunspot() {
let q = quiet_sun().field_strength_t;
let a = active_region().field_strength_t;
let s = sunspot_umbra().field_strength_t;
assert!(a > q && a < s);
}
#[test]
fn plasma_beta_quiet_sun_above_one() {
assert!(quiet_sun().plasma_beta() > 1.0, "Quiet sun β should be > 1");
}
#[test]
fn alfven_speed_positive() {
for r in [
quiet_sun(),
active_region(),
sunspot_umbra(),
coronal_field(),
] {
assert!(r.alfven_speed() > 0.0, "{} v_A <= 0", r.name);
}
}
#[test]
fn differential_rotation_equator_faster() {
let omega_eq = differential_rotation_rate(0.0);
let omega_60 = differential_rotation_rate(60.0);
assert!(omega_eq > omega_60, "Equator should rotate faster");
}
#[test]
fn rotation_period_equator_about_25_days() {
let p = differential_rotation_period_days(0.0);
assert!((p - 25.0).abs() < 3.0, "P_eq={p}d");
}
#[test]
fn rotation_period_increases_with_latitude() {
let p_eq = differential_rotation_period_days(0.0);
let p_60 = differential_rotation_period_days(60.0);
assert!(p_60 > p_eq);
}
#[test]
fn omega_effect_shear_positive() {
assert!(omega_effect_shear() > 0.0);
}
#[test]
fn dynamo_number_estimate_nonzero() {
assert!(solar_dynamo_number_estimate().abs() > 0.0);
}
#[test]
fn magnetic_flux_proportional_to_field() {
let phi1 = magnetic_flux(1.0, 1e6);
let phi2 = magnetic_flux(2.0, 1e6);
assert!((phi2 / phi1 - 2.0).abs() < 0.01);
}
#[test]
fn magnetic_energy_proportional_to_b_squared() {
let e1 = magnetic_energy(1.0, 1e6);
let e2 = magnetic_energy(2.0, 1e6);
assert!((e2 / e1 - 4.0).abs() < 0.01);
}
#[test]
fn proton_larmor_larger_than_electron() {
let r_p = proton_larmor_radius(1e6, 0.01);
let r_e = electron_larmor_radius(1e6, 0.01);
assert!(r_p > r_e);
}
#[test]
fn electron_cyclotron_faster_than_proton() {
let f_p = proton_cyclotron_frequency(0.01);
let f_e = electron_cyclotron_frequency(0.01);
assert!(f_e > f_p);
}
#[test]
fn reconnection_rate_subluminal() {
let v = sweet_parker_reconnection_rate(1e6, 1e12);
assert!(v > 0.0 && v < 3e8);
}
#[test]
fn hale_cycle_about_22_years() {
let p = hale_cycle_period();
assert!((p - 22.0).abs() < 1.0, "Hale period={p}yr");
}
#[test]
fn core_mean_free_path_very_small() {
let mfp = core_radiative_transfer().mean_free_path();
assert!(mfp > 0.0 && mfp < 1.0, "core MFP={mfp}m");
}
#[test]
fn radiative_zone_mfp_larger_than_core() {
let core_mfp = core_radiative_transfer().mean_free_path();
let rz_mfp = radiative_zone_transfer().mean_free_path();
assert!(rz_mfp > core_mfp);
}
#[test]
fn photon_diffusion_speed_subluminal() {
let v = core_radiative_transfer().photon_diffusion_speed();
assert!(v > 0.0 && v < 3e8);
}
#[test]
fn convective_velocity_reasonable() {
let v = surface_convection().convective_velocity();
assert!(v > 10.0 && v < 1e5, "v_conv={v}");
}
#[test]
fn convective_flux_positive() {
assert!(surface_convection().convective_flux() > 0.0);
}
#[test]
fn photon_random_walk_time_thousands_of_years() {
let t = photon_random_walk_time() / 3.15576e7;
assert!(t > 1e3 && t < 1e7, "walk time={t}yr");
}
#[test]
fn schwarzschild_criterion_convective_when_grad_rad_large() {
assert!(schwarzschild_criterion(1.0, adiabatic_gradient()));
}
#[test]
fn schwarzschild_criterion_radiative_when_grad_rad_small() {
assert!(!schwarzschild_criterion(0.1, adiabatic_gradient()));
}
#[test]
fn kramers_opacity_decreases_with_temperature() {
let k1 = kramers_opacity(1e3, 1e6);
let k2 = kramers_opacity(1e3, 1e7);
assert!(k2 < k1);
}
#[test]
fn electron_scattering_opacity_proportional_to_x() {
let k1 = electron_scattering_opacity(0.5);
let k2 = electron_scattering_opacity(0.7);
assert!(k2 > k1);
}
#[test]
fn limb_darkening_center_brightest() {
let center = limb_darkening(0.0);
let edge = limb_darkening(1.4);
assert!(center > edge);
}
#[test]
fn limb_darkening_at_center_is_one() {
let i = limb_darkening(0.0);
assert!((i - 1.0).abs() < 0.01, "I(0)={i}");
}
#[test]
fn adiabatic_gradient_about_0_4() {
let g = adiabatic_gradient();
assert!((g - 0.4).abs() < 0.05, "∇_ad={g}");
}
#[test]
fn deep_convection_slower_velocity() {
let vs = surface_convection().convective_velocity();
let vd = deep_convection().convective_velocity();
assert!(vd > 0.0 && vs > 0.0);
}