如何使用Optim.jl的optimize函数求解带求和项的参数优化问题

实现方案 #

1. 额外参数传递：用闭包封装固定参数 #

Optim.jl的待优化函数不需要强制改参数格式，你可以通过闭包把μ、E、cutoff、目标数组$\langle n_c \rangle_i$、粒子总数$N$这些不参与优化的固定参数封装进去，只把待优化的$\beta$作为函数输入即可，你原本的pcalc函数不需要做任何修改：

using Optim, LinearAlgebra

# 你原有的pcalc函数直接保留即可
function pcalc(β, μ, E, cutoff)
    # 你原本的p_i计算逻辑
end

封装目标函数示例：

# 先定义所有固定参数，替换为你实际的取值
const μ = 1.0 
const E = rand(100) # 你的能级数组
const cutoff = 1e-3 
const avg_nc = rand(100) # 待匹配的目标数组$\langle n_c \rangle_i$
const N = 1000.0 # 粒子总数

# 闭包封装，仅保留待优化的β作为输入（Optim默认输入为数组格式，因此取β[1]为实际参数值）
function objective(β)
    p = pcalc(β[1], μ, E, cutoff)
    return sum((p .- avg_nc).^2)
end

2. 约束条件处理 #

针对$\sum_i p_i = N$的等式约束，使用Optim的IPNewton求解器即可直接处理，约束函数同样用闭包封装固定参数：

# 定义等式约束函数，要求sum(p_i) - N = 0
function constraint!(c, β)
    p = pcalc(β[1], μ, E, cutoff)
    c[1] = sum(p) - N
end

# 等式约束上下界均为0
lc = [0.0]
uc = [0.0]

# 定义优化问题，[1.0]为初始β猜测值，可替换为你的实际初始值
df = TwiceDifferentiable(objective, [1.0])
dfc = TwiceDifferentiableConstraints(constraint!, lc, uc, [1.0])

# 调用优化器求解
res = optimize(df, dfc, [1.0], IPNewton())

3. 加速优化的可选配置 #

如果你的pcalc是无副作用的纯数值计算，可以开启自动微分代替手动求导，大幅降低计算耗时：

using ForwardDiff
# 定义时指定用ForwardDiff自动计算梯度、黑塞矩阵
df = TwiceDifferentiable(objective, [1.0], autodiff=:forward)
dfc = TwiceDifferentiableConstraints(constraint!, lc, uc, [1.0], autodiff=:forward)

优化完成后通过Optim.minimizer(res)即可获取最优$\beta$取值，Optim.minimum(res)可获取最小化后的目标函数值。

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