# Difficulty in using JuMP with Complex Functions
**URL:** https://discourse.julialang.org/t/difficulty-in-using-jump-with-complex-functions/126045
**Category:** Optimization (Mathematical)
**Tags:** question, jump, symbolics
**Created:** [February 18, 2025, 9:04pm UTC](https://discourse.julialang.org/t/difficulty-in-using-jump-with-complex-functions/126045 "2025-02-18T21:04:21Z")
**Posts on this page:** 4
**Page:** 1
### Author:  [@MechRob](https://discourse.julialang.org/u/MechRob)
#### Post date: [February 18, 2025, 9:04pm UTC](https://discourse.julialang.org/t/difficulty-in-using-jump-with-complex-functions/126045/1 "2025-02-18T21:04:21Z")
I am trying to create an optimization problem with JuMP. In order to create the constraints I have a function, which I generated its elements with Symbolics. When I pass numerical values, I can calculate the terms in tenth of a mili second, but when I pass Symbolics or JuMP variables, even after waiting two hours I couldn’t make it to work. Is there any work around that I can solve this problem? I will put the pseudo code here:
Here how each function is generated with symbolics:
```julia
@eval @fastmath const T_func(xi, r, L0, q) = $T
println("T_func_generated")
@eval @fastmath const TE1_diff_1_func(xi, r, L0, q) = $TE1_diff_1
println("TE1_diff_1_func_generated")
@eval @fastmath const TE1_diff_2_func(xi, r, L0, q) = $TE1_diff_2
println("TE1_diff_2_func_generated")
@eval @fastmath const TE1_diff_3_func(xi, r, L0, q) = $TE1_diff_3
println("TE1_diff_3_func_generated")
@eval @fastmath const p1_int_func(xi, r, L0, q) = $p1_int
println("p1_int_func_generated")
@eval @fastmath const p1_p1_transpose_int_func(xi, r, L0, q) = $p1_p1_transpose_int
println("p1_p1_transpose_int_func_generated")
@eval @fastmath const p1_diff_1_int_func(xi, r, L0, q) = $p1_diff_1_int
println("p1_diff_1_int_func_generated")
@eval @fastmath const p1_diff_2_int_func(xi, r, L0, q) = $p1_diff_2_int
println("p1_diff_2_int_func_generated")
@eval @fastmath const p1_diff_3_int_func(xi, r, L0, q) = $p1_diff_3_int
println("p1_diff_3_int_func_generated")
@eval @fastmath const p1_p1_diff_1_t_int_func(xi, r, L0, q) = $p1_p1_diff_1_t_int
println("p1_p1_diff_1_t_int_func_generated")
@eval @fastmath const p1_p1_diff_2_t_int_func(xi, r, L0, q) = $p1_p1_diff_2_t_int
println("p1_p1_diff_2_t_int_func_generated")
@eval @fastmath const p1_p1_diff_3_t_int_func(xi, r, L0, q) = $p1_p1_diff_3_t_int
println("p1_p1_diff_3_t_int_func_generated")
@eval @fastmath const M_func(xi,r,L0,q) = $M
```
This is the rest of the code that I wrap these elements to get my terms:
```julia
function Kin(L0,r,q)
T = T_func(1,r,L0,q)
TE1_diff_1 = TE1_diff_1_func(1,r,L0,q)
TE1_diff_2 = TE1_diff_2_func(1,r,L0,q)
TE1_diff_3 = TE1_diff_3_func(1,r,L0,q)
p1_int = p1_int_func(1,r,L0,q)
p1_p1_transpose_int = p1_p1_transpose_int_func(1,r,L0,q)
p1_diff_1_int = p1_diff_1_int_func(1,r,L0,q)
p1_diff_2_int = p1_diff_2_int_func(1,r,L0,q)
p1_diff_3_int = p1_diff_3_int_func(1,r,L0,q)
p1_p1_diff_1_t_int = p1_p1_diff_1_t_int_func(1,r,L0,q)
p1_p1_diff_2_t_int = p1_p1_diff_2_t_int_func(1,r,L0,q)
p1_p1_diff_3_t_int = p1_p1_diff_3_t_int_func(1,r,L0,q)
M = M_func(1,r,L0,q)
return T,TE1_diff_1,TE1_diff_2,TE1_diff_3,p1_int,p1_p1_transpose_int,p1_diff_1_int,p1_diff_2_int,p1_diff_3_int,p1_p1_diff_1_t_int,p1_p1_diff_2_t_int,p1_p1_diff_3_t_int, M
end
include("JC_func.jl")
include("DynamicsTerm.jl")
include("BiasT.jl")
include("RP.jl")
function Terms(q)
L0 = 0.15
beta1 = pi/3
beta2 = -pi/3
b1 = 40e-3
b2 = 40e-3
m1, m2, m3 = 1.17, 0.54, 0.265
r1, r2, r3 = 46.32e-3, 32e-3, 24e-3
g = [0,0,-9.81]
q1 = q[1:3]
q2 = q[4:6]
q3 = q[7:9]
TE1, TE1_diff_1, TE1_diff_2, TE1_diff_3, P1_int, p1_p1_int, p1_diff_1_int, p1_diff_2_int, p1_diff_3_int, p1_p1_diff_1_int, p1_p1_diff_2_int, p1_p1_diff_3_int, JJ1 = Kin(L0, r1, q1)
TE2, TE2_diff_1, TE2_diff_2, TE2_diff_3, P2_int, p2_p2_int, p2_diff_1_int, p2_diff_2_int, p2_diff_3_int, p2_p2_diff_1_int, p2_p2_diff_2_int, p2_p2_diff_3_int, JJ2 = Kin(L0, r2, q2)
TE3, TE3_diff_1, TE3_diff_2, TE3_diff_3, P3_int, p3_p3_int, p3_diff_1_int, p3_diff_2_int, p3_diff_3_int, p3_p3_diff_1_int, p3_p3_diff_2_int, p3_p3_diff_3_int, JJ3 = Kin(L0, r3, q3)
R1B, P1B = BiasT(beta1,b1)
R2B, P2B = BiasT(beta2,b2)
R1, P1 = RP(TE1);
R2, P2 = RP(TE2);
R3, P3 = RP(TE3);
R1_diff_1, p1_diff_1 = RP(TE1_diff_1)
R1_diff_2, p1_diff_2 = RP(TE1_diff_2)
R1_diff_3, p1_diff_3 = RP(TE1_diff_3)
R2_diff_1, p2_diff_1 = RP(TE2_diff_1)
R2_diff_2, p2_diff_2 = RP(TE2_diff_2)
R2_diff_3, p2_diff_3 = RP(TE2_diff_3)
R3_diff_1, p3_diff_1 = RP(TE3_diff_1)
R3_diff_2, p3_diff_2 = RP(TE3_diff_2)
R3_diff_3, p3_diff_3 = RP(TE3_diff_3)
(Jv1, Jw1, Jv2, Jw2, Jv3, Jw3, Jv1b, Jw1b, Jv2b, Jw2b, Jv1_base, Jw1_base, Jv2_base, Jw2_base, Jv3_base, Jw3_base, P_out, R_out) = JC(R1, R1B, R2, R2B, R3, P1, P1B, P2, P2B, P3,
p1_diff_1, p1_diff_2, p1_diff_3, p2_diff_1, p2_diff_2, p2_diff_3, p3_diff_1, p3_diff_2, p3_diff_3,
R1_diff_1, R1_diff_2, R1_diff_3, R2_diff_1, R2_diff_2, R2_diff_3, R3_diff_1, R3_diff_2, R3_diff_3)
(M,G) = Dyn(JJ1, JJ2, JJ3,
Jv1, Jw1, Jv2, Jw2, Jv3, Jw3,
P1_int, P2_int, P3_int,
p1_p1_int, p2_p2_int, p3_p3_int,
p1_diff_1_int, p1_diff_2_int, p1_diff_3_int,
p2_diff_1_int, p2_diff_2_int, p2_diff_3_int,
p3_diff_1_int, p3_diff_2_int, p3_diff_3_int,
p1_p1_diff_1_int, p1_p1_diff_2_int, p1_p1_diff_3_int,
p2_p2_diff_1_int, p2_p2_diff_2_int, p2_p2_diff_3_int,
p3_p3_diff_1_int, p3_p3_diff_2_int, p3_p3_diff_3_int,
m1, m2, m3,
R1, R1B, R2, R2B, R3,
g)
return P_out, R_out, Jv3_base, Jw3_base, M, G
end
```
when I do
`@time Terms(numerica_variables)`
I get sub,mili second time
but when I do
```julia
@variables qq[1:9]
@time Terms(qq)
```
It just get stuck unfortunately.
---
### Author:  [@odow](https://discourse.julialang.org/u/odow)
#### Post date: [February 18, 2025, 9:28pm UTC](https://discourse.julialang.org/t/difficulty-in-using-jump-with-complex-functions/126045/2 "2025-02-18T21:28:23Z")
Can you provide a reproducible example? I don’t see anything related to JuMP in your example.
The issue is almost certainly that the expression returned by `Terms` is too complicated. JuMP likes many variables and many constraints that are sparse, not a single complicated expression.
As a different approach, you might be able to write this as a user-defined operator:
> **[Nonlinear Modeling · JuMP](https://jump.dev/JuMP.jl/stable/manual/nonlinear/#jump_user_defined_operators)**
>
> Documentation for JuMP.
I will also caution you that using `@eval` is almost certainly the wrong thing to do. How and why are you constructing functions like this?
---
### Author:  [@MechRob](https://discourse.julialang.org/u/MechRob)
#### Post date: [February 18, 2025, 10:00pm UTC](https://discourse.julialang.org/t/difficulty-in-using-jump-with-complex-functions/126045/3 "2025-02-18T22:00:11Z")
Thank you for your reply. So this is a reduced version of the code that is reproducible. I will explain why I am writing the code in this way.
I had a very complicated symbolic expression that needed to be integrated. I couldn’t get good accuracy with Julia packages, so I used Sympy to calculate the integration. Then, I saved the Julia code generated by sympy as `M_julia2.jl` [M\_julia2.jl](https://discourse.julialang.org/uploads/short-url/zUAugcPIoIS0hz6l05DFYWtjQpI.jl) (954.1 KB), which is also attached to this reply. Julia parser had a lot of problem parsing its own generated code [Parsing Problem In Julia](https://discourse.julialang.org/t/symbolic-integration-for-hard-expressions-and-generate-function/124744/8). That is why I am loading the file in this way.
I’m working on an optimization problem, where I need to enforce the dynamic constraint of my system as equality constraints + some other constraint. In order to do that I need to calculate the inertia matrix + some other terms. The M function here calculates only an element of inertia matrix [the first 3-3 elements of a 9-9 matrix], which is a function of state of the system or decision variables which I defined here as qq.
I started learning Julia and Jump with the hope that I can use it to speed up my optimization and leverage some of the solvers that developed in Julia such as complementarity here. But unfortunately, I am stuck generating my terms to use them in the optimization. I am continuing my work in MATLAB again, but I am hoping that by solving this problem I can use Julia and JuMP.
Regarding eval, I agree that it might not be the best option, I tried to use build\_function too without evaluating it, and then use jump but still cannot make the code output anything when a symbolic vector is passed.
```julia
using Symbolics
Symbolics.@variables xi r L0 q[1:3]
M_julia = read("M_julia2.jl", String)
function parsed_data(expr)
expr_fixed = replace(expr, r"\.([/*+^])" => s"\1")
expr_parsed = Meta.parse(expr_fixed)
return expr_parsed
end
M = parsed_data(M_julia)
@eval @fastmath const M_func(xi,r,L0,q) = $M
function TermsReduced(q)
L0 = 0.15
xi = 1
r = 24e-3
M = M_func(xi, r, L0, q)
return M
end
q1 = reshape([0.05, 0.03, 0.01], (3, 1))
M = TermsReduced(q1)
@time TermsReduced(q1)
println("numerical calculation is finished")
using Complementarity
using JuMP
m = MCPModel()
@variable(m, qq[1:3,1:10], start = 0.0)
println("Using JuMP variables in the function")
@time TermsReduced(qq[:,1])
```
---
### Author:  [@odow](https://discourse.julialang.org/u/odow)
#### Post date: [February 18, 2025, 11:42pm UTC](https://discourse.julialang.org/t/difficulty-in-using-jump-with-complex-functions/126045/4 "2025-02-18T23:42:04Z")
I didn’t run your code, but consider something like:
```Julia
using JuMP
model = Model()
@variable(model, qq[1:9])
@operator(model, op_Terms, 9, (x...) -> Terms(collect(x)))
@expression(model, op_Terms(qq...))
```
Note that JuMP now supports complementarity problems natively. You don’t need to use Complementarity.jl. See [Example: mixed complementarity problems · JuMP](https://jump.dev/JuMP.jl/stable/tutorials/nonlinear/complementarity/)