HESSIAN

=HESSIAN(expression, variables, [variable_names])

=HESSIAN("x**2 + 3*x*y + y**2", {{1.0, 2.0}}, {{"x", "y"}})

=HESSIAN("x**4 + y**4 + 2*x**2*y**2", {{1.0, 1.0}}, {{"x", "y"}})

=HESSIAN("x**2", {{2.0, 0.0}}, {{"x", "y"}})

=HESSIAN("y**2", {{0.0, 3.0}}, {{"x", "y"}})

import casadi as ca
 
# Note: scipy.differentiate module is not available in either Excel PY or Pyodide
 
def hessian(expression, variables, variable_names=None):
    """Compute the Hessian matrix (second derivatives) of a scalar function using CasADi symbolic differentiation.
 
    Computes the matrix of second-order partial derivatives for accurate analysis of function curvature,
    optimization, and sensitivity analysis. See [CasADi documentation](https://web.casadi.org/).
 
    Args:
        expression (str): Scalar function as a string expression. Example: "x**2 + 3*x*y + y**2".
        variables (list[list]): 2D list of numeric values representing the point at which to evaluate
            the Hessian. Example: [[1.0, 2.0]].
        variable_names (list[list], optional): 2D list of variable names as strings. If not provided,
            names are inferred from the expression. Example: [["x", "y"]].
 
    Returns:
        list[list[float]]: The Hessian matrix evaluated at the given point, as a 2D list of floats.
        str: Error message if calculation fails.
 
    This example function is provided as-is without any representation of accuracy.
    """
 
    def to2d(x):
        return [[x]] if not isinstance(x, list) else x
 
    if not isinstance(expression, str) or not expression.strip():
        return "expression must be a non-empty string."
 
    variables = to2d(variables)
 
    if not variables or len(variables[0]) == 0:
        return "variables must be a 2D list of numeric values."
 
    try:
        var_vals = [float(v) for v in variables[0]]
    except (TypeError, ValueError):
        return "variables must contain numeric values."
 
    if len(var_vals) == 0:
        return "variables must include at least one value."
 
    if variable_names is not None:
        variable_names = to2d(variable_names)
        if variable_names and len(variable_names[0]) > 0:
            names = variable_names[0]
        else:
            return "variable_names must be a 2D list of strings."
 
        if len(names) != len(var_vals):
            return "Number of variable names must equal number of variable values."
 
        if not all(isinstance(name, str) for name in names):
            return "variable_names must contain string values."
    else:
        import re
        candidates = re.findall(r"[A-Za-z_]\w*", expression)
        names = []
        common_functions = {
            'sin', 'cos', 'tan', 'exp', 'log', 'sqrt', 'abs', 'pi', 'e',
            'ca', 'fabs', 'floor', 'ceil', 'sinh', 'cosh', 'tanh',
            'asin', 'acos', 'atan'
        }
        for candidate in candidates:
            if candidate not in names and candidate not in common_functions:
                names.append(candidate)
 
        if len(names) != len(var_vals):
            return "Unable to infer variable names matching the number of values."
 
    if len(names) == 0:
        return "No variable names found."
 
    try:
        sym_vars = [ca.SX.sym(name) for name in names]
 
        evaluation_context = {name: sym_vars[idx] for idx, name in enumerate(names)}
        evaluation_context["ca"] = ca
        evaluation_context["sin"] = ca.sin
        evaluation_context["cos"] = ca.cos
        evaluation_context["tan"] = ca.tan
        evaluation_context["exp"] = ca.exp
        evaluation_context["log"] = ca.log
        evaluation_context["sqrt"] = ca.sqrt
        evaluation_context["atan"] = ca.atan
        evaluation_context["asin"] = ca.asin
        evaluation_context["acos"] = ca.acos
        evaluation_context["sinh"] = ca.sinh
        evaluation_context["cosh"] = ca.cosh
        evaluation_context["tanh"] = ca.tanh
        evaluation_context["fabs"] = ca.fabs
        evaluation_context["floor"] = ca.floor
        evaluation_context["ceil"] = ca.ceil
 
        expr_obj = eval(expression, evaluation_context)
 
        H = ca.hessian(expr_obj, ca.vertcat(*sym_vars))[0]
        hess_func = ca.Function('hess_func', sym_vars, [H])
        result = hess_func(*var_vals)
 
        if isinstance(result, ca.DM):
            return result.full().tolist()
        else:
            return "Error during CasADi calculation: Unexpected result type."
 
    except Exception as e:
        return f"Invalid expression: {e}"