Error when attempting to calculate the Hessian of a vector-valued function

[dobos]

I'm not on the bleeding edge version but I'm getting this error when trying to calculate the Hessian of a function that returns a 1D array:

TypeError: only size-1 arrays can be converted to Python scalars

The above exception was the direct cause of the following exception:

Traceback (most recent call last):
  File "/home/dobos/project/ga_pfsspec_all/python/test/pfs/ga/pfsspec/stellar/rvfit/test_modelgridrvfit.py", line 155, in test_fit_rv
    rv, rv_err, params, params_err = rvfit.fit_rv({'mr': spec}, rv_0=100.0, rv_bounds=(0, 200), params_0=params_0, params_fixed=params_fixed, params_bounds=params_bounds)
  File "/home/dobos/project/ga_pfsspec_all/python/pfs/ga/pfsspec/stellar/rvfit/modelgridrvfit.py", line 155, in fit_rv
    self.calculate_F(spectra, rv, rv_bounds, None, params, params_bounds, None, params_free, params_fixed)
  File "/home/dobos/project/ga_pfsspec_all/python/pfs/ga/pfsspec/stellar/rvfit/modelgridrvfit.py", line 58, in calculate_F
    ddpc = ddphichi(self.params_to_array(rv_0, params_free, **params_0))
  File "/datascope/slurm/miniconda3/envs/astro-dnn/lib/python3.7/site-packages/numdifftools/core.py", line 849, in __call__
    return super(Hessdiag, self).__call__(np.atleast_1d(x), *args, **kwds)
  File "/datascope/slurm/miniconda3/envs/astro-dnn/lib/python3.7/site-packages/numdifftools/core.py", line 376, in __call__
    results = self._derivative(xi, args, kwds)
  File "/datascope/slurm/miniconda3/envs/astro-dnn/lib/python3.7/site-packages/numdifftools/core.py", line 260, in _derivative_nonzero_order
    results = [diff(f, fxi, xi, h) for h in steps]
  File "/datascope/slurm/miniconda3/envs/astro-dnn/lib/python3.7/site-packages/numdifftools/core.py", line 260, in <listcomp>
    results = [diff(f, fxi, xi, h) for h in steps]
  File "/datascope/slurm/miniconda3/envs/astro-dnn/lib/python3.7/site-packages/numdifftools/core.py", line 894, in _central_even
    hess[i, i] = (f(x + 2 * ee[i, :]) - 2 * fx + f(x - 2 * ee[i, :])) / (4. * hess[i, i])
ValueError: setting an array element with a sequence.

It works with the Gradient and the Jacobian. It would also be great if these accepted Nd arrays instead of just 1D arrays and the resulting new dimension would be the last (right now it's the second regardless of the number of input dimensions).

numdifftools 0.9.39 py_0 conda-forge