@@ -536,14 +536,15 @@ function init_optimization!(
536536 # Test new experimental feature:
537537
538538 model = mpc. estim. model
539- grad, jac = mpc . gradient, mpc. jacobian
539+ jac = mpc. jacobian
540540 nu, ny, nx̂, nϵ, nk = model. nu, model. ny, mpc. estim. nx̂, mpc. nϵ, model. nk
541541 Hp, Hc = mpc. Hp, mpc. Hc
542542 ng, nc, neq = length (mpc. con. i_g), mpc. con. nc, mpc. con. neq
543543 nZ̃, nU, nŶ, nX̂, nK = length (mpc. Z̃), Hp* nu, Hp* ny, Hp* nx̂, Hp* nk
544544 nΔŨ, nUe, nŶe = nu* Hc + nϵ, nU + nu, nŶ + ny
545545 strict = Val (true )
546546 myNaN = convert (JNT, NaN )
547+ myInf = convert (JNT, Inf )
547548
548549 ΔŨ:: Vector{JNT} = zeros (JNT, nΔŨ)
549550 x̂0end:: Vector{JNT} = zeros (JNT, nx̂)
@@ -554,9 +555,47 @@ function init_optimization!(
554555 gc:: Vector{JNT} , g:: Vector{JNT} = zeros (JNT, nc), zeros (JNT, ng)
555556 geq:: Vector{JNT} = zeros (JNT, neq)
556557
557- geq_min = zeros (JNT, mpc. con. neq)
558- geq_max = zeros (JNT, mpc. con. neq)
558+ function gfunc! (g, Z̃, ΔŨ, x̂0end, Ue, Ŷe, U0, Ŷ0, Û0, K0, X̂0, gc, geq)
559+ update_predictions! (ΔŨ, x̂0end, Ue, Ŷe, U0, Ŷ0, Û0, K0, X̂0, gc, g, geq, mpc, Z̃)
560+ return nothing
561+ end
562+ Z̃_∇g = fill (myNaN, nZ̃) # NaN to force update_predictions! at first call
563+ ∇g_context = (
564+ Cache (ΔŨ), Cache (x̂0end), Cache (Ue), Cache (Ŷe), Cache (U0), Cache (Ŷ0),
565+ Cache (Û0), Cache (K0), Cache (X̂0),
566+ Cache (gc), Cache (geq),
567+ )
568+ # temporarily enable all the inequality constraints for sparsity detection:
569+ mpc. con. i_g[1 : end - nc] .= true
570+ ∇g_prep = prepare_jacobian (gfunc!, g, jac, Z̃_∇g, ∇g_context... ; strict)
571+ mpc. con. i_g[1 : end - nc] .= false
572+ ∇g = init_diffmat (JNT, jac, ∇g_prep, nZ̃, ng)
573+
574+
575+ function gfunc_set! (g, Z̃)
576+ return gfunc! (g, Z̃, ΔŨ, x̂0end, Ue, Ŷe, U0, Ŷ0, Û0, K0, X̂0, gc, geq)
577+ end
578+ function ∇gfunc_set! (∇g_vec, Z̃)
579+ value_and_jacobian! (gfunc!, g, ∇g, ∇g_prep, jac, Z̃, ∇g_context... )
580+ ∇g_vec .= nonzeros (∇g)
581+ return nothing
582+ end
583+
584+ g_min = fill (- myInf, ng)
585+ g_max = zeros (JNT, ng)
586+
587+ I_∇g, J_∇g = SparseArrays. findnz (∇g)
588+ ∇g_structure = collect (zip (I_∇g, J_∇g))
559589
590+ g_set = Ipopt. _VectorNonlinearOracle (;
591+ dimension = nZ̃,
592+ l = g_min,
593+ u = g_max,
594+ eval_f = gfunc_set!,
595+ jacobian_structure = ∇g_structure,
596+ eval_jacobian = ∇gfunc_set!
597+ )
598+ @constraint (optim, Z̃var in g_set)
560599
561600 function geqfunc! (geq, Z̃, ΔŨ, x̂0end, Ue, Ŷe, U0, Ŷ0, Û0, K0, X̂0, gc, g)
562601 update_predictions! (ΔŨ, x̂0end, Ue, Ŷe, U0, Ŷ0, Û0, K0, X̂0, gc, g, geq, mpc, Z̃)
@@ -580,6 +619,10 @@ function init_optimization!(
580619 return nothing
581620 end
582621
622+
623+ geq_min = zeros (JNT, mpc. con. neq)
624+ geq_max = zeros (JNT, mpc. con. neq)
625+
583626 I_∇geq, J_∇geq = SparseArrays. findnz (∇geq)
584627 ∇geq_structure = collect (zip (I_∇geq, J_∇geq))
585628
@@ -594,15 +637,16 @@ function init_optimization!(
594637 =#
595638
596639
597- set = Ipopt. _VectorNonlinearOracle (;
640+ geq_set = Ipopt. _VectorNonlinearOracle (;
598641 dimension = nZ̃,
599642 l = geq_min,
600643 u = geq_max,
601644 eval_f = geqfunc_set!,
602645 jacobian_structure = ∇geq_structure,
603646 eval_jacobian = ∇geqfunc_set!
604647 )
605- @constraint (optim, Z̃var in set)
648+ @constraint (optim, Z̃var in geq_set)
649+
606650 end
607651 return nothing
608652end
0 commit comments