0

I've been recently working on a JS program that models this scholarly article. The problem I've been recently having is this:

"Surface tension is usually caused by attraction between particles. Inside the fluid, this attraction, cancels out, but for the molecules near the surface, asymmetry in neighbor distribution results in a non-zero net force towards the fluid."

I'm assuming that this happens naturally in the double density relaxation method, but my program seems to show no sign of it.

Here's an excerpt from the double density relaxation method:

function doubleDensityRelaxation() {
  for(var i = 0; i<particles.length; i++) {
  var pi = particles[i];
  var p = 0;
  var pNear = 0;
  for(var j = i+1; j<particles.length; j++) {
    var pj = particles[j];
    var rij = sub(pi.r,pj.r);
    var q = length(rij)/h;
    if(q < 1) {
      p += (1-q)*(1-q);
      pNear += (1-q)*(1-q)*(1-q);
    }
  }
  var P = k*(p-p0);
  var PNear = kNear*pNear;
  var dx = [0,0];
  for(var j = i+1; j<particles.length; j++) {
    var pj = particles[j];
    var rij = sub(pi.r,pj.r);
    var q = length(rij)/h;
    if(q < 1) {

      var D = scale(normalize(rij), dtsq*(P*(1-q)+PNear*(1-q)*(1-q)));

      pj.r = add(pj.r,scale(D,0.5));
      dx = sub(dx,scale(D,0.5));
    }
  }
  pi.r = add(pi.r,dx);

  }
}
WestLangley
  • 102,557
  • 10
  • 276
  • 276

1 Answers1

0

Here's one that behaves more like a recognizable low-viscosity fluid. This sort of simulation seems to go numerically unstable at the drop of a hat, and more so the lower the frame rate.

I found what I think are a few problems with your code.

  • As I read the paper, Lij are elements of a triangular matrix of spring lengths retained across time steps; gamma and alpha control the rate at which it changes. This matrix wasn't being retained in your original code, so viscoelasticity just didn't work.

  • Inter-particle collisions are unnecessary and actually destabilize the simulation further. Beta in the viscosity pass provides that repulsion.

  • The ratio of p0/k can't be too low or you'll get instability in the form of rapidly oscillating collapsed clusters. When the clusters do break, you get individual particles bouncing off the walls like out-of-control Flubber balls.

  • The walls still have some problems: even with the tweaks I've made, particles tend to get lined up along them in columns, with new particles forced in on the bottom as old ones fountain off the top, creating a sort of twin vortex in the fluid.

  • I was experimenting with setting a lower bound on spring length, which can sometimes stabilize an otherwise unstable parameter set. My off-the-cuff guess as to the physical property it models is dilatant behavior.

Jeffrey Hantin
  • 35,734
  • 7
  • 75
  • 94