  The Fokker-Planck (FP) equation is a foundational PDE in stochastic
processes. However, curse of dimensionality (CoD) poses challenge when dealing
with high-dimensional FP PDEs. Although Monte Carlo and vanilla
Physics-Informed Neural Networks (PINNs) have shown the potential to tackle
CoD, both methods exhibit numerical errors in high dimensions when dealing with
the probability density function (PDF) associated with Brownian motion. The
point-wise PDF values tend to decrease exponentially as dimension increases,
surpassing the precision of numerical simulations and resulting in substantial
errors. Moreover, due to its massive sampling, Monte Carlo fails to offer fast
sampling. Modeling the logarithm likelihood (LL) via vanilla PINNs transforms
the FP equation into a difficult HJB equation, whose error grows rapidly with
dimension. To this end, we propose a novel approach utilizing a score-based
solver to fit the score function in SDEs. The score function, defined as the
gradient of the LL, plays a fundamental role in inferring LL and PDF and
enables fast SDE sampling. Three fitting methods, Score Matching (SM), Sliced
SM (SSM), and Score-PINN, are introduced. The proposed score-based SDE solver
operates in two stages: first, employing SM, SSM, or Score-PINN to acquire the
score; and second, solving the LL via an ODE using the obtained score.
Comparative evaluations across these methods showcase varying trade-offs. The
proposed method is evaluated across diverse SDEs, including anisotropic OU
processes, geometric Brownian, and Brownian with varying eigenspace. We also
test various distributions, including Gaussian, Log-normal, Laplace, and
Cauchy. The numerical results demonstrate the score-based SDE solver's
stability, speed, and performance across different settings, solidifying its
potential as a solution to CoD for high-dimensional FP equations.
