跳跃:期望float64,typeError的问题:在TypeSert中,预期的float64,forwarddiff.dual autodiff = true和exp()的问题
因此,我试图以最少的示例来根据我编写的更复杂的代码提出问题:
- 我遇到的一个巨大的常见错误是float64,而是forwarddiff.dual-偶然 - 有人可以给我一个提示。我总是确保避免使用此错误。我觉得每当我做一个新的优化问题时,我都必须重新发明方向盘以使其
- 显然无法自动化Julia Exp()函数?有人知道如何使它起作用吗?
- 我做了一个有限的总和以通过泰勒系列近似它。在我的一个功能中,如果我有20个术语,自动驾驶仪的工作还不够准确 - 所以我去了40个任期Autodiff现在不起作用 - 有人为此解决了吗?
任何建议将不胜感激!
using Cubature
using Juniper
using Ipopt
using JuMP
using LinearAlgebra
using Base.Threads
using Cbc
using DifferentialEquations
using Trapz
function mat_exp(x::AbstractVector{T},dim,num_terms,A) where T
sum = zeros(Complex{T},(dim,dim))
A[1,1] = A[1,1]*x[1]
A[2,2] = A[2,2]*x[2]
return exp(A)-1
end
function exp_approx_no_big(x::AbstractVector{T},dim,num_terms,A) where T
sum = zeros(Complex{T},(dim,dim))
A[1,1] = A[1,1]*x[1]
A[2,2] = A[2,2]*x[2]
for k=0:num_terms-1
sum = sum + (1.0/factorial(k))*A^k
end
return norm(sum)-1
end
function exp_approx_big(x::AbstractVector{T},dim,num_terms,A) where T
sum = zeros(Complex{T},(dim,dim))
A[1,1] = A[1,1]*x[1]
A[2,2] = A[2,2]*x[2]
for k=0:num_terms-1
sum = sum + (1.0/factorial(big(k)))*A^k
end
return norm(sum)-1
end
optimizer = Juniper.Optimizer
nl_solver= optimizer_with_attributes(Ipopt.Optimizer, "print_level" => 0)
mip_solver = optimizer_with_attributes(Cbc.Optimizer, "logLevel" => 0, "threads"=>nthreads())
m = Model(optimizer_with_attributes(optimizer, "nl_solver"=>nl_solver, "mip_solver"=>mip_solver))
@variable(m, 0.0<=x[1:2]<=1.0)
dim=5
A=zeros(Complex,(dim,dim))
for k=1:dim
A[k,k]=1.0
end
println(A)
f(x...) = exp_approx_no_big(collect(x),dim,20,A)
g(x...) = exp_approx_big(collect(x),dim,40,A)
h(x...) = mat_exp(collect(x),dim,20,A)
register(m, :f, 2, f; autodiff = true)
@NLobjective(m, Min, f(x...))
optimize!(m)
println(JuMP.value.(x))
println(JuMP.objective_value(m))
println(JuMP.termination_status(m))
So I tried to make a minimum example to ask questions based on a more complicated piece of code I have written:
- A HUGE common error I'm getting is expecting float64 and instead got ForwardDiff.Dual - can someone give me a tip how in general I always make sure I avoid this bug. I feel like every time I do a new optimization problem I have to reinvent the wheel to try to make this go away
- Apparently you cannot autodiff the julia exp() function? Does anyone know how to make it work?
- A workaround is I did a finite sum to approximate it via the taylor series. In my one function if I had 20 terms the autodiff worked, but it wasn't accurate enough - so I went to 40 terms but then julia told me to do factorial(big(k)) and then when I try to do that with autodiff it doesn't work now - anyone have a fix for this?
Any advice would be greatly appreciated!
using Cubature
using Juniper
using Ipopt
using JuMP
using LinearAlgebra
using Base.Threads
using Cbc
using DifferentialEquations
using Trapz
function mat_exp(x::AbstractVector{T},dim,num_terms,A) where T
sum = zeros(Complex{T},(dim,dim))
A[1,1] = A[1,1]*x[1]
A[2,2] = A[2,2]*x[2]
return exp(A)-1
end
function exp_approx_no_big(x::AbstractVector{T},dim,num_terms,A) where T
sum = zeros(Complex{T},(dim,dim))
A[1,1] = A[1,1]*x[1]
A[2,2] = A[2,2]*x[2]
for k=0:num_terms-1
sum = sum + (1.0/factorial(k))*A^k
end
return norm(sum)-1
end
function exp_approx_big(x::AbstractVector{T},dim,num_terms,A) where T
sum = zeros(Complex{T},(dim,dim))
A[1,1] = A[1,1]*x[1]
A[2,2] = A[2,2]*x[2]
for k=0:num_terms-1
sum = sum + (1.0/factorial(big(k)))*A^k
end
return norm(sum)-1
end
optimizer = Juniper.Optimizer
nl_solver= optimizer_with_attributes(Ipopt.Optimizer, "print_level" => 0)
mip_solver = optimizer_with_attributes(Cbc.Optimizer, "logLevel" => 0, "threads"=>nthreads())
m = Model(optimizer_with_attributes(optimizer, "nl_solver"=>nl_solver, "mip_solver"=>mip_solver))
@variable(m, 0.0<=x[1:2]<=1.0)
dim=5
A=zeros(Complex,(dim,dim))
for k=1:dim
A[k,k]=1.0
end
println(A)
f(x...) = exp_approx_no_big(collect(x),dim,20,A)
g(x...) = exp_approx_big(collect(x),dim,40,A)
h(x...) = mat_exp(collect(x),dim,20,A)
register(m, :f, 2, f; autodiff = true)
@NLobjective(m, Min, f(x...))
optimize!(m)
println(JuMP.value.(x))
println(JuMP.objective_value(m))
println(JuMP.termination_status(m))
如果你对这篇内容有疑问,欢迎到本站社区发帖提问 参与讨论,获取更多帮助,或者扫码二维码加入 Web 技术交流群。

绑定邮箱获取回复消息
由于您还没有绑定你的真实邮箱,如果其他用户或者作者回复了您的评论,将不能在第一时间通知您!
发布评论
评论(1)
您的
mat_exp
函数有很多问题:a
就位,因此重复的调用不会执行您认为exp(x) -1
,是矩阵。跳转仅支持标量调用norm(EXP(x))-1
exp
进行区分,我也不知道为什么使用杜松,或安装了许多其他软件包。
如果您想对此进行讨论,请加入社区论坛: https:// usise.julialang.org/c/domain/opt/13 。 (对于来回而言,这比StackoverFlow要好得多。)有人可能有建议,但我不知道朱莉娅(Julia)中的广告工具可以通过矩阵指数分化。
There are quite a few problems with your
mat_exp
function:A
in-place, so repeated calls will not do what you thinkexp(x) - 1
, which is a matrix. JuMP only supports scalar callsnorm(exp(x)) - 1
exp
I also don't know why you're using Juniper, or that you have a bunch of other packages installed.
If you want to have a discussion on this, come join the community forum: https://discourse.julialang.org/c/domain/opt/13. (It's much better for back-and-forth than stackoverflow.) Someone might have suggestions, but I don't know of an AD tool in Julia that can differentiate through a matrix exponential.