use nalgebra::{DMatrix, DVector};
use crate::bytecode_tape::{BtapeThreadLocal, BytecodeTape};
use crate::float::Float;
use crate::BReverse;
pub fn grad_nalgebra<F: Float + BtapeThreadLocal>(
f: impl FnOnce(&[BReverse<F>]) -> BReverse<F>,
x: &DVector<F>,
) -> DVector<F> {
let (mut tape, _) = crate::api::record(f, x.as_slice());
tape_gradient_nalgebra(&mut tape, x)
}
pub fn grad_nalgebra_val<F: Float + BtapeThreadLocal>(
f: impl FnOnce(&[BReverse<F>]) -> BReverse<F>,
x: &DVector<F>,
) -> (F, DVector<F>) {
let (mut tape, val) = crate::api::record(f, x.as_slice());
(val, tape_gradient_nalgebra(&mut tape, x))
}
pub fn hessian_nalgebra<F: Float + BtapeThreadLocal>(
f: impl FnOnce(&[BReverse<F>]) -> BReverse<F>,
x: &DVector<F>,
) -> (F, DVector<F>, DMatrix<F>) {
let (tape, _) = crate::api::record(f, x.as_slice());
tape_hessian_nalgebra(&tape, x)
}
pub fn jacobian_nalgebra<F: Float + BtapeThreadLocal>(
f: impl FnOnce(&[BReverse<F>]) -> Vec<BReverse<F>>,
x: &DVector<F>,
) -> DMatrix<F> {
let xs = x.as_slice();
let (mut tape, _) = crate::api::record_multi(f, xs);
let jac = tape.jacobian(xs);
let m = jac.len();
let n = if m > 0 { jac[0].len() } else { xs.len() };
DMatrix::from_fn(m, n, |i, j| jac[i][j])
}
pub fn tape_gradient_nalgebra<F: Float>(tape: &mut BytecodeTape<F>, x: &DVector<F>) -> DVector<F> {
let g = tape.gradient(x.as_slice());
DVector::from_vec(g)
}
#[must_use]
pub fn tape_hessian_nalgebra<F: Float>(
tape: &BytecodeTape<F>,
x: &DVector<F>,
) -> (F, DVector<F>, DMatrix<F>) {
let xs = x.as_slice();
let (val, grad, hess) = tape.hessian(xs);
let n = xs.len();
(
val,
DVector::from_vec(grad),
DMatrix::from_fn(n, n, |i, j| hess[i][j]),
)
}