Angular Velocity Matrix demo

Er, JavaScript doesn't seem to be running, so basically all you're going to see is the bare code...

// Make a bunch of spheres
var g=new THREE.SphereGeometry(1, 3, 3);
var m=new THREE.MeshLambertMaterial(
    {color: 0xffffff})
sim.nspheres=5;
sim.spheres=[];
for (var i=0;i<sim.nspheres;i++) {
  var s=new THREE.Mesh(g,m);
  s.position.set(0,0,0);
  s.V=new vec3(+20.0*Math.random(),Math.random(),+15.0*Math.random());
  s.P=new vec3(-3.0,0,0);
  s.castShadow=true;
  s.receiveShadow=true;
  scene.add(s); // so it gets rendered
  sim.spheres[i]=s; // so we can animate it
}

// Ground
var groundTex=THREE.ImageUtils.loadTexture( "textures/checkerboard_noisy.jpg" );
groundTex.wrapS=groundTex.wrapT=THREE.RepeatWrapping;
groundTex.repeat=new vec2(25.0,25.0);

sim.ground = new THREE.Mesh(
  new THREE.CubeGeometry(50,50,0.0001),
  new THREE.MeshLambertMaterial(
    {color: 0xffccaa, opacity: 1, map:groundTex})
);
sim.ground.position.z=-1;  /* hack for radius 1 */
sim.ground.receiveShadow=true;
scene.add(sim.ground);

// Sun-like spotlight
var l=new THREE.SpotLight();
sim.light = l;
l.position.set( -50, -100, 100 );
l.castShadow=true;
l.shadowCameraNear = 50; 
l.shadowCameraFar = 500; 
l.shadowCameraFov = 10; // degrees field of view
scene.add(l);

camera.lookAt(scene.position);

/** Matrix preliminaries **/
// Assign this column of the matrix from this vector.
  THREE.Matrix4.prototype.setCol=function(col,v) {
     this.elements[col*4+0]=v.x;
     this.elements[col*4+1]=v.y;
     this.elements[col*4+2]=v.z;
  }

// Extract this column of the matrix, and return it.
  THREE.Matrix4.prototype.getCol=function(col) {
     return new vec3(this.elements[col*4+0], this.elements[col*4+1], this.elements[col*4+2]);
  }

// Normalize the vec3 making up this column.
  THREE.Matrix4.prototype.normCol=function(col) {
     var v=this.getCol(col);
     v.normalize();
     this.setCol(col,v);
  }

// Normalize the XYZ axes, so they have unit length.
  THREE.Matrix4.prototype.normalize=function() {
     this.normCol(0);
     this.normCol(1);
     this.normCol(2);
  }

// Do time integration: add v*dt to each element.
  THREE.Matrix4.prototype.integrate=function(v,dt) {
      var te=this.elements; var ve=v.elements;
      // loop ends early, so we don't mess up W.
      for (var i=0;i<15;i++) 
          te[i]+=ve[i]*dt;
  };

// User Interface and camera control
var rot=new vec3(); // user-applied torques

// X control via J and L
if (lib.key['j']) rot.pe(new vec3(-1,0,0));
if (lib.key['l']) rot.pe(new vec3(+1,0,0));
// Y control via I and K
if (lib.key['i']) rot.pe(new vec3(0,+1,0));
if (lib.key['k']) rot.pe(new vec3(0,-1,0));
// Z control via U and M
if (lib.key['u']) rot.pe(new vec3(0,0,+1));
if (lib.key['m']) rot.pe(new vec3(0,0,-1));

// Per-object rotation physics
for (var i=0;i<sim.nspheres;i++) {
  var s=sim.spheres[i];
  s.matrixAutoUpdate=false; // we will write matrix
  s.matrixWorldNeedsUpdate=true; // recompute world
  var m=s.matrix;

// Angular velocity, radians/sec per axis
//   This is in body-local coordinates.
  if (!s.W) { s.W=new vec3(0,1.0,0); }
  s.W.pe(rot.t(lib.dt*1.0)); // keyboard torque in

// Compute skew matrix, which is dm/dt (change in m)
  var skew=new THREE.Matrix4().identity();
  skew.elements[1]-=s.W.z;
  skew.elements[2]+=s.W.y;
  
  skew.elements[4]+=s.W.z;
  skew.elements[6]-=s.W.x;
 
  skew.elements[8]-=s.W.y;
  skew.elements[9]+=s.W.x;
  
// Rotate skew matrix by our current rotation
// (puts it in body coordinates--else rotation stops!)
  skew.multiply(m);

// Time integrate skew matrix (m+=skew*dt)
  m.integrate(skew,lib.dt);

// Normalize our matrix (or else object keeps growing!)
  m.normalize();

// Drop in the position column
  m.setCol(3,s.P);
}

/*********** Collisions and Linear Motion ********/
/**
  s: sphere object, with s.P position and s.V velocity.
  N: surface normal, facing out of region.
  P: point on surface of boundary region.
*/
function boundaryBounce(s,N,P) {
  N.normalize();
  var rP=s.P.m(P); // sphere position, relative to boundary
  var d=rP.dot(N);
  if (d>0.0) { // we hit the boundary
    s.P.pe(N.t(-d)); // move inside
    var dV=s.V.dot(N); 
    if (dV>0.0) { // flip velocity
      s.V.pe(N.t(-2.0*dV ));
    }
    s.V.te(0.6); // bounce friction (inelastic)

}
}

// Linear bounce physics
for (var i=0;i<sim.nspheres;i++) {
  var s=sim.spheres[i];
  var A=new vec3(0,0,-9.8); 
  s.V.pe(A.t(lib.dt));
  s.P.pe(s.V.t(lib.dt));

boundaryBounce(s,new vec3(0,0,-1), new vec3()); // floor
  boundaryBounce(s,new vec3(+1,0,0), new vec3(5,0,0));
  boundaryBounce(s,new vec3(-1,0,0), new vec3(-4,0,0));
  boundaryBounce(s,new vec3(0,0,+1), new vec3(0,0,5)); // lid

s.position=s.P; // copy to 3D display
}