6
\$\begingroup\$

Inspired by this question, I modified the script to run a simple simulation of Taylor-Green vortex using LBM and post-processes it using P5*js (an official Javascript port of the Processing API). Unfortunately, it is not performing well, but as I am new to Javascript programming, I lack the insight to really optimize the code. Getting some advice and suggestions for improvement would be greatly appreciated.

A working codepen version can be found here.

// 2D vector class
function vec2(x,y){
    this.x = x;
    this.y = y;
}
vec2.prototype = {
    scale: function(s){
        return new vec2(this.x*s,this.y*s);
    },
    add: function(v){
        return new vec2(this.x+v.x,this.y+v.y);
    },
    subtract: function(v){
        return new vec2(this.x-v.x,this.y-v.y);
    },
    dot: function(v){
        return this.x*v.x+this.y+v.y;
    },
    length: function(){
        return Math.sqrt(this.x*this.x+this.y*this.y);
    },
    normalise: function(){
        var l = 1/this.length();
        return new vec2(this.x*l, this.y*l);
    },
};

// domain class
function domain(nx, ny) {
    this.nx = nx; // domain width
    this.ny = ny; // domain height
    this.f = []; // distribution array
    this.ftmp = []; // temporary array
    this.dens = []; // density array
    this.vel = []; // velocity array
    this.omega = 1; // relaxation frequency
    this.e = [ // discrete velocity set
        new vec2(0,0),
        new vec2(0,1), new vec2(1,0), new vec2(0,-1), new vec2(-1,0),
        new vec2(1,1), new vec2(1,-1), new vec2(-1,-1), new vec2(-1,1)
    ];
    this.w  = [ // weights
        4/9, 
        1/9, 1/9, 1/9, 1/9, 
        1/36, 1/36, 1/36, 1/36
    ];
    // Arrays initialization
    for (var x=0; x<this.nx; x++) {
        this.f[x] = [];
        this.ftmp[x] = [];
        this.dens[x] = [];
        this.vel[x] = [];
        for (var y=0; y<this.ny; y++) {
            this.f[x][y] = [];
            this.ftmp[x][y] = [];
        }
    }
}
domain.prototype = {
    init: function(){
        // Initializes Taylor-Green vortex
        var kx = 2*Math.PI/this.nx;
        var ky = 2*Math.PI/this.ny;
        var kxkx = kx*kx;
        var kyky = ky*ky;
        var ksq = kxkx + kyky;
        var k = Math.sqrt(ksq);
        var dens0 = 1;
        this.umax = 0.1;
        var u0 = 4*this.umax;
        this.densmax = dens0 + 3*dens0*u0*u0/4; 
        for (var x=0; x<this.nx; x++){
            for (var y=0; y<this.ny; y++){
                var u = u0*ky/k*Math.cos(kx*x)*Math.sin(ky*y);
                var v = -u0*kx/k*Math.sin(kx*x)*Math.cos(ky*y);
                this.dens[x][y] = dens0 + 3*dens0*u0*u0/4*(kyky/ksq*Math.cos(2*kx*x)+kxkx/ksq*Math.sin(2*ky*y));
                this.vel[x][y] = new vec2(u,v);
                for(var i=0; i<9; i++){
                    // Initialize using equilibrium distribution
                    var uu = this.vel[x][y].x*this.vel[x][y].x + this.vel[x][y].y*this.vel[x][y].y;
                    var eu = this.e[i].x*this.vel[x][y].x + this.e[i].y*this.vel[x][y].y;
                    this.f[x][y][i] = this.w[i]*this.dens[x][y]*(1+3*eu+4.5*eu*eu-1.5*uu);
                }
            }
        }
    },
    collide: function(){
        for(var x=1; x<this.nx-1; x++){
            for(var y=1; y<this.ny-1; y++){
                // calculate density
                var rho = 0;
                for(var i=0; i<9; i++){
                    rho += this.f[x][y][i];
                }
                this.dens[x][y] = rho;

                // calculate velocity
                var u = new vec2(0,0);
                for(var i=1; i<9; i++){
                    u  = u.add( this.e[i].scale( this.f[x][y][i] ) );
                }
                u = u.scale( 1/rho );
                this.vel[x][y] = u;

                // Perform collision step and save to temp array ftmp
                var uu = u.x*u.x + u.y*u.y;
                for(var i=0; i<9; i++){
                    var eu = u.x*this.e[i].x + u.y*this.e[i].y;
                    var fiEq = this.w[i]*rho*(1+3*eu+4.5*eu*eu-1.5*uu);
                    var fiCol = -this.omega*(this.f[x][y][i]-fiEq); // bgk
                    this.ftmp[x][y][i] = this.f[x][y][i] + fiCol;
                }
            }
        }
    },
    periodic: function(){
        // Apply periodic boundary conditions on ftmp
        // x-periodic
        for(var y=1; y<this.ny-1; y++){ 
            for(var i=0; i<9; i++){
                this.ftmp[0][y][i] = this.ftmp[this.nx-2][y][i];
                this.ftmp[this.nx-1][y][i] = this.ftmp[1][y][i];
            }
        }
         // y-periodic
        for(var x=1; x<this.nx-1; x++){ 
            for(var i=0; i<9; i++){
                this.ftmp[x][0][i] = this.ftmp[x][this.ny-2][i];
                this.ftmp[x][this.ny-1][i] = this.ftmp[x][1][i];
            }                        
        }
        // corner treatment
        for(var i=0; i<9; i++){
            this.ftmp[0][0][i] = this.ftmp[this.nx-2][this.ny-2][i];
            this.ftmp[this.nx-1][this.ny-1][i] = this.ftmp[1][1][i];
            this.ftmp[this.nx-1][0][i] = this.ftmp[1][this.ny-2][i];
            this.ftmp[0][this.ny-1][i] = this.ftmp[this.nx-2][1][i];
        }
    },
    stream: function(){
        // Perform streaming step ftmp -> f
        for(var x=1; x<this.nx-1; x++){
            for(var y=1; y<this.ny-1; y++){
                for(var i=0; i<9; i++){
                    this.f[x][y][i] = this.ftmp[x-this.e[i].x][y-this.e[i].y][i];
                }
            }
        }
    }
}
function simulation(){
    var sim = function(p) {
        var nx = 200, ny = 200;
        var myDomain = new domain(nx, ny);
        p.setup = function() {
            p.createCanvas(nx, ny)
             .parent('sim');
            //p.frameRate(30);
            myDomain.init();
            p.noStroke();
            p.colorMode(p.RGB, 1);
        }
        p.draw = function() {
            myDomain.collide();
            myDomain.periodic();
            myDomain.stream();
            var dens = myDomain.dens;
            var vel = myDomain.vel;
            var densmax = myDomain.densmax;
            var umax = myDomain.umax;
            for (var x = 0; x < p.width; x++) {
                for (var y = 0; y < p.height; y++ ) {
                    //var v = dens[x][y]/densmax;
                    var v = vel[x][y].length()/umax;
                    p.stroke(v, 0, 1-v);
                    p.point(x, y);
                }
            }
        }
    }
    return new p5(sim);
}
$( document ).ready(simulation());

Edit 1: Running the profiler in Chrome results in:

enter image description here

This shows that most of the CPU is being used by P5*JS rather than the code.

\$\endgroup\$
1
  • 2
    \$\begingroup\$ Hi! Welcome to Code Review. We only review working code; I know your code is working, but we cannot really answer "Anyone any idea why dot method goes wrong?" \$\endgroup\$ Commented Oct 9, 2015 at 1:57

1 Answer 1

2
\$\begingroup\$

I found that this specific part of the processing code was performing badly:

for (var x = 0; x < p.width; x++) {
  for (var y = 0; y < p.height; y++ ) {
    // removed irrelevant code
    p.stroke(v, 0, 1-v);
    p.point(x, y);
  }
}

It turns out that this write to the screen at every point. It is much better to buffer the writes by set-ting the pixels to certain color and then writing the whole buffer using updatePixels() like:

for (var x = 0; x < p.width; x++) {
  for (var y = 0; y < p.height; y++ ) {
    // removed irrelevant code
    p.color(v, 0, 1-v);
    p.set(x, y, c);
  }
}
p.updatePixels();

Now the draw calls take up about as much time as the collide calls (about 35%). I'm still looking to reduce this further.

\$\endgroup\$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.