// Two deterministic experiments using packed unit-bearing arrays.
// Requires the updated native Aner development runtime; no packaged release is implied.
// All inputs are synthetic. Population spread is not a measurement uncertainty model.
import aner.units
import aner.math
import aner.tensor

fn close(actual: Float64, expected: Float64, tolerance: Float64 = 0.000000001) -> Bool {
    return math.abs(actual - expected) <= tolerance
}
fn verify(label: String, condition: Bool) -> Unit {
    print(label)
    print(condition)
}

// Temperature experiment: normalize Fahrenheit readings, then add a declared offset.
let source = DynamicArray<Point<degF>>(capacity: 4)
source.append(68.0@degF)
source.append(71.6@degF)
source.append(75.2@degF)
source.append(78.8@degF)
let celsius: DynamicArray<Point<degC>> = source.to("degC")
let correction = DynamicArray<Quantity<delta_degC>>(capacity: 4)
correction.resize(size: 4, fill: 0.5@delta_degC)
let calibrated: DynamicArray<Point<degC>> = celsius.add(correction)
verify("Four readings normalized and calibrated", calibrated.len() == 4)
verify("Mean is 23.5 degC", close(calibrated.mean().value(), 23.5))
verify("Calibration retains Point<degC>", type(calibrated.mean()) == "Point<degC>")
verify("Population spread is sqrt(5) delta_degC", close(calibrated.std().value(), math.sqrt(5.0)))
verify("Spread retains an interval unit", calibrated.std().unit() == "delta_degC")
verify("Constant correction preserves observed spread", close(celsius.std().value(), calibrated.std().value()))
let kelvin = calibrated.to("K")
verify("Mean is 296.65 K", close(kelvin.mean().value(), 296.65))
let reference = DynamicArray<Point<degC>>(capacity: 4)
reference.resize(size: 4, fill: 23.5@degC)
let anomalies: DynamicArray<Quantity<delta_degC>> = calibrated.sub(reference)
verify("Anomalies sum to zero", close(anomalies.sum().value(), 0.0))
verify("First anomaly is -3 delta_degC", close(anomalies.get(0).value(), -3.0))

// Travel experiment: rows identify two devices, columns three runs in fixed order.
let travel = NDArray<Quantity<cm>>(shape: List(2, 3), fill: 0.0@cm)
travel[0, 0] = 100.0@cm
travel[0, 1] = 120.0@cm
travel[0, 2] = 140.0@cm
travel[1, 0] = 80.0@cm
travel[1, 1] = 100.0@cm
travel[1, 2] = 120.0@cm
let elapsed = NDArray<Quantity<s>>(shape: List(1, 3), fill: 0.0@s)
elapsed[0, 0] = 2.0@s
elapsed[0, 1] = 4.0@s
elapsed[0, 2] = 5.0@s
let speeds: NDArray<Quantity<"m/s">> = travel.div(elapsed)
verify("Broadcast result has two rows and three runs", cast<Int64>(speeds.shape().get(0)) == 2 && cast<Int64>(speeds.shape().get(1)) == 3)
verify("Canonical speed unit is m/s", speeds.unit() == "m/s")
verify("First speed is 0.5 m/s", close(speeds[0, 0].value(), 0.5))
verify("Last speed is 0.24 m/s", close(speeds[1, 2].value(), 0.24))
verify("First speed converts to 1.8 km/h", close(speeds.to("km/h")[0, 0].value(), 1.8))

// Dimensionless weights retain a registered product unit in rank-two matmul.
let weights = NDArray<Quantity<"1">>(shape: List(3, 1), fill: units.quantity(1.0 / 3.0, "1"))
let average: NDArray<Quantity<"m/s">> = speeds.matmul(weights)
verify("Device one mean speed is 0.36 m/s", close(average[0, 0].value(), 0.36))
verify("Device two mean speed is 0.2966666666666667 m/s", close(average[1, 0].value(), 0.2966666666666667))

// Export is deliberate and detached. Retain the chosen unit beside raw values.
let export_unit = speeds.unit()
let raw = speeds.values()
let model_input = raw.to_tensor()
verify("Export identifies m/s and Float64 storage", export_unit == "m/s" && dtype(raw) == "Float64")
verify("Tensor export has the same six values", tensor.rows(model_input) == 2 && tensor.cols(model_input) == 3 && close(tensor.value(model_input, 1, 2), 0.24))
raw[0, 0] = -999.0
verify("Raw export does not alias the unit buffer", close(speeds[0, 0].value(), 0.5))
let same_unit_copy = speeds.to("m/s")
same_unit_copy[0, 0] = 99.0@"m/s"
verify("Same-unit conversion is also independent", close(speeds[0, 0].value(), 0.5))