A code for many-source simulation

In recent posts, I have described the creation and dynamics of "boson" in the proposed local-realistic model of quantum mechanics. In this post, we shall present a simulation code that implements these rules, allowing the local-realistic emergence of the QM outcome of a double slit experiment. The scenario is that introduced here.

Written in the Matlab language, the code would like like that:

%%% Simulate an ensemble of Np particles emitted at intervals Nti 
%%% from either of Ns distinct sources xs(1,...,Ns) having 
%%% probability Ps(1,...,Ns). Evaluate the frequency of arrivals at %%% a 'screen' after Nt iterations.
%%%
%%% Parameters: Np,Ns,Nt,Nti,xs,Ps.
%%%
%%% Evaluate the number of possible bosons for this scenario.

B = Ns*(Ns-1)/2;

%%% Attribute a label to each possible boson.
%%% Bosons generated from the same pair of sources 
%%% will have the same label.

lab = NaN(Ns,Ns);
for i = 1:Ns
    lab(i,i+1:Ns) = [(i-1)*(Ns-i/2)+1:i*(Ns-(i+1)/2)];
end
for i = 1:Ns
    lab(i,1:i-1) = lab(1:i-1,i);
end

%%% Pre-allocate lattice boson trace (lam), lifetime (h), and 
%%% momenta (ome). The nodes that are possibly visited are those 
%%% spanning from min(xs)-Nt to max(ns)+Nt.

lam = NaN(2*Nt+1+max(xs)-min(xs),Nt);
h = zeros(2*Nt+1+max(xs)-min(xs),Nt,B);
ome = zeros(2*Nt+1+max(xs)-min(xs),Nt,B);
ome0 = zeros(2*Nt+1+max(xs)-min(xs),Nt,B);

%%% Start simulating
%
for i = 1:Np
   
   %%% Pre-allocate particle boson momenta (w) and lifetime (k)
   
   k = zeros(1,B);
   w = zeros(1,B);
   
   %%% Attribute a random source momentum (v0)
   
   v0 = -1+2*rand;
   
   %%% Choose the source node according to probability vector
   
   x0 = xs(find(rand<cumsum(Ps),1,'first'));
   
   %%% Position the particle to be emitted at the source node x0
   
   x = x0;
   
   %%% Initialize span l
   
   l = 0;
   
   %%% Simulate Nt time steps 

   avegw = 0;
   occur = 0;
   for t = 1:Nt

      %%% Evaluate momentum and energy propensity

      vQ = sum(w);
      V = v0-2*vQ;
      e = (1+V^2)/2;

      %%% Attribute a momentum v

      v1 = [1,0,-1];
      v = v1(find(rand<[(e+V)/2,1-(e-V)/2,1],1,'first'));

      %%% Increase the particle's position s according to momentum

      x = x+v;

      %%% Increase the particle's span

      l = l+v;

      %%% Update PB lifetime and momenta

      k = k+1;
      w = w.*(1-1./(2*k));

      %%% Update LB lifetime and momenta

      h = h+1;
      ome = ome.*(1-(ome0./h).^2);

      %%% Detect Quantum Reset 

      ix = x+Nt+1-min(xs);
      if l ~= lam(ix,t)

          %%% First-time visited node

          if isnan(lam(ix,t))
              lam(ix,t) = l;

          %%% Detect boson type involved

          else
              iB = lab(xs==x-l,xs==x-lam(ix,t));
              if ~isnan(iB)

                  %%% Exchange span with trace
                  aux = l;
                  l = lam(ix,t);
                  lam(ix,t) = aux;

                  %%% Reset momentum and lifetime of PB of type iB

                  k(iB) = 0;
                  w(iB) = ome(ix,t,iB)/abs(l-lam(ix,t));

                  %%% Reset momentum and lifetime of LB of type iB

                  h(ix,t,iB) = 0;
                  ome(ix,t,iB) = abs(l-lam(ix,t))*V;%(l/t);
                  ome0(ix,t,iB) = abs(l-lam(ix,t))*V;%(l/t);
              end
          end
      end
   end

   %%% Record arrival node for further analysis

   X(i) = x;

   %%% Wait until the next emission and update lattice bosons

   for t = 1:Nti-Nt
      h = h+1;
      ome = ome.*(1-(ome0./h).^2);
   end
end

When run with the following parameters: Np = 300000, Ns = 2, xs = [-1,1], Ps = [0.5,0.5], Nt = 300, Nti = 0, the code takes about 3 hours to run on a scientific computer. Plotting the no. of arrivals per node by keeping only the last 5000 particles, we obtain a plot similar to the following:

Clearly, the interference pattern predicted by QM (in red) is about to being built. Arrivals are much more probable for the three nodes corresponding to the interference peaks. For two 'slits' (e.g., two source nodes) separated by a distance of 2 (in lattice units), these peaks occurs at positions 0 and ±Nt, as shown in the figure.

To obtain a result that is even closer to QM, it would take an unpractical computing time. Therefore, accelerated codes will be presented in future posts that will be based on the separate convergence of the lattice boson process and of the particle boson process, a thing to be studied next.