get linear terms

[baggepinnen]

The following code produces an expression that contains terms

0.5(u2[1] - u1[1])*(1.212532u2[1] - 1.212532u1[1])

from which I would like to be able to extract the coefficient of u[1]^1, which is alpha*u2[1] for some value of alpha. linear_expansion gives me

julia> Symbolics.linear_expansion(L, U[1]) # U[1] == u[1]
(0, 0, false)

using Symbolics
function vvariable(name, length)
    un = Symbol(name)
    u = @variables $un[1:length]
    collect(u[])
    # u[]
end
function loss(x, u, ulast, n)
    c = 0.5*(dot(x, Q1, x) + dot(u, Q2, u))
    if n > 1
        du = u - ulast
        c +=  0.5dot(du, Q3, du)
    end
    c
end
function final_cost(x)
    0.5x'QN3*x # TODO: replace by Riccati solution
end



nx = 1
nu = 1
Ts = 1
N = 3
Q1 = randn(nx,nx)
Q2 = randn(nu,nu)
Q3 = randn(nu,nu)
QN3 = randn(nx,nx)

X   = Num[] # variables
U   = Num[] # variables
L = 0
x = vvariable("x1", nx) # initial value variable

append!(X, x)

dynamics(x,u) = x

u0 = [0]
u = u0#vvariable("u0", nu)
n = 1
for n = 1:N # for whole time horizon N
    global x, L, u
    ulast = u
    u = vvariable("u$n", nu)
    append!(U, u)
    xp = dynamics(x, u)
    L += loss(x, u, ulast, n) # discrete integration of loss instead of rk4 integration, this makes the hessian constant
    
    x = vvariable("x$(n+1)", nx) # x in next time point
    append!(X, x)
end
L += final_cost(x)

Symbolics.linear_expansion(L, U[1])

[YingboMa]

- 0.7708555204507617u1[1]*u2[1] is not linear

[baggepinnen]

Symbolics.semilinear_form([L], [U[1]])

does what I want