Trouble debugging gradients of nonlinear model
还没有人认领这个 Issue。
评估
调研方向
从使用 DiffOpt.nonlinear_diff_model、reverse_differentiate! 和有限差分函数 solve_a_ϕ0_manual 的精简版 Julia 复现程序开始。比较 ϕ0 = 0.0 和 ϕ0 = 1.0 时报告的导数，然后确定差异是由数值行为还是由 DiffOpt 的非线性感度计算造成的。原因得到解释，并且任何所需的修正都由可复现的测试覆盖后，即视为完成。
由索引模型根据 Issue 内容生成。
描述
I have a nonlinear (nonconvex) program that I am solving with Ipopt and using DiffOpt to differentiate through. I have been having trouble getting correct gradients and am not sure if this is a numerical issue with my model formulation, or if there is something more fundamental going on. Here is a pared-down example that demonstrates my issue.
# define parameter value
ϕ0=[0.0];
# define constants used in model
Phigh = 0.4412123910369964;
Plow = 0.357171935601378;
c(Nsep) = (1 / (1 + 2*Nsep))^2
cols = 2;
reqwakes = 1;
Nwakes = 1;
aprev = 0.1 * ones(2, 1);
rhist = Phigh * ones(1);
normalized_rhist = [1];
predhorz = 2;
KP = 0.5 * π/4;
model = DiffOpt.nonlinear_diff_model(Ipopt.Optimizer)
# define decision variables
JuMP.@variable(model, 0.05 <= a[1:cols, -reqwakes:Nwakes] <= 1/3)
# fix certain values
for j=1:cols, n=1:reqwakes
JuMP.fix(a[j, n - cols], aprev[j, n]; force=true)
end
JuMP.fix.(a[end, 0:end], 0.333333; force=true)
# define parameter
JuMP.@variable(model, ϕ in JuMP.Parameter(ϕ0[1]))
# define parameterized signal
scale = (Phigh - Plow) / 2;
shift = (Phigh + Plow) / 2;
JuMP.@variables(model, begin
Pref_unscaled[t=1:predhorz]
Pref_var[t=1:predhorz]
end)
JuMP.fix.(Pref_unscaled[1], reverse(normalized_rhist))
JuMP.@constraints(model, begin
[t=2:predhorz], Pref_unscaled[t] .== ϕ * Pref_unscaled[t-1]
[t=1:predhorz], Pref_var[t] == (scale * Pref_unscaled[t] + shift)
end)
# define NL model expressions
JuMP.@expressions(model, begin
CP[j=1:cols, n=0:Nwakes], 4 * a[j,n] * (1 - a[j,n])^2
δV[j=1:cols, n=0:Nwakes], 2 * sqrt( sum( (a[prev,n-(j-prev)] * c(j-prev) )^2 for prev=1:j-1; init=0) )
u[j=1:cols, n=0:Nwakes], (1 - δV[j,n])
Pcap[j=1:cols, n=0:Nwakes], CP[j,n] * u[j,n]^3
Pgen[n=0:Nwakes], sum(Pcap[j,n] for j=1:cols)
end)
# define objective
JuMP.@objective(model, Min, sum( (KP * Pgen[n] - Pref_var[1+n])^2 for n=0:Nwakes))
With this model, I vary φ and find the optimal a, producing the following curve.
I then take the derivative of a w.r.t. φ using both DiffOpt
DiffOpt.empty_input_sensitivities!(model)
DiffOpt.set_reverse_variable.(model, model[:a][1, 0], 1)
DiffOpt.reverse_differentiate!(model)
dadϕ = DiffOpt.get_reverse_parameter.(model, model[:ϕ]);
and finite differencing.
function solve_a_ϕ0_manual(ϕ0)
JuMP.set_parameter_value.(model[:ϕ], ϕ0)
JuMP.optimize!(model)
JuMP.assert_is_solved_and_feasible(model)
asolved_new = JuMP.value.(model[:a]);
return asolved_new[1, 0]
end
dadϕ_fdm = FiniteDifferences.grad(FiniteDifferences.forward_fdm(2, 1; factor=1e6), solve_a_ϕ0_manual, ϕ0);
This returns problematic results. In the flat portion of the plot, e.g., φ0 = 0.0, both gradients are "close" to 0, but still differ significantly.
dadϕ = -9.871980439887741e-10
dadϕ_fdm = ([-5.972732875606716e-8],)
dadϕ_fdm[1] / dadϕ = [60.501871047818184]
In the changing portion of the plot, e.g., φ0 = 1.0, the gradients are closer but don't really match well.
dadϕ = -0.03988513330240084
dadϕ_fdm = ([-0.04102676532022559],)
dadϕ_fdm[1] / dadϕ = [1.02862299617176]
Plotting the gradients as tangents on the graphs shows that the finite differencing is giving me more accurate results here.
Any insight would be appreciated!
- 主要语言
- Julia
- 星标
- 146
- 派生
- 21
- 平均合并
- 1 小时 25 分钟
- 30 天内合并 PR
- 2
环境准备
这个项目没有提供开发容器、Dockerfile 或贡献指南，环境需要你自己搭建：先看它的 README，通用步骤见我们的新手贡献指南。
从这里开始
- 先读完整个 Issue，再读项目的贡献指南。
- 在 Issue 下留言说明你要接手 —— 这能避免两个人做同样的事。
- Fork 仓库，在一个分支上完成修改。
- 提交 Pull Request，并在描述里引用这个 Issue 编号。
jump-dev/DiffOpt.jl 的其他 Issue
-
难度 3/5 1-2 天 新手友好度 55/100
jump-dev/DiffOpt.jl#362 · 1 条评论 ·
-
难度 5/5 一周以上 新手友好度 25/100
jump-dev/DiffOpt.jl#358 · 2 条评论 ·
-
难度 3/5 1-2 天 新手友好度 55/100
jump-dev/DiffOpt.jl#348 ·
-
good first issue
难度 4/5 3-5 天 新手友好度 48/100
jump-dev/DiffOpt.jl#336 · 4 条评论 ·
-
难度 5/5 一周以上 新手友好度 35/100
jump-dev/DiffOpt.jl#331 · 1 条评论 ·
查看 jump-dev/DiffOpt.jl 的全部 Issue
相似的 Issue
-
难度 1/5 1 小时以内 新手友好度 83/100
JuliaMath/Polynomials.jl#647 ·
-
难度 2/5 1-3 小时 新手友好度 74/100
OpenSourceAWE/SymbolicAWEModels.jl#431 ·
维护者通常 1 天内回复
-
难度 2/5 1-3 小时 新手友好度 66/100
-
难度 1/5 1 小时以内 新手友好度 89/100
维护者通常 4 天内回复
-
难度 2/5 1-3 小时 新手友好度 76/100
维护者通常 1 天内回复