On 26 January 2015 at 19:42, Chris Barker <chris.barker@noaa.gov> wrote:
OK -- here is a real-life application
[details snipped]
The simplist iterative solution is something like this: recast the equation as:
k_2 = omega**2 / (g * tanh(k_1 * h))
- guess a k_1 - compute k_2 - check if k_2 is close to k_1 - if not, set k_1 to k_2 repeat.
[more details snipped]
In fact, you write a generic iteration routine:
def iterate(func, x_initial, *args): """ iterate to find a solution to the function passed in
func should be a function that takes x as a first argument, and computes an approximation to x as a result.
x_initial is an initial guess for the unknown value """ x_1 = x_initial while True: x_2 = func(x_1, *args) if is_close_to(x_2, x_1): break x_1 = x_2 return x_2
The point I was trying to make earlier (which came across as "nobody should be coding their own Newton iterations") is essentially this - surely iterating a relation until the value you need converges is a common operation, and libraries exist that do this for you? And those libraries would already have their own "is the result close enough" calculation (possibly tunable via arguments to the routine). Do people actually write the low-level algorithms by hand, so that a building block like is_close is worthwhile? The lack of people coming up with use cases suggests to me that maybe it isn't. FWIW, a quick google search came up with scipy.optimize.fixed_point, which is, AFAICT, the equivalent of your iterate(). (You probably know better than I do if that's the case). Paul