# GalactimOptim + AutoForwardDiff + IPNewton problem

**URL:** <https://discourse.julialang.org/t/galactimoptim-autoforwarddiff-ipnewton-problem/72508>\
**Category:** Optimization (Mathematical)\
**Created:** [December 3, 2021, 9:10am UTC](https://discourse.julialang.org/t/galactimoptim-autoforwarddiff-ipnewton-problem/72508 "2021-12-03T09:10:02Z")\
**Posts on this page:** 7\
**Page:** 1

<div class="post-metadata">

**Author:** ![mottelet](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mottelet/32/31089_2.png) [@mottelet](https://discourse.julialang.org/u/mottelet)\
**Post date:** [December 3, 2021, 9:10am UTC](https://discourse.julialang.org/t/galactimoptim-autoforwarddiff-ipnewton-problem/72508/1 "2021-12-03T09:10:02Z")

</div>

Hi, I am trying to solve a parameter fitting problem in a stiff biological ODE. The problem solves with correct precision in Matlab (with fmincon/IP) using the complex-step technique to compute the gradient, but since I am migrating a lot of stuff to Julia this example is very important to me. I made a first port and used Flux but constraints are not handled and the choice of solvers is small. I adapted the code for GalacticOptim, but the Optim.IPNewton() solver does not seem to like by combination:

```julia
ERROR: LoadError: Use OptimizationFunction to pass the derivatives or automatically generate them with one of the autodiff backends

```

Optim.BFGS() likes it, but fails to converge correctly.

Really guys, this is very HARD to migrate stuff to Julia, even for motivated guys like me (Julia code below).

S.

```julia
function syst_AM2(state,theta,t) 
    Y_x1,Y_x2,mu1_max,mu2_max,Ks1,Ks2,KI,mu0,S00,S10,S20,coeff = theta;
    X1,X2,S0,S1,S2,CH4 = state
    r0=mu0*S0;
    r1=mu1_max*S1*X1/(S1+Ks1);
    r2=mu2_max*S2*X2/(S2+Ks2+(S2^2/KI));
    [Y_x1*r1, Y_x2*r2, -r0, -r1+r0, (1-Y_x1)*r1-r2, (1-Y_x2)*r2]
end

function loss(theta,p)
    state0 = [X10;X20;theta[9:11];0];
    coeff = theta[12];
    prob = ODEProblem(syst_AM2,state0,(t_CH4[1],t_CH4[end]),theta)
    state = solve(prob, Rodas5(),saveat=t_CH4)
    S2 = state[5,:]
    CH4 = state[6,:]
    deb_CH4 = CH4[2:end]-CH4[1:end-1]
    e1 = CH4[1:end-1]-Exp_CH4*1
    e2 = (deb_CH4-Exp_deb_CH4)*4
    e3 = (S2[1:end-1]-coeff*Exp_AGV)*1
    (sum(e1.^2)+sum(e2.^2)+sum(e3.^2))/length(t_CH4);
end

function mycallback(p, l)
    display(loss(p))
    return false
end

using Plots, DifferentialEquations, GalacticOptim, Optim

H=0.39;   
R=0.195; 
Surf=pi*R^2; 
Exp_CH4 = [0,0.146695751,0.452562719,0.85245185,1.390508363,1.88574014,2.333176518,2.79329662,3.417909766,4.299654989,6.515952412,10.29002283,13.99900971,16.46767139,18.77271153,21.23159231,23.82660177,26.81067853,30.00135498,32.59923757,34.18509462,35.5158513,36.74890559,37.86754592,38.79872368,39.69511074,40.56481051,41.3865614,42.13182338,42.84323204,43.49045536];

t_CH4=0.0:size(Exp_CH4)[1]; # temps (jours)
Exp_AGV = [0.6307493,1.8081079,3.0002554,3.4371393,4.1710245,4.306228,4.4414315,4.576635,5.0225553,5.2529314,4.6852442,3.6354049,3.1327115,2.6300181,2.1273247,1.9948686,1.8574543,1.5407059,1.2099629,1.1136956,1.0174283,0.9211824,0.9121245,0.9195297,0.8298519,0.9236508,0.9041827,0.8846932,0.8652251,0.8372571,0.808326];

Exp_deb_CH4 = 0*Exp_CH4;
Exp_deb_CH4[2:end]= Exp_CH4[2:end]-Exp_CH4[1:end-1];

X10 = 0.283861214/(H*Surf);   
X20 = 0.094620405/(H*Surf);
S00 = 5.180333444/(H*Surf);
S10 = 0.01/(H*Surf); 
S20 = 0.029385961/(H*Surf); 

theta_c0 = [0.25, 0.43, 0.5, 3.1, 21.0, 23.0, 4.0, 0.04, S00, S10, S20, 1.0]

cstr = [
0 1; #Y_X1
0 1; #Y_X2
0 Inf; #mu1max
0 Inf; #mu2max
1e-3 Inf; #Ks1
1e-3 Inf; #Ks2
1e-3 Inf; #KI
0 Inf; # mu0
0 Inf; # S00
0 Inf; # S10
0 Inf; # S20
0 Inf]; # AGV coefficient

fun = OptimizationFunction(loss, GalacticOptim.AutoForwardDiff())
prob = OptimizationProblem(fun, theta_c0, [], lb=cstr[:,1], ub=cstr[:,2])
sol = solve(prob, Optim.IPNewton(), cb=mycallback)

```

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [December 4, 2021, 3:32pm UTC](https://discourse.julialang.org/t/galactimoptim-autoforwarddiff-ipnewton-problem/72508/2 "2021-12-04T15:32:12Z")

</div>

That was a bad error message. Fixed: [Better error message for constrained optimizers by ChrisRackauckas · Pull Request #198 · SciML/Optimization.jl · GitHub](https://github.com/SciML/GalacticOptim.jl/pull/198) .

> [@mottelet](#):
>
> Optim.BFGS() likes it, but fails to converge correctly.

NLopt BFGS acts similarly, so I don’t think it’s a problem of Optim as much of BFGS here. I found that `NLopt.LD_TNEWTON_PRECOND_RESTART` works best here.

```julia
using NLopt
sol = solve(prob, NLopt.LD_LBFGS(), cb=mycallback, maxiters = 10000).minimum # 4.156774187867866
sol = solve(prob, NLopt.LD_TNEWTON_PRECOND_RESTART(), cb=mycallback, maxiters = 10000).minimum # 3.2253551599824806

```

Note that one thing you can do to help optimizations is help it exit correctly when the solver fails, i.e.

```julia
    if state.retcode != :Success
        return Inf
    end

```

---

<div class="post-metadata">

**Author:** ![mottelet](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mottelet/32/31089_2.png) [@mottelet](https://discourse.julialang.org/u/mottelet)\
**Post date:** [December 6, 2021, 8:21am UTC](https://discourse.julialang.org/t/galactimoptim-autoforwarddiff-ipnewton-problem/72508/3 "2021-12-06T08:21:13Z")

</div>

Hi Chris. Thanks for the fix and advices. BTW, the display callback silently fails (there is no display), is the keyword different for NLopt solvers ? Another question: how do you get the final gradient ? Maybe using GalacticOptim is too high level for me (what would you advise as an alternative to it, still having easy access to AD) ?

S.

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [December 6, 2021, 3:12pm UTC](https://discourse.julialang.org/t/galactimoptim-autoforwarddiff-ipnewton-problem/72508/4 "2021-12-06T15:12:35Z")

</div>

Open an issue so it shows up in my emails. It’s like a one line fix to fix that but I just haven’t gotten to it and right now it’ll get lost in the v1.7 update stuff, but an email could probably put it in queue for tonight.

---

<div class="post-metadata">

**Author:** ![mottelet](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mottelet/32/31089_2.png) [@mottelet](https://discourse.julialang.org/u/mottelet)\
**Post date:** [December 6, 2021, 3:26pm UTC](https://discourse.julialang.org/t/galactimoptim-autoforwarddiff-ipnewton-problem/72508/5 "2021-12-06T15:26:13Z")

</div>

OK, where do I file the issue ? NLopt.jl ?

S.

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [December 6, 2021, 3:31pm UTC](https://discourse.julialang.org/t/galactimoptim-autoforwarddiff-ipnewton-problem/72508/6 "2021-12-06T15:31:45Z")

</div>

GalacticOptim.jl

---

<div class="post-metadata">

**Author:** ![mottelet](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mottelet/32/31089_2.png) [@mottelet](https://discourse.julialang.org/u/mottelet)\
**Post date:** [December 6, 2021, 3:44pm UTC](https://discourse.julialang.org/t/galactimoptim-autoforwarddiff-ipnewton-problem/72508/7 "2021-12-06T15:44:41Z")

</div>

[https://github.com/SciML/GalacticOptim.jl/issues/199](https://github.com/SciML/GalacticOptim.jl/issues/199)