Global Optimization of a well-defined function with gradient information
One alternative is to use Newton's method.  Newton's method uses information from the first and second derivative (i.e., the Hessian matrix), so it can converge in fewer steps than gradient descent (which uses only information from the first derivative) -- but the cost is that each single step is slower, because it requires more complex computations.

If Newton's method is too slow, there are faster variants, such as BFGS and other quasi-Newton methods or the conjugate gradient method.

Another approach is to continue using gradient descent, but make some adjustments.  One simple thing you can try is to apply it many times, each time with a different random starting point.  Also there are various tuning knobs you can tweak, such as adjusting the step size or using momentum.  A third possibility could be to try to make a "smarter" choice of starting point for gradient descent: you can approximate each cosine by a quadratic function (using the Taylor series expansion at a particular point), obtain a function $f^*$ that is a quadratic approximation to $f$, and then find the exact minimum for $f^*$ using least-squares method.  This basically amounts to applying one step of Newton's method and then continuing with gradient descent from there.

Trying several different methods shouldn't be too hard if you use the right software, as there are many software packages that let you specify the function and then make it easy to try multiple different optimization methods.