///|
/// Explain the numerical implications of a simulation configuration.
pub(all) struct ConfigAudit {
  valid : Bool
  model_name : String
  viscosity : Double
  recommended_dt : Double
  relaxation_margin : Double
  max_mach : Double
} derive(Debug)

///|
/// Return a stable model name.
pub fn collision_model_name(model : CollisionModel) -> String {
  match model {
    Bgk(_) => "bgk"
    Regularized(_) => "regularized"
    Mrt(_) => "mrt"
  }
}

///|
/// Audit options before a long numerical run.
pub fn audit_options(options : SimulationOptions) -> ConfigAudit {
  let omega = options.model.omega()
  let margin = omega.min(2.0 - omega)
  let valid = omega > 0.0 &&
    omega < 2.0 &&
    options.max_mach > 0.0 &&
    options.density_floor > 0.0
  {
    valid,
    model_name: collision_model_name(options.model),
    viscosity: kinematic_viscosity(options.model),
    recommended_dt: if valid {
      1.0
    } else {
      0.0
    },
    relaxation_margin: margin.max(0.0),
    max_mach: options.max_mach,
  }
}

///|
/// Return a human-readable configuration audit.
pub fn ConfigAudit::to_string(self : ConfigAudit) -> String {
  "model=\{self.model_name}\nvalid=\{self.valid}\nviscosity=\{self.viscosity}\nrecommended_dt=\{self.recommended_dt}\nrelaxation_margin=\{self.relaxation_margin}\nmax_mach=\{self.max_mach}\n"
}

///|
/// Return the relaxation range that passes the ordinary BGK stability criterion.
pub fn stable_relaxation_range() -> (Double, Double) {
  (0.0000001, 1.9999999)
}

///|
/// Return whether a relaxation value is in the open stable interval.
pub fn relaxation_is_stable(omega : Double) -> Bool {
  omega > 0.0 && omega < 2.0
}

///|
/// Estimate a relaxation parameter from a target viscosity.
pub fn relaxation_from_viscosity(viscosity : Double) -> Double {
  if viscosity <= 0.0 {
    1.0
  } else {
    1.0 / (3.0 * viscosity + 0.5)
  }
}

///|
/// Estimate the BGK viscosity for an omega value.
pub fn viscosity_from_relaxation(omega : Double) -> Double {
  if relaxation_is_stable(omega) {
    (1.0 / omega - 0.5) / 3.0
  } else {
    0.0
  }
}

///|
/// Return a conservative inlet speed bound for a Mach limit.
pub fn speed_limit_from_mach(max_mach~ : Double) -> Double {
  max_mach.max(0.0) * lattice_sound_speed()
}

///|
/// Return a stable force magnitude bound.
pub fn force_limit_from_density(
  density~ : Double,
  max_mach~ : Double,
) -> Double {
  density.max(0.0) * speed_limit_from_mach(max_mach~) * 0.1
}