Of course, generally the optimization is NP-hard. However there is a couple tricks that can be played in non-convex case.
First, if $d$ is domain dimension, the probability of, say, gradient descent getting stuck in local minimum decreases exponentially as $d \rightarrow \infty$
Thus, the main obstacle are, in this case, saddle points.

To escape the saddle points one may use the modified Newton method, that, instead of taking steps $-\alpha H^{-1}\nabla f$ takes steps of $-\alpha |H|^{-1}\nabla f.$

This is so-called Saddle-free Newton, which requires $O(d^3)$ time to take a step (matrix inversion is slow). One can use a version that approximates the subspace containing eigenvectors corresponding to largest (by absolute value) eigenvalues with technique similar to Lanzchos procedure. Time constraint becomes $O(kd)$ for some constant $k.$

Both methods were intorduced and proved in this article:
https://arxiv.org/abs/1406.2572

As for public implementation, there is one here:
https://github.com/smdrozdov/saddle_free_newton