mercurys 0.0.3

Mercury celestial simulation crate for the MilkyWay SolarSystem workspace
Documentation
use mercurys::{ORBITAL_PERIOD_S, SECONDS_PER_YEAR, SEMI_MAJOR_AXIS_M, SOLAR_DAY_S};
use mercurys::{
    atmosphere, geodata, geology, hydrology, lighting, physics, rendering, satellites, temporal,
    terrain,
};

const DEFAULT_DT: f64 = 3600.0;
const STEP_COUNT: usize = 240;

struct MercuryDiagnostics {
    elapsed_years: f64,
    subsolar_longitude_deg: f64,
    solar_flux_w_m2: f64,
    daytime: bool,
    surface_temp_k: f64,
    thermal_peak_wavelength_m: f64,
    thermal_emission: [f32; 3],
    solar_tide_bulge_m: f64,
    tidal_dissipation_w: f64,
    solar_wind_pressure_pa: f64,
    sodium_tail_m: f64,
    cme_enhancement: f64,
    gardening_depth_m: f64,
    darkening_factor: f64,
    sputtering_rate_kg_s: f64,
    thermal_crack_growth: f64,
    volcanic_coverage: f64,
    hollows_count: usize,
    total_hollows_depth_m: f64,
    polar_ice_mass_kg: f64,
    polar_deposit_count: usize,
    psr_temperature_k: f64,
    leaf_count: usize,
    orbiters_mean_period_h: f64,
    observer_altitude_m: f64,
    checksum: f64,
}

impl Default for MercuryDiagnostics {
    fn default() -> Self {
        Self {
            elapsed_years: 0.0,
            subsolar_longitude_deg: 0.0,
            solar_flux_w_m2: 0.0,
            daytime: false,
            surface_temp_k: 0.0,
            thermal_peak_wavelength_m: 0.0,
            thermal_emission: [0.0; 3],
            solar_tide_bulge_m: 0.0,
            tidal_dissipation_w: 0.0,
            solar_wind_pressure_pa: 0.0,
            sodium_tail_m: 0.0,
            cme_enhancement: 0.0,
            gardening_depth_m: 0.0,
            darkening_factor: 0.0,
            sputtering_rate_kg_s: 0.0,
            thermal_crack_growth: 0.0,
            volcanic_coverage: 0.0,
            hollows_count: 0,
            total_hollows_depth_m: 0.0,
            polar_ice_mass_kg: 0.0,
            polar_deposit_count: 0,
            psr_temperature_k: 0.0,
            leaf_count: 0,
            orbiters_mean_period_h: 0.0,
            observer_altitude_m: 0.0,
            checksum: 0.0,
        }
    }
}

impl MercuryDiagnostics {
    fn recompute_checksum(&mut self) {
        self.checksum = self.elapsed_years
            + self.subsolar_longitude_deg * 1e-2
            + self.solar_flux_w_m2 * 1e-4
            + if self.daytime { 1.0 } else { 0.0 }
            + self.surface_temp_k * 1e-2
            + self.thermal_peak_wavelength_m * 1e6
            + self.thermal_emission.iter().map(|v| *v as f64).sum::<f64>()
            + self.solar_tide_bulge_m * 1e8
            + self.tidal_dissipation_w * 1e-10
            + self.solar_wind_pressure_pa * 1e10
            + self.sodium_tail_m * 1e-8
            + self.cme_enhancement
            + self.gardening_depth_m
            + self.darkening_factor
            + self.sputtering_rate_kg_s * 1e-8
            + self.thermal_crack_growth * 1e12
            + self.volcanic_coverage
            + self.hollows_count as f64 * 1e-3
            + self.total_hollows_depth_m * 1e-2
            + self.polar_ice_mass_kg * 1e-15
            + self.polar_deposit_count as f64 * 1e-2
            + self.psr_temperature_k * 1e-2
            + self.leaf_count as f64 * 1e-3
            + self.orbiters_mean_period_h * 1e-2
            + self.observer_altitude_m * 1e-3;
    }
}

pub struct MercurySimulation {
    orbit: physics::orbit::MercuryOrbit,
    rotation: physics::rotation::MercuryRotation,
    thermal: atmosphere::climate::MercuryThermalState,
    time_scale: temporal::time_scale::TimeScale,
    elapsed_s: f64,
    day_night: lighting::day_night::DayNightCycle,
    observer: geodata::coordinates::LatLon,
    heightmap: terrain::heightmap::Heightmap,
    lod: terrain::lod::LodTerrain,
    weathering: geology::erosion::SpaceWeathering,
    regolith_material: rendering::materials::PbrMaterial,
    exosphere: rendering::atmosphere_scattering::ExosphereEndpoint,
    surface_shader: rendering::shaders::ShaderEndpoint,
    exosphere_shader: rendering::shaders::ShaderEndpoint,
    missions: Vec<satellites::artificial::ArtificialSatellite>,
    diagnostics: MercuryDiagnostics,
}

