Rutherford散乱

int NUM=24;
int dia=10;
float dt=1;
float[][] pos=new float[NUM][2];
float[][] tmp_pos=new float[NUM][2];
float[][] velocity=new float[NUM][2];
float[][] dotvelocity=new float[NUM][2];

void setup(){
  size(400,400);
  background(255);
  for(int i=0;i<NUM;i++){
    pos[i][0]=-width/2;
    pos[i][1]=-height/2+dia/2+i*17;
    tmp_pos[i][0]=pos[i][0];
    tmp_pos[i][1]=pos[i][1];
    velocity[i][0]=1;
    velocity[i][1]=0; 
    dotvelocity[i][0]=0;
    dotvelocity[i][1]=0;
  }
  fill(255,150,0);
  ellipse(width/2,height/2,50,50);
  fill(0);
  textSize(30);
  textAlign(CENTER);
  text("Au",width/2+60,height/2+10);
}

void draw(){
  Calculate();
  Disp();
}

void Calculate(){
   for(int i=0;i<NUM;i++){
    float r=mag(tmp_pos[i][0],tmp_pos[i][1]);
    dotvelocity[i][0]=tmp_pos[i][0]*40*dt/pow(r,3);
    dotvelocity[i][1]=tmp_pos[i][1]*40*dt/pow(r,3);
    velocity[i][0]+=dotvelocity[i][0];
    velocity[i][1]+=dotvelocity[i][1];
    pos[i][0]+=velocity[i][0]*dt;
    pos[i][1]+=velocity[i][1]*dt;
    tmp_pos[i][0]=pos[i][0];
    tmp_pos[i][1]=pos[i][1];
  }
}

void Disp(){
  for(int i=0;i<NUM;i++){
    strokeWeight(0.01);
    fill(255,0,0,100);
    ellipse(width/2+pos[i][0],height/2+pos[i][1],dia,dia);
  }
}