The response of a structural system to white noise excitation (deltacorrelated) constitutes a Markov vector process whose transitional probability density function (TPDF) is governed by both the forward Fokker-Planck and backward Kolmogorov equations. Numerical solution of these equations by finite element and finite difference methods for dynamical systems of engineering interest has been hindered by the problem of dimensionality. In this paper numerical solution of the stationary and transient form of the Fokker-Planck (FP) equation corresponding to two state nonlinear systems is obtained by standard sequential finite element method (FEM) using C⁰ shape function and Crank-Nicholson time integration scheme. The method is applied to Van-der-Pol and Duffing oscillators providing good agreement between results obtained by it and exact results. An extension of the finite difference discretization scheme developed by Spencer, Bergman and Wojtkiewicz is also presented. This paper presents an extension of the finite difference method for the solution of FP equation up to four dimensions. The difficulties associated in extending these methods to higher dimensional systems are discussed.