Add try_pass_hessian to the highs API

[GermanHeim]

This PR exposes Highs_passHessian from the C API with the following methods:
try_pass_hessian and pass_hessian.

I've also added some tests. I thought about adding an example, but it doesn't seem we have an examples folder, and I am not sure if adding to the README is fine.

Closes #43

[lovasoa]

I think the rust api can still be improved

[GermanHeim]

I think the suggestions of IntoIterator and the double allocations overlap. I don't think we can fully do both. Since Highs_passHessian takes three separate contiguous pointers (start, index, value), this makes it so that the two suggestions can't be implemented together:

• Pairs iterator: It is clean and matches add_column's &[(Row, N)] that is already in the API, and makes a length mismatch unrepresentable. But to pass it to the C API, we have to unzip the pairs into two Vecs, which adds a copy of value.
• Preallocated typed slices: this is good because it is zero copies, but it keeps the two parallel arrays, and pushes the HighsInt/i32 type onto the caller.

I would do the pairs iterator since the crate is mostly an ergonomic wrapper, and the API being better is worth one value copy, in my opinion.

[lovasoa]

Yes, we need either an iterator OR an array of the right type as input. But an array of the wrong type that has to be converted is the worst of both worlds :p

I think the IntoIterator version is cleaner, even if it may indeed often result in double allocation in many cases. At least performance-minded devs have the option to stream the values from whatever their source is if necessary.

[GermanHeim]

I coded an API like:

pub fn try_pass_hessian<S, E, I>(&mut self, format: HessianFormat, start: S, entries: E)
    -> Result<(), HighsStatus>
where
    S: IntoIterator, S::Item: TryInto<HighsInt>,
    E: IntoIterator<Item = (I, f64)>, I: TryInto<HighsInt>;

However, pass_hessian panic can't use dim/nnz, since iterators consume once. If we want a better panic message, we would need to collect first, which I don't think is worth it. Is this fine? I wanted to check before committing the code.

[lovasoa]

Yes the api looks great !

In the future, don't hesitate to commit and update the pr even when you are not sure. Github allows setting the pr state to "draft" while you are iterating on api ideas.

[lovasoa]

undefined behavior is very scary in rust. does that mean "memory unsafe" (in which case we must either mark the function unsafe in rust, or verify the inputs before calling highs) ? Or does it just mean highs does not give any guarantee as to what the result will be (but guarantees it won't expose uninitialized memory, read out of bounds, or anything nasty) ?

[GermanHeim]

What I meant by this is that HiGHS doesn't detect some ill-defined problems (in my experience).
For example:
min x^2 + 4xy + y^2 s.t. x + y = 1 (variables free).

use highs::{HessianFormat, RowProblem, Sense};

fn main() {
    let mut pb = RowProblem::new();
    let x = pb.add_column(0.0, f64::NEG_INFINITY..=f64::INFINITY);
    let y = pb.add_column(0.0, f64::NEG_INFINITY..=f64::INFINITY);
    pb.add_row(1.0..=1.0, [(x, 1.0), (y, 1.0)]);

    let mut model = pb.optimise(Sense::Minimise);

    // Lower-tri CSC of Q = [[2, 4], [4, 2]]:
    //   col 0: rows (0, 2.0) and (1, 4.0)
    //   col 1: row  (1, 2.0)
    let start = [0, 2];
    let index = [0, 1, 1];
    let value = [2.0, 4.0, 2.0];

    let status = model.try_pass_hessian(
        HessianFormat::Triangular,
        start,
        index.into_iter().zip(value),
    );
    println!("try_pass_hessian = {status:?}");

    let solved = model.solve();
    println!("status    = {:?}", solved.status());
    println!("objective = {}", solved.objective_value());
    println!("primal    = {:?}", solved.get_solution().columns());
}

This returns:

try_pass_hessian = Ok(())
status    = Optimal
objective = 1
primal    = [0.0, 1.0]

But should be -inf.
But a negative diagonal is detected, and HiGHS returns an error status at solve().

Both are memory safe. Maybe we could be clearer in the comments?

[jajhall]

Indeed, HiGHS checks only for negative diagonal values on entry, and it's arguable that it shouldn't, as convexity of a QP is dependent of the feasible region.

[lovasoa]

@jajhall can passing invalid problems result in true undefined behavior, or is it just unspecified behavior (highs returns well defined but arbitrary values) ?

[jajhall]

If HiGHS detects negative curvature, then it stops with an error message, since non-convexity has been identified. I'd call that well-defined behaviour

[lovasoa]

Good. And when highs does not detect non-convexity, like in the example German found above, the results are still well-defined (but incorrect), right ?

@GermanHeim it may be worth updating the note in the function's doc, to make that clear.

[jajhall]

I'll code up the example and see precisely what happens in the QP solver

[jajhall]

So, the HiGHS QP solver had a bug in it, leading to the incorrect assertion of optimality at [0, 1] for

min x^2 + 4xy + y^2 s.t. x + y = 1 (variables free).

It now identifies that the Hessian has negative curvature in the feasible direction at [0, 1], logs that the problem is non-convex, and returns an error.

[GermanHeim]

Amazing work @jajhall! @lovasoa maybe we should update the docs for the bindings? I am not sure if there are any other "undefined behavior" (in the mathematical sense)

[lovasoa]

The API looks good now I think, but two points still bother me:

• start yields the offset into entries : this still makes it possible to represent impossible states. If instead we had an iterator of columns, with each column an iterator of (row index, coefficient), we would make it impossible to represent invalid inputs. Better ?
• the function returns a single bit of information. all these .map_err(|_| HighsStatus::Error)) erase the error and make it unavailable to the caller. and panic!("{e:?}") cannot display anything actually useful. Maybe we should have a custom error enum that allows having actually helpful error messages (like "the dimension of the quadratic objective matrix (hessian Q) is too large: got XXX but expected at most YYY").

[GermanHeim]

Yes, I agree, those are good ideas. I will work on it tomorrow.

[GermanHeim]

Hello, I added the HessianError enum. I also added the columns iterator, and I tried to improve the message about the UB.

[lovasoa]

Published as v2.2.0.

Thank you very much @GermanHeim , this is a great contribution.

[GermanHeim]

Thanks for the quick response and the great reviews @lovasoa!