This paper is based on a formulation of the Navier-Stokes equations developed by P. Constantin and the first author (texttt{arxiv:math.PR/0511067}, to appear), where the velocity field of a viscous incompressible fluid is written as the expected value of a stochastic process. In this paper, we take $N$ copies of the above process (each based on independent Wiener processes), and replace the expected value with $frac{1}{N}$ times the sum over these $N$ copies. (We remark that our formulation requires one to keep track of $N$ stochastic flows of diffeomorphisms, and not just the motion of $N$ particles.) We prove that in two dimensions, this system of interacting diffeomorphisms has (time) global solutions with initial data in the space $holderspace{1}{alpha}$ which consists of differentiable functions whose first derivative is $alpha$ Holder continuous (see Section ref{sGexist} for the precise definition). Further, we show that as $N to infty$ the system converges to the solution of Navier-Stokes equations on any finite interval $[0,T]$. However for fixed $N$, we prove that this system retains roughly $O(frac{1}{N})$ times its original energy as $t to infty$. Hence the limit $N to infty$ and $Tto infty$ do not commute. For general flows, we only provide a lower bound to this effect. In the special case of shear flows, we compute the behaviour as $t to infty$ explicitly.