///|
/// Auditable property correlations used for preliminary process screening.
pub(all) struct AntoineCorrelation {
  component : String
  a : Double
  b : Double
  c : Double
  pressure_unit : String
} derive(Debug, Eq, ToJson)

///|
pub(all) struct DensityCorrelation {
  component : String
  reference_temperature : Double
  reference_density : Double
  expansion_coefficient : Double
} derive(Debug, Eq, ToJson)

///|
pub(all) struct ViscosityCorrelation {
  component : String
  reference_temperature : Double
  reference_viscosity : Double
  temperature_sensitivity : Double
} derive(Debug, Eq, ToJson)

///|
pub(all) struct PropertyPoint {
  temperature : Double
  pressure : Double
  vapor_pressure : Double
  density : Double
  viscosity : Double
  surface_tension : Double
} derive(Debug, Eq, ToJson)

///|
pub fn AntoineCorrelation::vapor_pressure(
  self : AntoineCorrelation,
  temperature : Double,
) -> Double raise @domain.ExtractionError {
  if temperature + self.c <= 0.0 {
    raise @domain.ExtractionError::InvalidArgument(
      message="Antoine temperature is invalid",
    )
  }
  let exponent = self.a - self.b / (temperature + self.c)
  approximate_ten_power(exponent)
}

