# Use ModelingToolkit.jl to solve a maximization problem

**URL:** https://discourse.julialang.org/t/use-modelingtoolkit-jl-to-solve-a-maximization-problem/61042
**Category:** New to Julia
**Tags:** question
**Created:** [May 12, 2021, 7:53pm UTC](https://discourse.julialang.org/t/use-modelingtoolkit-jl-to-solve-a-maximization-problem/61042 "2021-05-12T19:53:11Z")
**Posts on this page:** 7
**Page:** 1

<div class="post-metadata">

### Author: ![Chen](https://avatars.discourse-cdn.com/v4/letter/c/3d9bf3/32.png) [@Chen](https://discourse.julialang.org/u/Chen)
#### Post date: [May 12, 2021, 7:53pm UTC](https://discourse.julialang.org/t/use-modelingtoolkit-jl-to-solve-a-maximization-problem/61042/1 "2021-05-12T19:53:11Z")

</div>

I browsed the website of ModelingToolkit.jl and did not find an example that covers a very simple math problem.

Let’s say that I’d like to maximize y by choosing a value of x:

y = a - b x^2

where a = c + dx - e.

The math is simple. The solution is x^\* = \frac{d}{2b} where y is maximized.

I understand that using SymPy.jl is probably the more common approach. However, I’d like to explore ModelingToolkit.jl for its potentially better speed. I’d appreciate if someone could provide an example code to the simple optimization problem above.

---

<div class="post-metadata">

### Author: ![j\_verzani](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/j_verzani/32/8551_2.png) [@j\_verzani](https://discourse.julialang.org/u/j_verzani)
#### Post date: [May 12, 2021, 8:14pm UTC](https://discourse.julialang.org/t/use-modelingtoolkit-jl-to-solve-a-maximization-problem/61042/2 "2021-05-12T20:14:47Z")

</div>

For this problem, there isn’t much difference in the approach between Symbolics and Sympy, as the derivative is a linear function in the unknown, and `Symbolics.solve_for` (used below handles those with ease.)

```julia
julia> using Symbolics

julia> @variables a b c x
4-element Vector{Num}:
 a
 b
 c
 x

julia> ex = a*x^2 + b*x + c
c + a*(x^2) + b*x

julia> ex′ = Symbolics.derivative(ex, x)
b + 2a*x

julia> Symbolics.solve_for(ex′ ~ 0, x)
(-1//2)*b*(a^-1)

```

Whereas, you might see:

```julia
julia> using SymPy

julia> @syms a b c x
(a, b, c, x)

julia> ex = a*x^2 + b*x + c
   2          
a⋅x + b⋅x + c

julia> ex′ = diff(ex, x)
2⋅a⋅x + b

julia> solve(ex′, x)
1-element Vector{Sym}:
 -b/(2*a)

```

or using `solveset`:

```julia
julia> u = solveset(ex′, x)
⎧-b ⎫
⎨───⎬
⎩2⋅a⎭

julia> elements(u)
1-element Vector{Sym}:
 -b/(2*a)

```

I’m not aware of a direct optimizer in either, but would be happy to pointed to one.

---

<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: [May 13, 2021, 1:17am UTC](https://discourse.julialang.org/t/use-modelingtoolkit-jl-to-solve-a-maximization-problem/61042/3 "2021-05-13T01:17:10Z")

</div>

That’s the symbolic way, and then symbolic-numeric:

```julia
using ModelingToolkit, GalacticOptim

@variables x y
@parameters b c d e
a = c + d*x - e
loss = a - b*x^2
sys = OptimizationSystem(-loss,[x],[b,c,d,e])

u0 = [
    x=>0.0
]
p = [
    b => 1.2
    c => 0.2
    d => 0.3
    e => 0.5
]

prob = OptimizationProblem(sys,u0,p,grad=true)
sol = solve(prob,BFGS())
sol.u # [0.12499999999999999]

```

---

<div class="post-metadata">

### Author: ![Chen](https://avatars.discourse-cdn.com/v4/letter/c/3d9bf3/32.png) [@Chen](https://discourse.julialang.org/u/Chen)
#### Post date: [May 13, 2021, 9:34pm UTC](https://discourse.julialang.org/t/use-modelingtoolkit-jl-to-solve-a-maximization-problem/61042/4 "2021-05-13T21:34:34Z")

</div>

Thank you very much! At the moment I focus more on the symbolic math. The plan is to run numeric simulation later on, probably in a few weeks (or a few months or years, if I turn out to progress very slowly). Your example codes will be very helpful once I started numerical simulation, I think!

---

<div class="post-metadata">

### Author: ![Chen](https://avatars.discourse-cdn.com/v4/letter/c/3d9bf3/32.png) [@Chen](https://discourse.julialang.org/u/Chen)
#### Post date: [May 13, 2021, 9:40pm UTC](https://discourse.julialang.org/t/use-modelingtoolkit-jl-to-solve-a-maximization-problem/61042/5 "2021-05-13T21:40:48Z")

</div>

Thank you very much!

I see that your codes, there was one equation, `ex = a*x^2 + b*x + c`, for example.

In the example of y=a−bx^2, where a=c+dx−e, I wonder if it is possible to input two equations in Julia, and ask Julia to solve it?

Of course I could substitute a into y myself, and feed the resulted equation to Julia, using your example codes. My interest is to ask Julia to solve more complex problems where there are more equations involved. That’s why I’m interested in learning how I can feed more than one equations to Julia, and ask Julia to solve them.

---

<div class="post-metadata">

### Author: ![j\_verzani](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/j_verzani/32/8551_2.png) [@j\_verzani](https://discourse.julialang.org/u/j_verzani)
#### Post date: [May 13, 2021, 10:07pm UTC](https://discourse.julialang.org/t/use-modelingtoolkit-jl-to-solve-a-maximization-problem/61042/6 "2021-05-13T22:07:55Z")

</div>

Sure, both Symbolics and SymPy allow for vectors of equations. Symbolics only has linear equations for now.

An example that would work for either is:

```julia
julia> using SymPy

julia> @syms a[1:2,1:2] b[1:2] x y
(Sym[a₁_₁ a₁_₂; a₂_₁ a₂_₂], Sym[b₁, b₂], x, y)

julia> eqs = [Eq(sum(a[i,:] .* (x,y)), b[i]) for i ∈ 1:2]
2-element Vector{Sym}:
 a₁_₁⋅x + a₁_₂⋅y = b₁
 a₂_₁⋅x + a₂_₂⋅y = b₂

julia> solve(eqs, [x,y])
Dict{Any, Any} with 2 entries:
  y => (a₁_₁*b₂ - a₂_₁*b₁)/(a₁_₁*a₂_₂ - a₁_₂*a₂_₁)
  x => (-a₁_₂*b₂ + a₂_₂*b₁)/(a₁_₁*a₂_₂ - a₁_₂*a₂_₁)

```

Or

```julia
julia> using Symbolics

julia> @variables a[1:2,1:2] b[1:2] x y
4-element Vector{Any}:
  Num[a₁ˏ₁ a₁ˏ₂; a₂ˏ₁ a₂ˏ₂]
  Num[b₁, b₂]
 x
 y

julia> eqs = [sum(a[i,:] .* (x,y)) ~ b[i] for i ∈ 1:2]
2-element Vector{Equation}:
 a₁ˏ₁*x + a₁ˏ₂*y ~ b₁
 a₂ˏ₁*x + a₂ˏ₂*y ~ b₂

julia> Symbolics.solve_for(eqs, [x,y])
2-element Vector{Num}:
 -*(a₁ˏ₁^-1)*(a₁ˏ₂*(a₂ˏ₁*b₁*(a₁ˏ₁^-1) - b₂)*((a₁ˏ₂*a₂ˏ₁*(a₁ˏ₁^-1) - a₂ˏ₂)^-1) - b₁)
                         (a₂ˏ₁*b₁*(a₁ˏ₁^-1) - b₂)*((a₁ˏ₂*a₂ˏ₁*(a₁ˏ₁^-1) - a₂ˏ₂)^-1)

```

---

<div class="post-metadata">

### Author: ![Chen](https://avatars.discourse-cdn.com/v4/letter/c/3d9bf3/32.png) [@Chen](https://discourse.julialang.org/u/Chen)
#### Post date: [May 15, 2021, 7:23pm UTC](https://discourse.julialang.org/t/use-modelingtoolkit-jl-to-solve-a-maximization-problem/61042/7 "2021-05-15T19:23:08Z")

</div>

Thank you!

In general, I have been more interested in packages that are written in Julia, because I think such packages are going to have more speed advantages. However, does it mean that SymPy is at the moment a more capable package than Symbolics, for example, because Symbolics only has linear equations for now?

A side note is that the results from SymPy are much more readable than from Symbolics.
