We develop a new method to solve the Fokker-Planck or Kolmogorov's forward equation that governs the time evolution of the joint probability density function of a continuous-time stochastic nonlinear system. Numerical solution of this equation is fundamental for propagating the effect of initial condition, parametric and forcing uncertainties through a nonlinear dynamical system, and has applications encompassing but not limited to forecasting, risk assessment, nonlinear filtering and stochastic control. Our methodology breaks away from the traditional approach of spatial discretization for solving this second-order partial differential equation (PDE), which in general, suffers from the "curse-of-dimensionality". Instead, we numerically solve an infinite dimensional proximal recursion in the space of probability density functions, which is theoretically equivalent to solving the Fokker-Planck-Kolmogorov PDE. We show that the dual formulation along with the introduction of an entropic regularization, leads to a smooth convex optimization problem that can be implemented via suitable block co-ordinate iteration and has fast convergence due to certain contraction property that we establish. This approach enables meshless implementation leading to remarkably fast computation.