I implemented this in Python: multivar_horner
You can look at the approach used there and port it to Fortran.
NOTE: ToIn contrast to the best of my knowledgeone dimensional case there is no known algorithm for finding the optimal horner factorisationare multiple possible Horner factorisations of a multivariate polynomialpolynomials. The best oneOne can do isallow a search over the possible factorisations to find a minimal representation as described HERE. In most cases however it suffices to use a clever heuristic to find a single "good" factorisationsfactorisation. I implemented amultivar_horner implements the greedy heuristic similar to the one described in the paper "Greedy Algorithms for Optimizing Multivariate Horner Schemes".
The authors of some related publications (including the above mentioned) claim to have an implementation of their proposed algorithms, but I was not able to find any publicly available ones.