Stochastic PDE with conservative noise

I am trying to simulate a simple stochastic pde using the DifferentialEquations package, the equation I want to simulate is a simple conservation equation with conservative noise in 2 dimensions, i.e.

\partial_t \rho = D \nabla^2 \rho + \xi ; where the correlation of the noise \xi is
<\xi(x,t) \xi(x^\prime,t^\prime)> = -D \nabla^2 \delta(x-x^\prime) \delta (t-t^\prime).

Here the \nabla^2 in the noise correlation function stems from the conservation of the total density.

As a starting point I decided to adopt a very naive strategy, namely define \rho on a fixed N*N lattice, hence (\rho(x,y) is a matrix \rho_{ij}), and compute the Laplacian numerically. Defining an appropriate matrix K, the discrete Laplacian can be written as \rho.K+ K.\rho.

Though, I do not now how include the Laplacian in the correlation of the noise, probably I have to define a new noise dWc = dW.K+K.dW, where dW is matrix of N*N Wiener Processes. However I am not sure how to do that. From the documentation of the DifferentialEquations package it is clear how to implement colored noise or non diagonal noise (namely implementing something like M.dW); but I didn’t find a way to implement something like dWc = dW.K+K.dW. In general what I would like to do is to define a new noise which is a function of the old one dWnew = F(dW), is that possible?

You’d have to define a new AbstractNoiseProcess for this. Maybe @Ricardo_Rosa might have a trick in mind.

I am not sure I follow what you want. Do you want to take the Laplacian, in space, of a Dirac delta ? That would be pretty rough. Can you explain how that stems from the conservation? I think I need some background. I imagine there is a paper you can direct me to?

A good reference for this topic is: “Critical dynamics: a field theory approach to equilibrium and non equilibrium scaling behaviour” sec 3.2 by Uwe Tauber, hope that is not too technical.

I will try to explain what I mean. What I want is for the total density to be conserved in time, which means that
\int \rho(x,t) dx = cost
or in the discrete case

\sum_{i,j=1}^N \rho_{ij}(t) = const.

The SPDE for such a field looks like

\partial_t \rho = D\nabla^2 \rho + \xi

The first term in the RHS is a diffusion term and it does not change the density. The second term is a noise term and depending on its correlation may or not change the total density. If the noise is diagonal and delta correlated the total integral of density changes in time and this is not what I want. Conversely, if the correlation of the noise has a laplacian density is conserved. This can be understood by considering the noise correlation in the discrete case (I apologize for the heavy notation) ,

<\xi_{ij} \xi_{nm}> = D[ 4 \delta_{in} \delta_{jm} - \delta_{i+1,n} \delta_{jm} -\delta_{i-1,n} \delta_{jm} -\delta_{in}\delta_{j+1,m} - \delta_{in} \delta_{j-1,m}]

(I am omitting the temporal structure)
This equation is telling us that if there is a stochastic fluctuation of \rho on site (i,j) an opposite fluctuation must be happening on the neighboring sites, and these fluctuations sum up to zero.
This is the reason why having the laplacian in the correlation of the noise implies conservation.

A possible solution:
A possible way of doing this would be to define a noise term which is the gradient of a vectorial noise. This means that \xi(x,y) = \nabla \cdot \theta(x,y), where \theta is a delta correlated vectorial noise (namely a N\times N\times 2 object, where N is the linear size of the system I am simulating). However, I am not sure how I could implement that.

Sorry, I had not seen you had replied. I didn’t get any message from Discourse.

Anyway, yeah, your “possible solution” makes more sense than the previous one, to me. One variation of that is to multiply the noise by the density, to have something of the form

\partial \rho + \nabla \cdot (\rho \xi) = D\Delta \rho,

where the noise \xi is interpreted as random velocity fluctuations moving the mass around.

This seems like something for which we could embed support in MOL, does MTK have a symbolic model for random processes?

I don’t think it has an RODESystem yet, but it needs to get added soon. It needs an @noise and all of that, so we’re on the case.