///|
/// Amount units used for molar balances.
pub(all) enum AmountUnit {
  Mole
  KiloMole
  PoundMole
} derive(Debug, Eq, ToJson)

///|
pub(all) struct Amount {
  value : Double
  unit : AmountUnit
} derive(Debug, Eq, ToJson)

///|
pub fn Amount::new(
  value~ : Double,
  unit~ : AmountUnit,
) -> Amount raise @domain.ExtractionError {
  @domain.require_non_negative("amount", value)
  { value, unit }
}

///|
pub fn Amount::mole(self : Amount) -> Double {
  match self.unit {
    Mole => self.value
    KiloMole => self.value * 1000.0
    PoundMole => self.value * 453.59237
  }
}

///|
pub fn Amount::kilomole(self : Amount) -> Double {
  self.mole() / 1000.0
}

///|
pub fn Amount::pound_mole(self : Amount) -> Double {
  self.mole() / 453.59237
}

///|
pub fn Amount::convert(self : Amount, unit : AmountUnit) -> Amount {
  let value = match unit {
    Mole => self.mole()
    KiloMole => self.kilomole()
    PoundMole => self.pound_mole()
  }
  { value, unit }
}

///|
pub fn Amount::add(self : Amount, other : Amount) -> Amount {
  { value: self.value + other.convert(self.unit).value, unit: self.unit }
}

///|
pub(all) struct MolarFlow {
  amount : Amount
  duration : Duration
} derive(Debug, Eq, ToJson)

///|
pub fn MolarFlow::new(amount~ : Amount, duration~ : Duration) -> MolarFlow {
  { amount, duration }
}

///|
pub fn MolarFlow::mole_per_second(self : MolarFlow) -> Double {
  self.amount.mole() / self.duration.per_second()
}

///|
pub fn MolarFlow::kilomole_per_hour(self : MolarFlow) -> Double {
  self.mole_per_second() * 3.6
}

///|
pub fn MolarFlow::scale(
  self : MolarFlow,
  factor : Double,
) -> MolarFlow raise @domain.ExtractionError {
  @domain.require_non_negative("molar flow scale", factor)
  {
    amount: Amount::{
      value: self.amount.value * factor,
      unit: self.amount.unit,
    },
    duration: self.duration,
  }
}

///|
pub(all) struct Density {
  mass_kg_per_cubic_meter : Double
} derive(Debug, Eq, ToJson)

///|
pub fn Density::new(
  mass_kg_per_cubic_meter~ : Double,
) -> Density raise @domain.ExtractionError {
  @domain.require_positive("density", mass_kg_per_cubic_meter)
  { mass_kg_per_cubic_meter, }
}

///|
pub fn Density::mass(
  self : Density,
  volume : Volume,
) -> Mass raise @domain.ExtractionError {
  Mass::new(
    value=self.mass_kg_per_cubic_meter * volume.cubic_meter(),
    unit=Kilogram,
  )
}

///|
pub fn Density::volume(
  self : Density,
  mass : Mass,
) -> Volume raise @domain.ExtractionError {
  Volume::new(
    value=mass.convert(Kilogram).value() / self.mass_kg_per_cubic_meter,
    unit=CubicMeter,
  )
}

///|
pub(all) struct DynamicViscosity {
  pascal_second : Double
} derive(Debug, Eq, ToJson)

///|
pub fn DynamicViscosity::new(
  pascal_second~ : Double,
) -> DynamicViscosity raise @domain.ExtractionError {
  @domain.require_positive("dynamic viscosity", pascal_second)
  { pascal_second, }
}

///|
pub fn DynamicViscosity::centipoise(self : DynamicViscosity) -> Double {
  self.pascal_second * 1000.0
}

///|
pub fn DynamicViscosity::kinematic(
  self : DynamicViscosity,
  density : Density,
) -> Double {
  self.pascal_second / density.mass_kg_per_cubic_meter
}

///|
pub(all) struct UnitInterval {
  lower : Double
  upper : Double
} derive(Debug, Eq, ToJson)

///|
pub fn UnitInterval::new(
  lower~ : Double,
  upper~ : Double,
) -> UnitInterval raise @domain.ExtractionError {
  @domain.require_fraction("interval lower", lower)
  @domain.require_fraction("interval upper", upper)
  if upper < lower {
    raise @domain.ExtractionError::InvalidArgument(
      message="unit interval upper bound is below lower bound",
    )
  }
  { lower, upper }
}

