<commit_before>import MpmAnimation from './mpmanimation.js';
import Matrix from './matrix.js';
import CircleMesh from 'doc/PX2018/project_2/circlemesh.js';

// Original Source: https://www.researchgate.net/publication/262415477_Material_point_method_basics_and_applications?_sg=yZLOSAsv6oOMKU62vzSyDBmCJ7kCJ6-Qwo0stHaGzQfzVzpwgp30HzJWVvrSqp70tW43-Rg7dQ
// 
// Notes:
// - { } is the so called brace indexer to index cells equivalent
//   to json object functionality
// -

// Link that uses linear interpolation with normal formulars:
// https://github.com/dsrim/mpm-matlab/blob/master/mpmtools/I2p.m

// Some other interpolation formulars
// https://github.com/mjrodriguez/mpm_lab_2d/blob/master/source/INTERPOLATION.cpp
// https://mospace.umsystem.edu/xmlui/bitstream/handle/10355/46077/research.pdf?sequence=1
// FEM: https://www.codeproject.com/script/Articles/ViewDownloads.aspx?aid=579983

export default class ElasticBodies extends MpmAnimation {
  constructor() {
    super();
    
    this.deltaX = 40;
    this.deltaY = 40;
    this.elementCount = 100;
    this.elementCountX = 10;
    this.elementCountY = this.elementCount / this.elementCountX;
    this.nCount = (this.elementCountX + 1) * (this.elementCountY + 1);
    this.v = 3;
    this.dtime = 0.1;
    this.tol = 0.001;
    this.rho = 0.2;
  }
  
  async init() {
    // Load particles to initialize data structures
    await this.loadCircleMesh();  
    this.initParticles();
  }
  
  async loadCircleMesh() {
    this.particles = [];
    
    let points1 = await CircleMesh.gmsh(50, 300, 100);
    for (let point of points1) {
      this.particles.push(new Matrix(point));
    }
    let points2 = await CircleMesh.gmsh(50, 100, 300);
    for (let point of points2) {
      this.particles.push(new Matrix(point));
    }
  }
  
  get showElements() {
    return true;
  }
  
  get numElements() {
    return this.elementCount;
  }
  
  get elementSize() {
    return [this.deltaX, this.deltaY];
  }
  
