# Large-scale least squares with ModelingToolkit.jl?

**URL:** <https://discourse.julialang.org/t/large-scale-least-squares-with-modelingtoolkit-jl/67957>\
**Category:** Modelling & Simulations\
**Tags:** sciml, modelingtoolkit\
**Created:** [September 9, 2021, 8:28pm UTC](https://discourse.julialang.org/t/large-scale-least-squares-with-modelingtoolkit-jl/67957 "2021-09-09T20:28:46Z")\
**Posts on this page:** 2\
**Page:** 1

<div class="post-metadata">

**Author:** ![twadleigh](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/twadleigh/32/22778_2.png) [@twadleigh](https://discourse.julialang.org/u/twadleigh)\
**Post date:** [September 9, 2021, 8:28pm UTC](https://discourse.julialang.org/t/large-scale-least-squares-with-modelingtoolkit-jl/67957/1 "2021-09-09T20:28:46Z")

</div>

I’m starting to dig deep into `ModelingToolkit.jl` with the intention of using it to model a version of the [SLAM](https://en.wikipedia.org/wiki/Simultaneous_localization_and_mapping) problem. Ultimately, SLAM problems end up looking like large-scale nonlinear least squares problems. Instead of working with a single scalar objective function, large-scale least squares solvers exploit the extra structure in the problem by working with the residual function and its Jacobian directly.

`ModelingToolkit.jl` provides `OptimizationSystem`, but it is designed to wrap a single scalar objective function, and doesn’t provide a method for `generate_jacobian`.

`NonlinearSystem` seems like almost the right thing to use to build up the residual function by components (including having a `generate_jacobian` method) but requires the system be expressed as a vector of `Equation`s, rather than simple expressions.

My current thought is to just make use of `NonlinearSystem` for now, prepending my residual expressions with `0 ~` to turn them into equations.

A few questions, though:

1. Immediate term: will this hack work? I think so, but I haven’t tried yet.
2. Slightly longer term, what is the right way to make things less hacky for this use case? A little generalization of `NonlinearSystem`? A new concrete sub-type of `AbstractTimeIndependentSystem`?
3. Even longer term, where might large-scale least squares solvers for problems of this sort find their way into the ecosystem? Once I’ve got my residual function modeled, I’m perfectly comfortable cobbling together the solver bits I require from elsewhere in the julia ecosystem and beyond, but it would be nice to button things up with a more complete story fully contemplated by the reach of SciML.

I opened an [issue on github](https://github.com/SciML/ModelingToolkit.jl/issues/1256) about this, figuring it would ultimately lead to some opportunity to bring relevant extensions to `ModelingToolkit.jl` and/or packages in its larger ecosystem, but figured it would probably be best to cross post here, since at this point it is more about questions and discussions than specific software changes.

---

<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:** [September 11, 2021, 2:15am UTC](https://discourse.julialang.org/t/large-scale-least-squares-with-modelingtoolkit-jl/67957/2 "2021-09-11T02:15:07Z")

</div>

> [@twadleigh](#):
>
> `ModelingToolkit.jl` provides `OptimizationSystem` , but it is designed to wrap a single scalar objective function, and doesn’t provide a method for `generate_jacobian` .

Well the gradient is the one-row Jacobian of a scalar equation.

> [@twadleigh](#):
>
> My current thought is to just make use of `NonlinearSystem` for now, prepending my residual expressions with `0 ~` to turn them into equations.

Residual equations are the same thing, and sometimes can be a little bit more naturally written as an equation! Indeed, NonlinearSystem automatically moves them to the right and drops the zeros internally.

> [@twadleigh](#):
>
> Slightly longer term, what is the right way to make things less hacky for this use case? A little generalization of `NonlinearSystem` ? A new concrete sub-type of `AbstractTimeIndependentSystem` ?

If there are solvers which specialize on this, then it would make sense to do. Right now there are no solvers in the SciML ecosystem to target. See [https://scimlbase.sciml.ai/dev/](https://scimlbase.sciml.ai/dev/) . So it really needs to happen the other way: first there needs to be a set of specific least squares solvers and a LeastSquaresProblem, and then a symbolic form would make sense.