impl MercurySimulation {
    pub fn new() -> Self {
        let mut heightmap = terrain::heightmap::Heightmap::new(256);
        heightmap.generate_mercury_topography();

        Self {
            orbit: physics::orbit::MercuryOrbit::new(),
            rotation: physics::rotation::MercuryRotation::new(),
            thermal: atmosphere::climate::MercuryThermalState::at_distance(
                SEMI_MAJOR_AXIS_M / sciforge::hub::domain::common::constants::AU,
            ),
            time_scale: temporal::time_scale::TimeScale::new(SOLAR_DAY_S / 240.0),
            elapsed_s: 0.0,
            day_night: lighting::day_night::DayNightCycle::new(),
            observer: geodata::coordinates::LatLon::new(85.0, 30.0),
            heightmap,
            lod: terrain::lod::LodTerrain::new(terrain::lod::LodConfig::default_mercury()),
            weathering: geology::erosion::SpaceWeathering::new(),
            regolith_material: rendering::materials::PbrMaterial::default_regolith(),
            exosphere: rendering::atmosphere_scattering::ExosphereEndpoint::mercury(),
            surface_shader: rendering::shaders::ShaderEndpoint::surface(),
            exosphere_shader: rendering::shaders::ShaderEndpoint::exosphere(),
            missions: satellites::artificial::all_missions(),
            diagnostics: MercuryDiagnostics::default(),
        }
    }

    pub fn tick(&mut self, dt: f64) {
        let scaled_dt = self.time_scale.apply(dt);
        self.elapsed_s += scaled_dt;
        self.orbit.step(scaled_dt);
        self.thermal = atmosphere::climate::MercuryThermalState::at_distance(
            self.orbit.current_radius() / sciforge::hub::domain::common::constants::AU,
        );

        let observer_alt = self
            .heightmap
            .sample(self.observer.lat_deg, self.observer.lon_deg);
        let observer_ecef = self.observer.to_ecef(observer_alt + 1000.0);
        let subsolar_lon = self
            .rotation
            .subsolar_longitude_deg(self.orbit.true_anomaly_rad, self.elapsed_s);
        let sun = lighting::solar_position::SolarPosition::compute(
            self.observer.lat_deg,
            self.observer.lon_deg,
            self.orbit.true_anomaly_rad,
            subsolar_lon,
        );
        let zenith_deg = 90.0 - sun.elevation_deg();
        let surface_temp = self.thermal.temperature_at_zenith_k(zenith_deg.max(0.0));
        let thermal_emission = self.regolith_material.thermal_emission(surface_temp);
        let sun_distance_au =
            self.orbit.current_radius() / sciforge::hub::domain::common::constants::AU;
        let na_glow = self.exosphere.sodium_glow_intensity(sun_distance_au);
        self.surface_shader.uniforms[4].value =
            rendering::shaders::UniformValue::Float(surface_temp as f32);
        self.exosphere_shader.uniforms[0].value =
            rendering::shaders::UniformValue::Float(na_glow as f32);
        self.exosphere_shader.uniforms[2].value =
            rendering::shaders::UniformValue::Float(sun_distance_au as f32);
        let solar_tide = physics::tides::SolarTide::at_distance(self.orbit.current_radius());
        let cme = atmosphere::storms::CmeImpact::typical();
        let hollows = geology::volcanism::known_hollow_fields();
        let deposits = hydrology::glaciers::known_ice_deposits();

        self.diagnostics = MercuryDiagnostics {
            elapsed_years: self.elapsed_s / ORBITAL_PERIOD_S,
            subsolar_longitude_deg: subsolar_lon,
            solar_flux_w_m2: sun.irradiance(),
            daytime: self
                .day_night
                .is_daytime(self.observer.lon_deg, subsolar_lon),
            surface_temp_k: surface_temp,
            thermal_peak_wavelength_m: atmosphere::climate::wien_peak_wavelength(surface_temp),
            thermal_emission,
            solar_tide_bulge_m: solar_tide.tidal_bulge_height(),
            tidal_dissipation_w: physics::tides::tidal_dissipation_rate(
                self.orbit.current_radius(),
            ),
            solar_wind_pressure_pa: atmosphere::winds::solar_wind_dynamic_pressure(),
            sodium_tail_m: atmosphere::winds::sodium_tail_length_estimate_m(),
            cme_enhancement: cme.exosphere_enhancement(),
            gardening_depth_m: self
                .weathering
                .regolith_gardening_depth_m(self.elapsed_s / SECONDS_PER_YEAR),
            darkening_factor: self
                .weathering
                .surface_darkening_factor(self.elapsed_s / SECONDS_PER_YEAR / 1.0e6),
            sputtering_rate_kg_s: self.weathering.global_sputtering_rate_kg_s(),
            thermal_crack_growth: geology::erosion::thermal_fatigue_crack_growth_rate(
                geology::erosion::equatorial_temp_swing(),
                1.0,
            ),
            volcanic_coverage: geology::volcanism::MercuryVolcanism::new()
                .total_volcanic_coverage(),
            hollows_count: hollows.len(),
            total_hollows_depth_m: hollows.iter().map(|h| h.depth_m).sum::<f64>(),
            polar_ice_mass_kg: hydrology::glaciers::total_polar_ice_mass(),
            polar_deposit_count: deposits.len(),
            psr_temperature_k: hydrology::glaciers::psr_equilibrium_temp(self.observer.lat_deg),
            leaf_count: self.lod.total_leaves(),
            orbiters_mean_period_h: self
                .missions
                .iter()
                .map(|m| m.orbital_period_s())
                .sum::<f64>()
                / self.missions.len() as f64
                / 3600.0,
            observer_altitude_m: observer_ecef.altitude(),
            checksum: 0.0,
        };
        self.diagnostics.recompute_checksum();
    }
}

impl Default for MercurySimulation {
    fn default() -> Self {
        Self::new()
    }
}

fn main() {
    let mut sim = MercurySimulation::new();

    for _ in 0..STEP_COUNT {
        sim.tick(DEFAULT_DT);
    }
}