求解多维函数Hessian矩阵：如何将10维凸函数转为数组输入

你当前实现的hessian函数是基于有限差分逻辑，对已在网格上完成采样的函数值数组计算海森矩阵，不支持直接输入函数表达式。要实现构造10维凸函数的需求，按以下步骤操作：

• 第一步：定义10维凸函数表达式，最简单的严格凸函数可选正定二次函数：
  公式为 f(x) = Σ(i=1到10) i * x_i²，该函数的海森矩阵为对角矩阵，对角元素为2、4、…、20，所有特征值均为正，天然满足凸函数要求。
• 第二步：生成10维采样网格，为10个自变量分别设置采样区间与步长，比如每个自变量取[-2,2]区间内的5个均匀采样点，最终生成形状为(5,5,5,5,5,5,5,5,5,5)的10维坐标网格。
• 第三步：将定义好的凸函数在10维网格上逐点计算，得到所有网格点的函数值数组，该数组就是hessian函数可接收的输入参数。

现有测试代码的问题说明 #

你当前测试代码中的x = np.random.randn(100,100)是二维随机数组，对应的是二元函数的采样值，因此计算得到的海森矩阵也对应二元函数，和你需要的10维场景不匹配。

简化使用示例（3维场景，可类推到10维） #

import numpy as np

def hessian(x):
    """
    用有限差分计算海森矩阵
    参数:
       - x : ndarray 网格上的函数值数组
    返回:
       形状为 (x.ndim, x.ndim) + x.shape 的数组
       数组元素 array[i, j, ...] 对应二阶导数 x_ij
    """
    x_grad = np.gradient(x) 
    hessian = np.empty((x.ndim, x.ndim) + x.shape, dtype=x.dtype) 
    for k, grad_k in enumerate(x_grad):
        tmp_grad = np.gradient(grad_k) 
        for l, grad_kl in enumerate(tmp_grad):
            hessian[k, l, :, :] = grad_kl
    return hessian

# 生成3维网格，每个维度取10个采样点
x1 = np.linspace(-2, 2, 10)
x2 = np.linspace(-2, 2, 10)
x3 = np.linspace(-2, 2, 10)
X1, X2, X3 = np.meshgrid(x1, x2, x3, indexing='ij')
# 定义凸函数 f = x1² + 2x2² +3x3²
f_val = X1**2 + 2 * X2**2 + 3 * X3**2
# 计算海森矩阵
hess = hessian(f_val)

更高效的凸函数构造方案 #

如果你的核心目标只是构造10维凸函数，不需要用有限差分计算海森，可直接构造10阶对称正定矩阵作为海森：对二次型函数 f(x) = 0.5 * x^T A x，只要A是对称正定矩阵，f(x)就是严格凸函数，可直接通过构造正定矩阵A的方式快速得到符合要求的凸函数，不需要做采样差分计算。

内容的提问来源于stack exchange，提问作者grushenka