I may be mistaken, but it seems to me that you want to initialize your pde equation/system by a stochastic process (1D) or a stochastic field (nD) ?
If it is the case there exist a lot of methods for numerically simulating such objects (most commonly used and easier to understand are based on Cholesky decomposition and Fast Fourier Transform)
Here is just one link among others : http://www.icms.org.uk/downloads/PMPM2014/PowellPart2.pdf
If you are only interested in implementing a simulation method, a very quick overview of the FFT simulation method is provided here http://www.siam.org/students/g2s3/2011/docs/alea.pdf
For more details, see the reference site for gaussian processes http://www.gaussianprocess.org