use sciforge::hub::domain::common::constants::{R_GAS, SIGMA_SB};
use sciforge::hub::domain::physics::thermodynamics::clausius_clapeyron;
use sciforge::hub::prelude::constants::elements::atomic_mass;

pub const INTERIOROCEANDEPTHM: f64 = 5000000.0;
pub const IONICOCEANDENSITY: f64 = 1400.0;
pub const MEANSALINITYPSU: f64 = 0.0;

pub struct InteriorOceanLayer {
    pub name: &'static str,
    pub depthrangem: (f64, f64),
    pub meantemperaturek: f64,
    pub compositionh2ofraction: f64,
}

impl InteriorOceanLayer {
    pub fn densitykgm3(&self) -> f64 {
        let t = self.meantemperaturek;
        let basedensity = IONICOCEANDENSITY + 0.1 * (2000.0 - t).max(0.0);
        basedensity * (0.9 + 0.1 * self.compositionh2ofraction)
    }

    pub fn pressureatdepth(&self, depthm: f64) -> f64 {
        self.densitykgm3() * *crate::SURFACEGRAVITY * depthm
    }

    pub fn soundspeedms(&self) -> f64 {
        let t = self.meantemperaturek;
        let d = (self.depthrangem.0 + self.depthrangem.1) / 2.0;
        1500.0 + 2.0 * (t - 273.0) + 0.016 * d * 1e-3
    }
}

pub fn superioniclayer() -> InteriorOceanLayer {
    InteriorOceanLayer {
        name: "Superionic Water Layer",
        depthrangem: (7000000.0, 15000000.0),
        meantemperaturek: 5000.0,
        compositionh2ofraction: 0.7,
    }
}

pub fn ioniclayer() -> InteriorOceanLayer {
    InteriorOceanLayer {
        name: "Ionic Water Layer",
        depthrangem: (5000000.0, 7000000.0),
        meantemperaturek: 3000.0,
        compositionh2ofraction: 0.8,
    }
}

pub fn molecularlayer() -> InteriorOceanLayer {
    InteriorOceanLayer {
        name: "Molecular Fluid Layer",
        depthrangem: (100000.0, 5000000.0),
        meantemperaturek: 1000.0,
        compositionh2ofraction: 0.6,
    }
}

pub fn interiorheatcontentj(tempanomalyk: f64, depthm: f64) -> f64 {
    let surfacearea =
        4.0 * std::f64::consts::PI * crate::NEPTUNEMEANRADIUS * crate::NEPTUNEMEANRADIUS;
    IONICOCEANDENSITY * 4200.0 * surfacearea * depthm * tempanomalyk
}

pub fn surfacelongwaveradiationwm2(tempk: f64, emissivity: f64) -> f64 {
    emissivity * SIGMA_SB * tempk.powi(4)
}

pub fn vaporpressureovermethane(tempk: f64) -> f64 {
    let pref = 101325.0;
    let tref = 111.7;
    let lvap = 510000.0;
    clausius_clapeyron(pref, tref, tempk, lvap)
}

pub fn dissolvedgashenry(partialpressurepa: f64, temperaturek: f64, henryconstmolm3pa: f64) -> f64 {
    let correction = (-2400.0 * (1.0 / temperaturek - 1.0 / 298.15) / R_GAS).exp();
    henryconstmolm3pa * partialpressurepa * correction
}

pub fn h2omolarmass() -> f64 {
    (2.0 * atomic_mass(1) + atomic_mass(8)) * 1e-3
}
