Adding noise in frequency domain while solving an SPDE

I’m implementing a numerical solver to a certain kind of a 2D stochastic PDE using Fourier Transforms. So, this solver is a Galerkin type method, which means that choosing the right basis functions I can apply Fourier Transforms to compute the solution.

Long story short, applying the explicit Euler-Maruyama scheme to my SPDE and discretising the domain with $$N^2$$ nodes, the discrete version of the Fourier coefficients of $$V(t_{i+1})$$ is given by: $$$$\tag{1} \hat{v}_{mn}\left(t_{i+1}\right) = \hat{v}_{mn}\left(t_{i}\right) + \frac{\Delta t}{\alpha}\Big(\hat{\imath}_{mn}\left(t_{i}\right) - \hat{v}_{mn}\left(t_{i}\right) + \hat{a}_{mn}(t_i)\Big) + \epsilon \sqrt{\Delta t} \lambda_{mn}W\left(t_i\right),$$$$ where the hat terms corresponds to the Fourier Coefficients of the equation's functions, $$\epsilon$$ is the level of additive noise, $$W$$ is a complex random matrix with distribution $$N(0,1)$$ and $$\lambda$$ is a spatial correlation matrix.

Notice that the noise is being added in the frequency domain!

However, when comparing results with other numerical methods the effect of the noise is not the same, actually I have to scale up $$\epsilon$$ by a factor.

I'm using the fftw library in julia to compute the Fast Fourier Transforms and it is known that this library do gives us the Fourier coefficients multiplied by $$N^2$$, i.e. $$fft(V) = N^2 \hat{v}_{mn}$$. So, the stochastic term of (1) has not the same magnitude as the rest of the coefficients. Am I right?

To overcome this issue, I simply have to multiply $$\epsilon$$ by $$N^2$$, is this easy? Or adding noise in the frequency domain is more complicated than this? Sorry for my question, but I’m a bit newbie in the Fourier world.

Actually, an experimental evidence that the number of nodes influences the noise effect in the solution, is if I increase $$N$$ I'll also have to increase $$\epsilon$$ to continue to have similar behaviours than with less nodes.

Thank you! Tiago