Where we left off, I stated that the goal is:
`d/(d\theta)L_"total" = 0`
Now all that's left is to find an expression for `L_"total"` and solve for `r(theta)`, right?
Well the first part, yes, so let's tackle that first.
I split up the problem into 3 separate expressions for each side of the winch: `L_"winch"`, `L_"middle"`, and `L_"arc"`. This means that `L_"total"` is equal to:
`L_"winch1" + L_"middle1" + L_"arc1" + L_"constant" + L_"arc2" + L_"middle2" + L_"winch2"`
`L_"winch"` Formulation
The easiest way that I found to formulate `L_"winch"` is to use parametric equations, where the parameter is theta (rotation of the winch). Here is the basic set of parametric equations describing the xyz location of a point on the winch helix:
`vec w = ((r(theta)cos(theta)), (r(theta)sin(theta)), (Ptheta))`
where:
`r(theta)` is a function describing a variable radius in terms of theta
`P` is the pitch of the helix.
The length of a parametric equation can be found with an arc integral:
`L_"winch" = int_0^(theta)||vec w'(theta)|| d\theta`
We will use this equation to verify the solution at the end, but our problem statement only requires the derivative of the winch length, so we can simplify the expression we need into the L2 norm of the derivative of the winch parametric equation vector:
`d/(d\theta)L_"winch" = ||vec w'(theta)||`
This form of the expression is critical, because it means we don't need to solve for the integral, which turns out to not have a closed-form solution for most forms of `r(theta)` (except maybe with elliptic integrals for some simpler forms of `r(theta)`).
I decided to use SymPy for the algebraic manipulation and analytical derivatives, to simplify the code and to allow for `r(theta)` to be changed easily.
`L_"middle"` Formulation
At least for me, formulating `L_"middle"` was a bit easier, since it's just a matter of doing some trigonometry. I'm sure there is a more elegant way of finding `L_"middle"`, but I decided to just draw some triangles:
Substituting in `Ptheta' for z, we get:
`L_"middle" = sqrt((Ptheta - r_"pulley")^2 + (H - r(theta))^2 + r_"pulley"^2)`
where:
`H` is a constant for the height between the center of the winch axis and the center of the pulley.
`r_"pulley"` is a constant for the radius of the pulley.
Again, I then used SymPy to take the analytical derivative after substituting in an expression for `r(theta)`.
`L_"arc"` Formulation
In the previous figures, I showed another angle, `phi`. Given this angle, it is trivial to find the length of the arc on the pulley.
`phi = tan^-1(z/(H - r(theta)))`
Substituting in `Ptheta` for z, we get:
`phi = tan^-1(Ptheta/(H - r(theta)))`
The length of the arc is just equal to `phi` times `r_"pulley"`:
`L_"arc = "r_"pulley"tan^-1(Ptheta/(H - r(theta)))`
Solving for `r(theta)`
Ok so that's the first part out of the way. All the equations are defined. The issue is that now we have a pretty complex differential equation that we need to solve. Most differential equations simply don't have closed form solutions. There's no way to know for sure, but I'm fairly certain there is no closed form solution to this cable length equation. That means we need to use nonlinear optimization to find the best equation for `r(theta)`, to get `L_"total"` as close to zero as possible.
I started out by using SciPy's optimize.least_squares and optimize.minimize methods. My first implementation had some math errors, so I couldn't get the solvers to converge. In the process of troubleshooting, I switched to CasADi as my modelling and optimization tool. I'm sure one of SciPy's optimization algorithms would've worked fine, but this was a good opportunity for me to learn how to use CasADi anyways. CasADi is very powerful compared to SciPy, since it does automatic differentiation, which drastically improves solve times. For nonlinear programs, CasADi uses the IPOPT solver, which can solve large scale sparse problems, as well as smaller problems.
The optimization problem I've described is to optimize the derivative of the winch length, for all theta within a range. This is an infinitely dimensional optimization problem, since `r(theta)` could be anything. It could be a function with a billion parameters. So in order to actually solve the problem, we need to define `r(theta)` a little better.
This type of optimization problem, where the parameter that can be changed is a function, is most commonly seen in trajectory optimization problems. You can think of this winch problem as a controls problem, where there are no time dynamics, but the control input u is the instantaneous radius of the winch. Instead the differential equation defining the "dynamics" is relative to `theta`. So `theta` is analogous to time. The most common approach to trajectory optimization problems is to use splines to define the function, with some finite number of control points. This website has a lot of great resources about trajectory optimization: http://www.matthewpeterkelly.com/tutorials/trajectoryOptimization/
Spline based optimization is definitely the best option, but I decided to go with a simpler option to start, which is a single polynomial function. Since the solution is likely to be smooth and monotonic, a relatively low-order polynomial seemed like it would work. In the end, I found that a 1st order polynomial (i.e. a line) provided a remarkably good solution. The code allows you to specify any order polynomial. 2nd and 3rd order polynomials also worked well.
So what does our final winch look like?
There is still some residual error, so I decided trying some other polynomials. This is the optimal solution with a 3rd order polynomial:
This time the radius change is much more dramatic.
There's another important parameter that isn't captured in the optimization problem, which is to reduce the amount of radius change over the winch. That parameter is important because as the radius changes, the effective force that a motor can put on the cable changes. That means the motor needs to be oversized to be able to handle all possible transmission ratios. This cubic radius solution is clearly worse than the linear solution, even though it is technically better at maintaining constant tension. Even with the linear radius model, we're well beyond manufacturing tolerances, so any improved performance is not necessary. But regardless, this does perform better at maintaining constant length, with a cable length change of only 0.003mm.
So the final question: Does this actually work in real life?
I don't know quite yet. I'm doing some testing with 3d printed models, which will hopefully be accurate enough to answer that. In the script, the verification portion that produces these length plots is somewhat separated from the optimization portion, but ultimately if my mathematical model is incorrect, the optimization and the verification will be incorrect. So the only real option is to try it on real hardware. Stay tuned for part 3.
Check out the code for this project on my Github: https://github.com/kylebme/winch-optimization