///|
pub fn UnitInterval::width(self : UnitInterval) -> Double {
  self.upper - self.lower
}

///|
pub fn UnitInterval::midpoint(self : UnitInterval) -> Double {
  (self.lower + self.upper) / 2.0
}

///|
pub fn UnitInterval::contains(self : UnitInterval, value : Double) -> Bool {
  value >= self.lower && value <= self.upper
}

///|
pub fn UnitInterval::map(
  self : UnitInterval,
  value : Double,
  target : UnitInterval,
) -> Double {
  target.lower + self.normalized(value) * target.width()
}

///|
pub fn UnitInterval::normalized(self : UnitInterval, value : Double) -> Double {
  if self.width() <= 0.0 {
    0.0
  } else {
    (value - self.lower) / self.width()
  }
}

///|
pub fn volume_flow_from_mass(
  mass_flow~ : MassFlow,
  density~ : Density,
) -> Volume {
  Volume::{
    value: mass_flow.kg_per_hour() / density.mass_kg_per_cubic_meter,
    unit: CubicMeter,
  }
}

///|
pub fn mass_flow_from_volume(
  volume~ : Volume,
  density~ : Density,
  duration~ : Duration,
) -> MassFlow raise @domain.ExtractionError {
  let mass = density.mass(volume)
  MassFlow::new(mass~, duration~)
}

///|
pub fn reynolds_number(
  density~ : Density,
  velocity_m_per_second~ : Double,
  diameter_meter~ : Double,
  viscosity~ : DynamicViscosity,
) -> Double raise @domain.ExtractionError {
  @domain.require_non_negative("velocity", velocity_m_per_second)
  @domain.require_positive("diameter", diameter_meter)
  density.mass_kg_per_cubic_meter *
  velocity_m_per_second *
  diameter_meter /
  viscosity.pascal_second
}

///|
pub fn schmidt_number(
  viscosity~ : DynamicViscosity,
  density~ : Density,
  diffusivity_square_meter_per_second~ : Double,
) -> Double raise @domain.ExtractionError {
  @domain.require_positive("diffusivity", diffusivity_square_meter_per_second)
  viscosity.pascal_second /
  density.mass_kg_per_cubic_meter /
  diffusivity_square_meter_per_second
}

///|
pub fn peclet_number(
  velocity_m_per_second~ : Double,
  length_meter~ : Double,
  diffusivity_square_meter_per_second~ : Double,
) -> Double raise @domain.ExtractionError {
  @domain.require_non_negative("velocity", velocity_m_per_second)
  @domain.require_positive("length", length_meter)
  @domain.require_positive("diffusivity", diffusivity_square_meter_per_second)
  velocity_m_per_second * length_meter / diffusivity_square_meter_per_second
}

///|
pub fn weber_number(
  density~ : Density,
  velocity_m_per_second~ : Double,
  length_meter~ : Double,
  surface_tension_newton_per_meter~ : Double,
) -> Double raise @domain.ExtractionError {
  @domain.require_non_negative("velocity", velocity_m_per_second)
  @domain.require_positive("length", length_meter)
  @domain.require_positive("surface tension", surface_tension_newton_per_meter)
  density.mass_kg_per_cubic_meter *
  velocity_m_per_second *
  velocity_m_per_second *
  length_meter /
  surface_tension_newton_per_meter
}

///|
pub fn froude_number(
  velocity_m_per_second~ : Double,
  length_meter~ : Double,
  gravity_m_per_second_squared? : Double = 9.80665,
) -> Double raise @domain.ExtractionError {
  @domain.require_non_negative("velocity", velocity_m_per_second)
  @domain.require_positive("length", length_meter)
  @domain.require_positive("gravity", gravity_m_per_second_squared)
  velocity_m_per_second / (gravity_m_per_second_squared * length_meter).sqrt()
}

///|
pub(all) struct DimensionlessGroups {
  reynolds : Double
  schmidt : Double
  peclet : Double
  weber : Double
  froude : Double
} derive(Debug, Eq, ToJson)

