code: def _hessian(x, grid):
      """find hessian of ``f(x)``, sampled on the coordinate grid defined by ``x``
  
      Args:
        x (ndarray): an array of shape ``(npts, dim)`` or ``(npts,)``
        grid (ndarray): irregular coordinate grid generated from ``z = f(x)``
  
      Returns:
        array of shape indicated in Notes, the hessian of the points on grid
  
      Notes:
        - Output will be of the form ``(dim,dim) + grid.shape``, where the output
          ``hess[i,j]`` corresponds to the second derivative ``z_{i,j}`` with
          ``i,j`` in ``range(dim)``. For a 1D array ``x``, the output will be a
          1D array of the same length. The hessian is calculated using finite
          differences.
      """  # XXX: can unique be used in this function?
      import numpy as np
  
      x = np.asarray(x)
      hess = np.empty(
          (grid.ndim, grid.ndim) + grid.shape, dtype=grid.dtype
      )
      if (
          grid.ndim == 1
      ):  # XXX: is (1,1,N) really desirable when x is (N,1)?
          hess[0, 0] = _gradient(x, _gradient(x, grid))
          return hess.ravel() if x.ndim == 1 else hess
      for k, grad_k in enumerate(_gradient(x, grid)):
          # apply gradient to every component of the first derivative
          for l, grad_kl in enumerate(_gradient(x, grad_k)):
              hess[k, l] = grad_kl
      return hess
  
imports: from mystic.math import _rbf
label_desc: Use image node_id 25 for calling a global function with example usage: _hessian(x, grid) and returns: hess, unknown