standupmaths / bubble-cosh Goto Github PK

So I just had to get involved here... I couldn't leave the chaos I saw untouched.

The reason I'm posting this issue (which really is just a message, do with it what you see fit) is because I needed to point out that you're coding in Python! Always assume something is done for you already! And that something is SciPy which has a huge number of different optimization functions (and yes I'm spelling that optimization despite being English... I spotted that maths, but programming is written in "American", no matter how much I disagree)

I'm not even 100% sure what the code is minimizing! I just recognised that it was trying to do that. So below is the same code modified to use the scipy library to solve the problem with minimal effort.

Note some subtle changes:

1. error now takes one argument which is a list of both values. That's just how optimize works.
2. The error function returns the sum of squares rather than the absolutes. You get best results when your function is differentiable within the bounds that you care about (because ScyPy tries to be smart rather than brute force things). There are ways to handle the other situation, but this works better.

import math as maths
from scipy.optimize import minimize

d = 1.068
l = 0.6

def total_area(a, l):
    return maths.pi * (a ** 2) * (maths.sinh(l / a) + (l / a))

def error(args):
    a, b = args
    try:
        x1 = 0
        y1 = d / 2
        e1 = a * maths.cosh((x1 - b) / a) - y1
        x2 = l
        y2 = d / 2
        e2 = a * maths.cosh((x2 - b) / a) - y2
        error = e1**2 + e2**2
        return error
    except:
        return float("inf")

# Use SciPy to minimize our error function with minimal effort
# Note: The (1, 1) represents the starting values for the search
results = minimize(error, (1, 1))

# Extract out the two final results
a, b = results.x

print("for diameters of {0} and length of {1}".format(d, l))
print("a={0}, b={1}".format(a, b))

mid_radius = a * maths.cosh(((l / 2) - b) / a)

print("Area of {0}".format(total_area(a, l)))
print("mid dip of {0}".format((d / 2) - mid_radius))
print("mid gap of {0}".format(mid_radius * 2))

For reference, this gives an area 99.999999% of the original code. So... it's totally more accurate.

And to be clear, I have other issues with coding style... but I can forgive those of mathematicians. I didn't want to change the code too much in a random issue being posted on GitHub.