I have a set of equations where there are some input variables and one output variable. I want to minimize the output variable, and thought to do this with an OptimizationProblem. In this example, distance is the only input, and I want to minimize total_acc. My question is: how can I find the equation of total_acc relative to distance programmatically, so I can pass it to the OptimizationSystem?
# Falling mass, attached to linear spring
using ModelingToolkit, Optimization, OptimizationOptimJL, LinearAlgebra
using ModelingToolkit: t_nounits as t, D_nounits as D
G_EARTH::Vector{Float64} = [0.0, 0.0, -9.81] # gravitational acceleration [m/s²]
L0::Float64 = -10.0 # initial spring length [m]
# defing the model, Z component upwards
@parameters mass=1.0 c_spring=50.0 damping=0.5 l0=L0
@variables pos(t)[1:3]
@variables vel(t)[1:3]
@variables acc(t)[1:3]
@variables unit_vector(t)[1:3]
@variables spring_force(t)[1:3]
@variables norm1(t) spring_vel(t)
@variables total_acc(t)
@variables distance(t)
eqs = [
pos ~ distance * [0.0, 0.0, 1.0]
vel ~ [0.0, 0.0, 0.0]
norm1 ~ norm(pos)
unit_vector ~ -pos / norm1 # direction from point mass to origin
spring_vel ~ -unit_vector ⋅ vel
spring_force ~ (c_spring * (norm1 - abs(l0)) + damping * spring_vel) * unit_vector
acc ~ G_EARTH + spring_force/mass
total_acc ~ norm(acc)
]
eqs = Symbolics.scalarize.(reduce(vcat, Symbolics.scalarize.(eqs)))
@mtkbuild sys = OptimizationSystem(total_acc, [distance], [])
prob = OptimizationProblem(sys, [1.0], [], grad = true, hess = true)
u_opt = solve(prob, GradientDescent())
# solution should be the distance where the acceleration is the smallest