///|
pub fn approximate_ten_power(exponent : Double) -> Double {
  let bounded = exponent.max(-6.0).min(6.0)
  let integer = bounded.to_int()
  let fraction = bounded - integer.to_double()
  let mut integer_power = 1.0
  if integer >= 0 {
    for _ in 0.. Double raise @domain.ExtractionError {
  if correlation.reference_density <= 0.0 ||
    correlation.expansion_coefficient < 0.0 ||
    pressure_factor <= 0.0 {
    raise @domain.ExtractionError::InvalidArgument(
      message="density correlation is invalid",
    )
  }
  let thermal = 1.0 -
    correlation.expansion_coefficient *
    (temperature - correlation.reference_temperature)
  if thermal <= 0.0 {
    raise @domain.ExtractionError::InvalidArgument(
      message="density correlation extrapolates to non-positive density",
    )
  }
  correlation.reference_density * thermal * pressure_factor
}

///|
pub fn viscosity_at(
  correlation : ViscosityCorrelation,
  temperature : Double,
) -> Double raise @domain.ExtractionError {
  if correlation.reference_viscosity <= 0.0 || temperature <= 0.0 {
    raise @domain.ExtractionError::InvalidArgument(
      message="viscosity correlation is invalid",
    )
  }
  let factor = 1.0 +
    correlation.temperature_sensitivity *
    (correlation.reference_temperature / temperature - 1.0)
  correlation.reference_viscosity * factor.max(0.05)
}

///|
pub fn surface_tension_at(
  reference_surface_tension : Double,
  reference_temperature : Double,
  temperature : Double,
  critical_temperature : Double,
) -> Double raise @domain.ExtractionError {
  if reference_surface_tension < 0.0 ||
    reference_temperature <= 0.0 ||
    temperature <= 0.0 ||
    critical_temperature <= temperature {
    raise @domain.ExtractionError::InvalidArgument(
      message="surface tension inputs are invalid",
    )
  }
  let reduced = ((critical_temperature - temperature) /
  (critical_temperature - reference_temperature)).max(0.0)
  reference_surface_tension * reduced.sqrt() * reduced.sqrt().sqrt()
}

///|
pub fn property_point(
  temperature : Double,
  pressure : Double,
  vapor_pressure_model : AntoineCorrelation,
  density_model : DensityCorrelation,
  viscosity_model : ViscosityCorrelation,
  reference_surface_tension : Double,
  critical_temperature : Double,
) -> PropertyPoint raise @domain.ExtractionError {
  {
    temperature,
    pressure,
    vapor_pressure: vapor_pressure_model.vapor_pressure(temperature),
    density: density_at(density_model, temperature, pressure),
    viscosity: viscosity_at(viscosity_model, temperature),
    surface_tension: surface_tension_at(
      reference_surface_tension,
      density_model.reference_temperature,
      temperature,
      critical_temperature,
    ),
  }
}

///|
pub fn property_grid(
  temperatures : ArrayView[Double],
  pressures : ArrayView[Double],
  vapor_pressure_model : AntoineCorrelation,
  density_model : DensityCorrelation,
  viscosity_model : ViscosityCorrelation,
  reference_surface_tension : Double,
  critical_temperature : Double,
) -> Array[PropertyPoint] raise @domain.ExtractionError {
  let points = []
  for temperature in temperatures {
    for pressure in pressures {
      points.push(
        property_point(
          temperature, pressure, vapor_pressure_model, density_model, viscosity_model,
          reference_surface_tension, critical_temperature,
        ),
      )
    }
  }
  points
}

///|
pub fn property_average(points : ArrayView[PropertyPoint]) -> PropertyPoint {
  if points.length() == 0 {
    return {
      temperature: 0.0,
      pressure: 0.0,
      vapor_pressure: 0.0,
      density: 0.0,
      viscosity: 0.0,
      surface_tension: 0.0,
    }
  }
  let mut temperature = 0.0
  let mut pressure = 0.0
  let mut vapor = 0.0
  let mut density = 0.0
  let mut viscosity = 0.0
  let mut surface = 0.0
  for point in points {
    temperature += point.temperature
    pressure += point.pressure
    vapor += point.vapor_pressure
    density += point.density
    viscosity += point.viscosity
    surface += point.surface_tension
  }
  let count = points.length().to_double()
  {
    temperature: temperature / count,
    pressure: pressure / count,
    vapor_pressure: vapor / count,
    density: density / count,
    viscosity: viscosity / count,
    surface_tension: surface / count,
  }
}

///|
pub fn property_minimum(points : ArrayView[PropertyPoint]) -> PropertyPoint {
  if points.length() == 0 {
    return {
      temperature: 0.0,
      pressure: 0.0,
      vapor_pressure: 0.0,
      density: 0.0,
      viscosity: 0.0,
      surface_tension: 0.0,
    }
  }
  let mut result = points[0]
  for point in points {
    result = {
      temperature: result.temperature.min(point.temperature),
      pressure: result.pressure.min(point.pressure),
      vapor_pressure: result.vapor_pressure.min(point.vapor_pressure),
      density: result.density.min(point.density),
      viscosity: result.viscosity.min(point.viscosity),
      surface_tension: result.surface_tension.min(point.surface_tension),
    }
  }
  result
}

///|
pub fn property_maximum(points : ArrayView[PropertyPoint]) -> PropertyPoint {
  if points.length() == 0 {
    return {
      temperature: 0.0,
      pressure: 0.0,
      vapor_pressure: 0.0,
      density: 0.0,
      viscosity: 0.0,
      surface_tension: 0.0,
    }
  }
  let mut result = points[0]
  for point in points {
    result = {
      temperature: result.temperature.max(point.temperature),
      pressure: result.pressure.max(point.pressure),
      vapor_pressure: result.vapor_pressure.max(point.vapor_pressure),
      density: result.density.max(point.density),
      viscosity: result.viscosity.max(point.viscosity),
      surface_tension: result.surface_tension.max(point.surface_tension),
    }
  }
  result
}

///|
pub fn reynolds_number(
  density : Double,
  velocity : Double,
  diameter : Double,
  viscosity : Double,
) -> Double raise @domain.ExtractionError {
  if density <= 0.0 || velocity < 0.0 || diameter <= 0.0 || viscosity <= 0.0 {
    raise @domain.ExtractionError::InvalidArgument(
      message="Reynolds inputs are invalid",
    )
  }
  density * velocity * diameter / viscosity
}

///|
pub fn capillary_number(
  viscosity : Double,
  velocity : Double,
  surface_tension : Double,
) -> Double raise @domain.ExtractionError {
  if viscosity < 0.0 || velocity < 0.0 || surface_tension <= 0.0 {
    raise @domain.ExtractionError::InvalidArgument(
      message="capillary inputs are invalid",
    )
  }
  viscosity * velocity / surface_tension
}

///|
pub fn property_csv(points : ArrayView[PropertyPoint]) -> String {
  let lines = [
    "temperature,pressure,vapor_pressure,density,viscosity,surface_tension",
  ]
  for point in points {
    lines.push(
      point.temperature.to_string() +
      "," +
      point.pressure.to_string() +
      "," +
      point.vapor_pressure.to_string() +
      "," +
      point.density.to_string() +
      "," +
      point.viscosity.to_string() +
      "," +
      point.surface_tension.to_string(),
    )
  }
  lines.join("\n")
}

///|
pub fn property_report(points : ArrayView[PropertyPoint]) -> String {
  let average = property_average(points)
  let minimum = property_minimum(points)
  let maximum = property_maximum(points)
  "samples=" +
  points.length().to_string() +
  ",average_density=" +
  average.density.to_string() +
  ",density_range=[" +
  minimum.density.to_string() +
  "," +
  maximum.density.to_string() +
  "]" +
  ",average_viscosity=" +
  average.viscosity.to_string()
}