NASA-AMMOS / aerie

A software framework for modeling spacecraft.
https://nasa-ammos.github.io/aerie-docs/
MIT License
73 stars 19 forks source link

Follow up on X-value bounds in the rootfinding algorithm #1490

Open dandelany opened 5 months ago

dandelany commented 5 months ago

Split off from #1488 - a fix to enforce bounds of initial X-values in the rootfinding algorithm. In that PR, there was some discussion about why the is necessary, and whether there is a deeper improvement to be made to the underlying algorithm or heuristics. Let's continue that discussion here so we can merge & release #1488.

from @adrienmaillard:

Clarification about why this happens as I have seen this is not clear. Let's say you want to find the root of f. Your x domain is [1,3], your y domain is [5,6]. You start at x0, provided as an argument of the rootfinding procedure. In our case, it's gonna be 1 (we take the lower bound of the x domain). The algorithm actually needs 2 starting values to really begin. we compute f(x0) = y0, let's say that f(1) = 8. Then we use a heuristic that assumes that the function is constant so we say x1 (the next x to be tested) is x1 = x0 - 8 (as f(x) is assumed to be x+8 and we want f(x) = 0). x1 = -7 which is outside of the x domain [1,3]. This x1 is passed to the nextValueAt as an init value. The domain of values found by nextValueAt was checked... except for the init value. So in the faulty case, it computes f(x1) and returns. f(x1) satisfies constraints on the y domain and rootfinding is done even though x1 is outside of the x domain. The x domain is determined by the application of constraints beforehand, of which mutex is one.

from @Mythicaeda:

Why does the algorithm want the second value to be f(x)=0 instead of f(x_max)=y? If it took f(x_min) and f(x_max) as its first two values, then it shouldn't ever escape the fn's x-domain (and would remove the need for x_init, assuming of course that the caller isn't optimizing this). Even if x_init was kept, having the second point be the x-domain bound it's further from would still give the alg two points while ensuring it's in the bounds (while avoiding issues should x_init be x_min or x_max.

from @adrienmaillard:

Our overall use of rootfinding is to find f(x) = y but we transpose this problem to a traditional rootfinding problem (there is one method with the y argument that transposes to another method without the y... because it's gonna be 0). The goal of traditional rootfinding is to find x such that f(x) = 0. Here we add bounds for x and f(x). The bounds for f(x) are just the wiggle room around 0 that we are willing to go with.

Mythicaeda commented 5 months ago

@adrienmaillard I'm moving the discussion to this issue:

I understand that the goal is to find f(x)=0, but the algorithm is currently assuming that f(x) is monotonically increasing (as in f(x) = x + C) instead of checking if f(x) is a decreasing function or what its slope is in general. If the second value it picks is another actual x val in the domain rather than assuming off a single point the shape of the fn, then nextValueAt shouldn't have to check if it was actually handed a good value.

adrienmaillard commented 5 months ago

No, the algorithm is not assuming anything in general. This x1 is the only one for which we make this assumption. We can't know if it is decreasing or not before having a second point. Actually we can't know anything about the function, that's why we use this specific rootfinding. Otherwise we would used closed form solutions. The check about the second point could be moved out of nextValueAt for the initValue but we need it anyway in nextValueAt in the general case (when it is testing several x and f(x)).