Can't figure out how to create multivariate approximation

[baggepinnen]

Hi guys, I'm trying to figure out how to create a multivariate approximation to a nonlinear function dfx, but can't quite understand the multivariate interface. Below is an attempt at creating an order 7 approximation, but the two nonlinear equations get approximated linearly only. What am I doing wrong?

using TaylorModels
function dxf(xs, us)
    [
    xs[3]
    xs[4]
       -0.8333333333333334((0.1(0.9810000000000001sin(xs[2]) - 0.08333333333333334(-us[1] - 0.1(xs[4]^2)*sin(xs[2]))*cos(xs[2]))*cos(xs[2])) / (0.008333333333333335(cos(xs[2])^2) - 0.05) - us[1] - 0.1(xs[4]^2)*sin(xs[2]))
       (0.9810000000000001sin(xs[2]) - 0.08333333333333334(-us[1] - 0.1(xs[4]^2)*sin(xs[2]))*cos(xs[2])) / (0.008333333333333335(cos(xs[2])^2) - 0.05)
    ]
end

nx = 4; nu = 1
taylor_expand(zeros(nx+nu); order=7) do xu
    xs = xu[1:4]
    us = xu[5]
    dxf(xs, u)
end

4-element Vector{TaylorN{Float64}}:
                          1.0 x₃ + 𝒪(‖x‖⁸)
                          1.0 x₄ + 𝒪(‖x‖⁸)
  1.9620000000000002 x₂ + 1.0 x₅ + 𝒪(‖x‖²)
            - 23.544 x₂ - 2.0 x₅ + 𝒪(‖x‖²)

[PerezHz]

I think you can make your example work if you use set_variables before calling taylor_expand:

set_variables("x", order=7, numvars=4) # numvars is just the number of variables wrt which you're expanding

Probably it makes sense to add a comment about this in the docs?

[baggepinnen]

Thanks for the help, the following does work

function dxf(xt, ut)
    [
    xt[3]
    xt[4]
       -0.8333333333333334((0.1(0.9810000000000001sin(xt[2]) - 0.08333333333333334(-ut[1] - 0.1(xt[4]^2)*sin(xt[2]))*cos(xt[2]))*cos(xt[2])) / (0.008333333333333335(cos(xt[2])^2) - 0.05) - ut[1] - 0.1(xt[4]^2)*sin(xt[2]))
       (0.9810000000000001sin(xt[2]) - 0.08333333333333334(-ut[1] - 0.1(xt[4]^2)*sin(xt[2]))*cos(xt[2])) / (0.008333333333333335(cos(xt[2])^2) - 0.05)
    ]
end

nx = 4; nu = 1
set_variables("x", order=7, numvars=nx+nu)
taylor_expand(zeros(nx+nu); order=7) do xu
    xs = xu[1:4]
    us = xu[5]
    dxf(xs, us)
end

Is there an interface available that does not mutate a global state, like I assume set_variables does?

---

EDIT: Actually I had a mistake in the code above, it still does not work:

function dxf(xt, ut)
    [
    xt[3]
    xt[4]
       -0.8333333333333334((0.1(0.9810000000000001sin(xt[2]) - 0.08333333333333334(-ut[1] - 0.1(xt[4]^2)*sin(xt[2]))*cos(xt[2]))*cos(xt[2])) / (0.008333333333333335(cos(xt[2])^2) - 0.05) - ut[1] - 0.1(xt[4]^2)*sin(xt[2]))
       (0.9810000000000001sin(xt[2]) - 0.08333333333333334(-ut[1] - 0.1(xt[4]^2)*sin(xt[2]))*cos(xt[2])) / (0.008333333333333335(cos(xt[2])^2) - 0.05)
    ]
end

nx = 4; nu = 1
xu = set_variables("x", order=7, numvars=nx+nu)
xs = xu[1:4]
us = xu[5]
dxf(xs, us)

julia> dxf(xs, us)
4-element Vector{TaylorN{Float64}}:
                          1.0 x₃ + 𝒪(‖x‖⁸)
                          1.0 x₄ + 𝒪(‖x‖⁸)
  1.9620000000000002 x₂ + 1.0 x₅ + 𝒪(‖x‖²)
            - 23.544 x₂ - 2.0 x₅ + 𝒪(‖x‖²)

[baggepinnen]

Actually, the problem in this case might be that the multivariate Taylor expansion is the wrong tool for this job, and what I really should be looking at is a least-squares approximation over an interval instead

[PerezHz]

I think the problem is in the definition of dxf you're doing ut[1] and since ut is a TaylorN, then you're actually taking the 1st-order coefficient of its Taylor expansion; so if you instead do:

nx = 4; nu = 1
xu = set_variables("x", order=7, numvars=nx+nu)
xs = xu[1:4]
us = xu[5:5] # construct us as a 1-element Vector{TaylorN{Float64}}
dxf(xs, us)

then your example should work. You can also redefine dxf so that you read ut instead of ut[1]. I don't know the specific problem that you're trying to solve, but if for example you're trying to minimize a loss function then in general, if you have the 2-nd order Taylor expansion wrt your variables, it should be straight-forward to perform Newton method iterations.

[PerezHz]

To answer your other question, TaylorSeries.jl mutates a global state (hash-tables) to construct TaylorN variables...

[baggepinnen]

I think the problem is in the definition of dxf you're doing ut[1] and since ut is a TaylorN

Thanks, I just figured this out as well (after a lot of trial and error). I opened #331 about this since it's kind-of a violation of the Number interface that can cause all sorts of unexpected bugs if you pass TaylorN into a function into which you have no insight

[lbenet]

Is there an interface available that does not mutate a global state, like I assume set_variables does?

The way TaylorSeries works, for many variables, is by mutating a global state. Getting to the details, the many-variables polynomial is stored as a vector of coefficients which correspond to the homogeneous polynomials. The order of each monomial within such a homogeneous polynomial is fixed and defined by the number of variables and the order (degree). Changing those then changes the ordering. That's the naive reason for having chosen this design.