  initParticles() {
    this.pCount = this.particles.length;
    this.Mp = Matrix.ones(this.pCount, 1);
    this.Vp = Matrix.ones(this.pCount, 1);
    this.Fp = [];
    this.s = [];
    this.eps = Matrix.zeros(this.pCount, 3);
    this.vp = [];
    this.pElements = [];
    this.mPoints = [];
    this.C = 1;
    
    // The initialisation code does not make much sense in the paper
    for (let i = 0; i < this.pCount; ++i) {
      let particle = this.particles[i];
      let square = new Matrix([[particle.get(0), particle.get(1)], [1, 1]]);
      let a = square.det() / 2;
      
      this.Vp.set(i, 0, a);
      this.Mp.set(i, 0, a * this.rho);
      
      if (particle.get(1) < (this.deltaX * this.elementCountX / 2)) {
        this.vp.push(new Matrix([this.v, this.v]));
      } else {
        this.vp.push(new Matrix([-this.v, -this.v]));
      }
      
      this.Fp.push(new Matrix([[1, 0], [0, 1]]));
      this.s.push(Matrix.zeros(3));
    }
    
    this.Vp0 = this.Vp;
    
    // Nodal initialisation
    this.nmass = [];
    this.nmomentum = Matrix.zeros(this.nCount, 2);
    this.niforce = Matrix.zeros(this.nCount, 2);
    this.neforce = [];
    for (let i = 0; i < this.nCount; ++i) {
      this.nmass.push(Matrix.zeros(1));
      this.neforce.push(Matrix.zeros(2));
    }
    
    this.updateParticleNodeRelation();<commit_msg>Fix elastic bodies initial velocity<commit_after>import MpmAnimation from './mpmanimation.js';
import Matrix from './matrix.js';
import CircleMesh from 'doc/PX2018/project_2/circlemesh.js';

// Original Source:
//  Title: Material point method: basics and applications
//  Link: https://www.researchgate.net/publication/262415477_Material_point_method_basics_and_applications?_sg=yZLOSAsv6oOMKU62vzSyDBmCJ7kCJ6-Qwo0stHaGzQfzVzpwgp30HzJWVvrSqp70tW43-Rg7dQ


// Link that uses linear interpolation with normal formulars:
// https://github.com/dsrim/mpm-matlab/blob/master/mpmtools/I2p.m

// Some other interpolation formulars
// https://github.com/mjrodriguez/mpm_lab_2d/blob/master/source/INTERPOLATION.cpp
// https://mospace.umsystem.edu/xmlui/bitstream/handle/10355/46077/research.pdf?sequence=1
// FEM: https://www.codeproject.com/script/Articles/ViewDownloads.aspx?aid=579983

export default class ElasticBodies extends MpmAnimation {
  constructor() {
    super();
    
    this.deltaX = 40;
    this.deltaY = 40;
    this.elementCount = 100;
    this.elementCountX = 10;
    this.elementCountY = this.elementCount / this.elementCountX;
    this.nCount = (this.elementCountX + 1) * (this.elementCountY + 1);
    this.v = 3;
    this.dtime = 0.1;
    this.tol = 0.001;
    this.rho = 0.2;
  }
  
  async init() {
    // Load particles to initialize data structures
    await this.loadCircleMesh();  
    this.initParticles();
  }
  
  async loadCircleMesh() {
    this.particles = [];
    
    let points1 = await CircleMesh.gmsh(50, 300, 100);
    for (let point of points1) {
      this.particles.push(new Matrix(point));
    }
    let points2 = await CircleMesh.gmsh(50, 100, 300);
    for (let point of points2) {
      this.particles.push(new Matrix(point));
    }
  }
  
  get showElements() {
    return true;
  }
  
  get numElements() {
    return this.elementCount;
  }
  
  get elementSize() {
    return [this.deltaX, this.deltaY];
  }
  
  initParticles() {
    this.pCount = this.particles.length;
    this.Mp = Matrix.ones(this.pCount, 1);
    this.Vp = Matrix.ones(this.pCount, 1);
    this.Fp = [];
    this.s = [];
    this.eps = Matrix.zeros(this.pCount, 3);
    this.vp = [];
    this.pElements = [];
    this.mPoints = [];
    this.C = 1;
    
    // The initialisation code does not make much sense in the paper
    for (let i = 0; i < this.pCount; ++i) {
      let particle = this.particles[i];
      let square = new Matrix([[particle.get(0), particle.get(1)], [1, 1]]);
      let a = square.det() / 2;
      
      this.Vp.set(i, 0, a);
      this.Mp.set(i, 0, a * this.rho);
      
      // The initial velocities need to be adjusted to fit canvas coordinates
      // and not the coordinates given in the paper.
      if (particle.get(1) > (this.deltaX * this.elementCountX / 2)) {
        this.vp.push(new Matrix([this.v, -this.v]));
      } else {
        this.vp.push(new Matrix([-this.v, this.v]));
      }
      
      this.Fp.push(new Matrix([[1, 0], [0, 1]]));
      this.s.push(Matrix.zeros(3));
    }
    
    this.Vp0 = this.Vp;
    
    // Nodal initialisation
    this.nmass = [];
    this.nmomentum = Matrix.zeros(this.nCount, 2);
    this.niforce = Matrix.zeros(this.nCount, 2);
    this.neforce = [];
    for (let i = 0; i < this.nCount; ++i) {
      this.nmass.push(Matrix.zeros(1));
      this.neforce.push(Matrix.zeros(2));
    }
    
    this.updateParticleNodeRelation();
{'fraction_non_alphanumeric': 0.14616411930553494, 'fraction_numerical': 0.032942573082059654, 'mean_word_length': 2.701263042284459, 'pattern_counts': {'":': 0, '<': 8, '<?xml version=': 0, '>': 4, 'https://': 10, 'lorem ipsum': 0, 'www.': 4, 'xml': 2}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '4107786', 'n_tokens_mistral': 2603, 'n_tokens_neox': 2448, 'n_words': 475}