///|
pub fn calculate_dimensionless_groups(
  density~ : Density,
  velocity_m_per_second~ : Double,
  diameter_meter~ : Double,
  viscosity~ : DynamicViscosity,
  diffusivity_square_meter_per_second~ : Double,
  surface_tension_newton_per_meter~ : Double,
) -> DimensionlessGroups raise @domain.ExtractionError {
  {
    reynolds: reynolds_number(
      density~,
      velocity_m_per_second~,
      diameter_meter~,
      viscosity~,
    ),
    schmidt: schmidt_number(
      viscosity~,
      density~,
      diffusivity_square_meter_per_second~,
    ),
    peclet: peclet_number(
      velocity_m_per_second~,
      length_meter=diameter_meter,
      diffusivity_square_meter_per_second~,
    ),
    weber: weber_number(
      density~,
      velocity_m_per_second~,
      length_meter=diameter_meter,
      surface_tension_newton_per_meter~,
    ),
    froude: froude_number(velocity_m_per_second~, length_meter=diameter_meter),
  }
}

///|
pub fn mixing_power_number(
  power_watt~ : Double,
  density~ : Density,
  rotational_speed_per_second~ : Double,
  impeller_diameter_meter~ : Double,
) -> Double raise @domain.ExtractionError {
  @domain.require_non_negative("mixing power", power_watt)
  @domain.require_positive("rotational speed", rotational_speed_per_second)
  @domain.require_positive("impeller diameter", impeller_diameter_meter)
  power_watt /
  (
    density.mass_kg_per_cubic_meter *
    rotational_speed_per_second *
    rotational_speed_per_second *
    rotational_speed_per_second *
    impeller_diameter_meter *
    impeller_diameter_meter *
    impeller_diameter_meter *
    impeller_diameter_meter *
    impeller_diameter_meter
  )
}

///|
pub fn heat_transfer_coefficient(
  nusselt~ : Double,
  thermal_conductivity_w_per_meter_k~ : Double,
  characteristic_length_meter~ : Double,
) -> Double raise @domain.ExtractionError {
  @domain.require_non_negative("Nusselt number", nusselt)
  @domain.require_positive(
    "thermal conductivity", thermal_conductivity_w_per_meter_k,
  )
  @domain.require_positive("characteristic length", characteristic_length_meter)
  nusselt * thermal_conductivity_w_per_meter_k / characteristic_length_meter
}

///|
pub fn nusselt_laminar_pipe(
  reynolds : Double,
  prandtl : Double,
) -> Double raise @domain.ExtractionError {
  @domain.require_non_negative("Reynolds number", reynolds)
  @domain.require_positive("Prandtl number", prandtl)
  if reynolds <= 2300.0 {
    3.66 +
    0.0668 *
    reynolds *
    prandtl /
    (1.0 + 0.04 * @math.pow(reynolds * prandtl, 2.0 / 3.0))
  } else {
    0.023 * @math.pow(reynolds, 0.8) * @math.pow(prandtl, 0.4)
  }
}

///|
pub fn pressure_drop_pipe(
  friction_factor~ : Double,
  length_meter~ : Double,
  diameter_meter~ : Double,
  density~ : Density,
  velocity_m_per_second~ : Double,
) -> Pressure raise @domain.ExtractionError {
  @domain.require_non_negative("friction factor", friction_factor)
  @domain.require_positive("pipe length", length_meter)
  @domain.require_positive("pipe diameter", diameter_meter)
  @domain.require_non_negative("velocity", velocity_m_per_second)
  let pascal = friction_factor *
    length_meter /
    diameter_meter *
    density.mass_kg_per_cubic_meter *
    velocity_m_per_second *
    velocity_m_per_second /
    2.0
  Pressure::new(value=pascal, unit=PressureUnit::Pascal)
}

///|
pub fn friction_factor_blasius(
  reynolds : Double,
) -> Double raise @domain.ExtractionError {
  @domain.require_positive("Reynolds number", reynolds)
  if reynolds < 2300.0 {
    64.0 / reynolds
  } else {
    0.3164 / @math.pow(reynolds, 0.25)
  }
}

///|
pub fn residence_time_from_geometry(
  vessel_volume~ : Volume,
  carrier_flow~ : Volume,
  solvent_flow~ : Volume,
  duration~ : Duration,
) -> Duration raise @domain.ExtractionError {
  residence_time_from_flow(
    vessel_volume~,
    volumetric_flow=carrier_flow.add(solvent_flow),
    duration~,
  )
}