My goal is to find the local minimums of a smooth multivariate function $f_t(\vec y)$ for multiple values of $t$. I have created an algorithm $foo(t,\vec x)$ which returns the results of Newton's Method (by "results", I mean that $foo(t,\vec x)$ returns $\vec a$, where $\vec a$ is a local minimum for $f_t(\vec y)$) on $f_t(\vec y)$, where the initial guess is $\vec x$. I have 1,000 values of $t$ to solve for, and for each $t$, I try half-million values of $\vec x$. I then push all the return values of $foo(t,\vec x)$ into an array.
Essentially, I need to run $foo(t,\vec x)$ about half-billion times for different $t$ and $\vec x$. Currently, I have split the computations into four threads, each running $foo(t,\vec x)$ about 125,000,000 times. On my laptop this takes around 2 hours. I would like to shorten the running time as much as possible. What can I possibly do to achieve this?
Here are some details about $f_t(\vec y)$ and $foo(t,\vec x)$:
$f_t: \mathbb R^4 \to \mathbb R$, so this is a very low dimensional optimization problem.
$f_t$ is a radical function.
$foo(t,\vec x)$ uses a standard Newton's Method, i.e. tries to find the roots of the gradient $\nabla f_t$. It finds the descent direction $\vec p_{n+1} = -H_n^{-1}\nabla f_t$, where $H_n$ is the hessian matrix, and using inexact line search finds descent length $\alpha$. Then it sets $x_{n+1}=x_n+\alpha p_n$.
I have heard that machine learning uses GPU. Would my program also be able to utilize the GPU in some way? Maybe run $foo(t,\vec x)$ on each core of a GPU?
