import {EVENT_OBJECT_CREATE, EVENT_OBJECT_MERGE, OBJECT_HISTORY_SIZE} from './config.js'; import {MassObject} from './object.js'; import { add, copy, cross, degrees, direction, div, dot, magnitude, mult, square, sub, weightedAvg, zero } from './vector.js'; export class System { objects = []; creatingObject = undefined; selectedObject = undefined; selectObjectStart = undefined; paused = false; constructor(sim) { this.sim = sim; } toJSON() { return { objects: this.objects.map(obj => obj.toJSON()), } } fromJSON({objects} = {}) { objects = objects ?? []; // Replace current state with the provided one. // Assumes a backup has already been saved if desired. this.objects = []; for (const objectJSON of objects) { const obj = new MassObject(this.sim, 0, 0); obj.fromJSON(objectJSON); this.objects.push(obj); } } handlePointerDown({x, y}) { // If pointer is touching an object, select the object const touchingObject = this.objectAtLocation(x, y); if (touchingObject !== undefined) { this.selectObject(touchingObject, {x, y}); } else { // Otherwise, create a new object this.createObject(x, y); } } handlePointerUp() { const obj = this.getSelectedOrCreating(); if (obj === undefined) return; this.doneCreatingObject(); this.deselect(); // Convert pointer velocity to simulation scale obj.velocity = div(this.sim.pointer.latestVelocity, this.sim.display.scale); // Including time scale - if time is slow, our motion is relatively faster if (this.sim.getOption('compensate.timeScale')) { obj.velocity = div(obj.velocity, this.sim.timeScale); } obj.velocity = add(obj.velocity, this.sim.panning.velocity); } handlePointerMove(r) { // If the cursor moves while creating an object, or while an object is selected, // update the position using the pointer motion but the velocity using the pointer velocity const obj = this.getSelectedOrCreating(); if (obj === undefined) return; const start = this.selectedObjectStart; const delta = sub(r, start.pointer); // TODO: Calculate work done by pointer here? // Either interpolate the acceleration and use m*a, or // measure the change in the object's kinetic energy // Or both! // In the frame loop, for the selected object we already calculate // the work done by the pointer in preventing the acceleration the object // would have experienced. // But that doesn't capture the work done by the pointer moving the object. // If the sim is paused, the dot(force, displacement) method loses accuracy. // If the sim is paused though, then channge in kinetic energy will be easy to measure. // If the sim is not paused, we can use dot(force, displacement) method. obj.position = add(start, delta); obj.velocity = zero; } frame(elapsedTime) { // If we're creating an object, increment its mass // with the mass creation rate accelerating over time if (this.creatingObject !== undefined) { const obj = this.objects[this.creatingObject]; let massCreationRate = this.sim.getOption('param.massCreationRate'); massCreationRate /= this.sim.display.scale; // Keep consistent time scale if (this.sim.getOption('compensate.timeScale')) { massCreationRate /= this.sim.timeScale; } obj.mass += massCreationRate * elapsedTime; } // Calculate forces due to gravity. this.computeForces(); if (this.sim.playing) { // Predict positions (Velocity verlet method) this.forEachObject(obj => { obj.currentAcceleration = copy(obj.acceleration); // If this object is being created/selected, we're not going to let it move... // but we can calculate the work being done by holding it in place. obj.currentPosition = copy(obj.position); obj.position = add(obj.position, mult( elapsedTime, add( obj.velocity, mult(obj.currentAcceleration, elapsedTime / 2) ), )); }); // Recompute forces this.computeForces(); // Predict velocities this.forEachObject(obj => { const acceleration = copy(obj.acceleration); obj.acceleration = div(add(obj.currentAcceleration, acceleration), 2); obj.velocity = add(obj.velocity, mult(obj.acceleration, elapsedTime)); // If the user is positioning this object, we'll leave its position unchanged; // but let's compute how much work we're doing to accomplish it! if (obj.id === this.getSelectedOrCreating()?.id) { const delta = sub(obj.currentPosition, obj.position); const netForce = mult(obj.acceleration, obj.mass); const work = dot(netForce, delta); obj.position = obj.currentPosition; obj.workDoneByPointer += work; } // Append to object history obj.history.push({position: {...obj.position}}); // TODO: store object color changes in history // Enforce object history length while (obj.history.length > OBJECT_HISTORY_SIZE) { obj.history.shift(); } }); // Collisions this.forEachObject((A, i) => { this.forEachObject((B, j) => { if (this.objectsOverlap(A, B)) { this.bounceObjects(i, j, elapsedTime); } }, {alive: true, startWith: i + 1}); }); } // Display total energy const totalKE = this.reduce((total, obj) => { const energy = obj.mass * square(obj.velocity) / 2; return total + energy; }, 0, {alive: true}); const G = this.sim.getOption('param.gravity'); const totalGPE = this.reduce((total, A, i) => { return total + this.reduce((objTotal, B) => { const energy = -G * A.mass * B.mass / magnitude(sub(A.position, B.position)); return objTotal + energy; }, 0, {alive: true, startWith: i + 1}); }, 0, {alive: true}); const netMomentum = this.reduce((net, obj) => add(net, obj.momentum), zero); this.sim.info['Total E'] = (totalKE + totalGPE).toExponential(2); this.sim.info['Total K'] = totalKE.toExponential(2); this.sim.info['Total Ug'] = totalGPE.toExponential(2); this.sim.info['Net Momentum'] = magnitude(netMomentum).toExponential(2); // Display objects info // First clear info from previous frame this.forEachObject((_obj, i) => { delete this.sim.info[`Object ${i}`]; }, {alive: null}); if (this.sim.getOption('debug.objectsInfo')) { const aliveOnly = this.sim.getOption('debug.aliveObjects'); this.forEachObject((obj, i) => { const speed = magnitude(obj.velocity); const accel = magnitude(obj.acceleration); // Give angle counterclockwise from horizontal const velocityDir = -1 * degrees(direction(obj.velocity)); const accelDir = -1 * degrees(direction(obj.acceleration)); const {r, g, b} = obj.color; this.sim.info[`Object ${i}`] = [ `  `, `${obj.position.x.toPrecision(4)}, `, `${obj.position.y.toPrecision(4)}, `, `${obj.mass.toPrecision(4)} kg, `, `${speed.toPrecision(2)} m/s, ${velocityDir.toPrecision(2)}°`, `${accel.toPrecision(2)} m/s2, ${accelDir.toPrecision(2)}°`, `Alive: ${obj.alive}`, ]; }, {alive: aliveOnly || null}); } // Render the objects this.drawObjects(); } // Pause and resume to enable automatic pause on object create/select // in this mode (mass generation) pause() { this.sim.pause(); this.paused = true; } resume() { if (this.paused) { this.paused = false; this.sim.play(); } } // Create an object with mass that grows as pointer is held down createObject(x, y) { const idx = this.objects.length; const obj = new MassObject(this.sim, x, y); this.creatingObject = idx; this.selectedObjectStart = {x, y, pointer: {x, y}}; this.objects.push(obj); // Pause the simulation during mass creation; this avoids some complex local dynamics if (this.sim.getOption('pauseDuring.creation')) { this.pause(); } obj.velocity = copy(this.sim.panning.velocity); const e = new CustomEvent(EVENT_OBJECT_CREATE, {detail: {obj}}); this.sim.div.dispatchEvent(e); } doneCreatingObject() { if (this.creatingObject !== undefined) { this.creatingObject = undefined; this.resume(); } } // cb: (obj) => undefined onCreate(cb) { this.sim.div.addEventListener(EVENT_OBJECT_CREATE, ({detail: {obj}}) => { cb(obj); }); } // cb: ({surviving, merged}) => undefined onMerge(cb) { this.sim.div.addEventListener(EVENT_OBJECT_MERGE, ({detail: {surviving, merged}}) => { cb({surviving, merged}); }); } object(i) { return this.objects[i]; } selectObject(i, pointer) { this.selectedObject = i; const {x, y} = this.object(i).position; this.selectedObjectStart = {x, y, pointer}; if (this.sim.getOption('pauseDuring.selection')) { this.pause(); } } deselect() { this.selectedObject = undefined; this.selectedObjectStart = undefined; this.resume(); } getSelectedOrCreating() { let i = this.creatingObject ?? this.selectedObject; if (i !== undefined) { return this.objects[i]; } } get length() { return this.objects.length; } getBoundingBox(objects = []) { const box = this.reduce(({start, end}, obj) => { if (objects.length && !objects.includes(obj)) return {start, end}; const lx = obj.position.x - obj.radius; const gx = obj.position.x + obj.radius; const ly = obj.position.y - obj.radius; const gy = obj.position.y + obj.radius; if (start.x === undefined) { return { start: {x: lx, y: ly}, end: {x: gx, y: gy}, }; } return { start: { x: Math.min(start.x, lx), y: Math.min(start.y, ly), }, end: { x: Math.max(end.x, gx), y: Math.max(end.y, gy), } }; }, { start: {x: undefined, y: undefined}, end: {x: undefined, y: undefined}, }); box.start.x = (box.start.x ?? 0); box.start.y = (box.start.y ?? 0); box.end.x = (box.end.x ?? 0); box.end.y = (box.end.y ?? 0); return box; } objectAtLocation(x, y) { let idx = undefined; this.selectedObjectStart = undefined; this.forEachObject((obj, i) => { // If distance to object is less than object's radius, we are touching the object const dist = magnitude(sub(obj.position, {x, y})); if (dist <= obj.radius) { idx = i; return null; } }); return idx; } // cb: (obj, idx) => {} forEachObject(cb, {alive, startWith} = {}) { if (alive === undefined) alive = true; for (let i = startWith ?? 0; i < this.objects.length; i++) { const obj = this.objects[i]; if (alive === null || alive == obj.alive) { const ret = cb(obj, i); if (ret === null) break; } } } drawObjects() { // Draw all paths, all objects, and then all arrows this.forEachObject(obj => obj.drawPath(this.sim.display), {alive: null}); this.forEachObject(obj => obj.drawObject(this.sim.display), {alive: true}); this.forEachObject(obj => obj.drawSelection(this.sim.display), {alive: true}); this.forEachObject(obj => obj.drawArrows(this.sim.display), {alive: true}); } // cb: (acc, obj, idx) => {} reduce(cb, initial, opts) { let acc = initial; this.forEachObject((obj, idx) => { const ret = cb(acc, obj, idx); if (ret !== undefined) { acc = ret; } }, opts); return acc; } // cb: (obj, idx) => boolean filter(cb, opts) { let objects = []; this.forEachObject((obj, idx) => { const ret = cb(obj, idx); if (ret) { objects.push(obj); } }, opts); return objects; } objectsOverlap(A, B) { const dx = (B.position.x - A.position.x); const dy = (B.position.y - A.position.y); const dSquared = dx ** 2 + dy ** 2; const d = Math.sqrt(dSquared); return d < A.radius + B.radius; } bounceObjects(i, j, elapsedTime) { const A = this.object(i); const B = this.object(j); const totalMass = A.mass + B.mass; // TODO: Handle scenario where an object overlaps more than one other object const autoMerge = true; const Z = this.sim.getOption('param.elasticity'); // Elasticity: 0 = inelastic, 1 = elastic const r = sub(A.position, B.position); const normal = div(r, magnitude(r)); const tangent = {x: -normal.y, y: normal.x}; const vAt = dot(A.velocity, tangent); // Require that objects are moving toward each other! // Normal is directed toward A, so they are moving toward each other if vBn > vAn const vAn = dot(A.velocity, normal); const vBn = dot(B.velocity, normal); if (vAn >= vBn) { // The objects are already moving away from each other return; } // TODO: Relative tangential velocity, and initial rotation of each object, should impart angular momentum // Also TODO: In partially elastic case, compute partial angular momentum transfer // Are these objects sticking together? // One way to determine this is, // does either object have enough KE to escape their gravitational potential? if (autoMerge) { // const vSquared = Math.abs(vAn - vBn) ** 2; const vSquared = square(sub(A.velocity, B.velocity)); const G = this.sim.getOption('param.gravity'); const mergeThreshold = this.sim.getOption('param.mergeThreshold'); if (vSquared < mergeThreshold * G * totalMass / magnitude(r)) { // Neither object looks like it can escape! this.mergeObjects(i, j); return; } } // Calculate the collision const Zscaled = Math.sqrt(Z); const vAnNew = (vAn * (A.mass - Zscaled * B.mass) + vBn * (1 + Zscaled) * B.mass) / totalMass; const vA = add(mult(tangent, vAt), mult(normal, vAnNew)); const vB = div(add(mult(A.velocity, A.mass), mult(B.velocity, B.mass), mult(vA, -A.mass)), B.mass); // What about position?? The objects are still overlapping. // If we fire and forget, we're hoping the object escapes within the next frame. // But it often doesn't. We catch some of these by verifing a positive relative radial velocity. // However, if the objects have sufficient tangential velocity, // or if we're imparting tangential velocity due to initial rotation in partially elastic cases, // then they can remain overlapping in the next frame despite having substantial velocity. // Then, in the next frame, since the masses are so close together, they accelerate toward each other // and may remain overlapping in the next frame; and so they can orbit in this overlapping manner. // However, this is an unstable orbit due to its sensitivity to parameters, including framerate! // Therefore it is desireable not only to mitigate, but to prevent this scenario by enforcing that // the objects never remain overlapping at the end of a frame computation! // The question then becomes, what is the best accuracy we can reasonably hope to achieve? // We are altering the objects gravitational potentials if we move them from their simulated positions. // This might be fine given the relative infrequence of collisions. // 1. So a first order solution would be to back the object to its initial position, and leave it there. // 2. An improvement could be to move it from its initial position, to the position it would have after the same duration // as the current frame had elapsed, but with the velocity we've computed for it after collision. // 3. A further improvement would be to linearly interpolate along its computed path, and then update its position // using the remaining time in the current interval. // 4. A yet further improvement would be to interpolate quadratically when estimating the point of collision! // The linear interpolation solution sounds likely to be sufficient, given that the risk of error is not even a change in // the resulting deflection angle, but only of the exact position of the object. // What are likely issues with the simplest solution? An object could perhaps get stuck in place if we fail to update its position. // That brings us to the second solution. It sounds likely fine; the risk of error is just a slight change in the object's position, // which is guaranteed to be less than a frame's worth of distance travelled. // TODO: Specify desired level of accuracy ^_^ // Well, by experimenting we found the problem with the second solution; because the objects are close together, // giving them an altitude boost greatly diminishes their gravitational potential, increasing their separation. // It turns out that's the problem with the first solution too! // But it depends on how far into the frame the objects collide. Which brings us to interpolation. // Determine time of collision, relative to start of frame. // We know the objects are accelerating toward each other, possibly rapidly. // Therefore we really need to interpolate quadratically. // But even that will not be accurate, as the force is increasing substantially during the interval as well. // // What we may really want to do is to run the whole frame computation multiple times, searching for a moment before impact, to whatever specified precision. // Then we re-run the whole frame from that point... // Another pathway here is to approximate the resulting position but then to adjust // the resulting velocity in order to preserve the object's final energy, // by subtracting for the decrease in gravitational potential energy. // The objects never should have gotten so close... we should probably do interpolation and re-run whole frame. A.position = add(A.currentPosition, mult(elapsedTime, vA)); B.position = add(B.currentPosition, mult(elapsedTime, vB)); const G = this.sim.getOption('param.gravity'); const dV = 2 * G * ( 1 / magnitude(sub(A.position, B.position)) - 1 / magnitude(sub(A.currentPosition, B.currentPosition)) ); const vAmag = Math.sqrt(square(A.velocity) - B.mass * dV); const vBmag = Math.sqrt(square(B.velocity) - A.mass * dV); A.velocity = mult(vA, div(vAmag, magnitude(vA))); B.velocity = mult(vB, div(vBmag, magnitude(vB))); } mergeObjects(i, j) { const A = this.object(i); const B = this.object(j); let S, T; // Merge the older into the newer, in order to provide mass creation rate continuity if (A.age > B.age) { // A merges into B; B survives S = B; T = A; // If A was selected or being created, select S instead if (this.creatingObject === i) this.creatingObject = j; if (this.selectedObject === i) this.selectedObject = j; } else { // B merges into A; A survives S = A; T = B; // If B was selected or being created, select S instead if (this.creatingObject === j) this.creatingObject = i; if (this.selectedObject === j) this.selectedObject = i; } // Merge T into S: // Set position = center of mass // Set velocity = total momentum / total mass // Combine forces // Add masses // Average color S.position = weightedAvg([[S.position, S.mass], [T.position, T.mass]]); S.velocity = weightedAvg([[S.velocity, S.mass], [T.velocity, T.mass]]); S.mass += T.mass; S.color = { r: (S.mass * S.color.r + T.mass * T.color.r) / (S.mass + T.mass), g: (S.mass * S.color.g + T.mass * T.color.g) / (S.mass + T.mass), b: (S.mass * S.color.b + T.mass * T.color.b) / (S.mass + T.mass), }; T.alive = false; T.forces = []; const e = new CustomEvent(EVENT_OBJECT_MERGE, {detail: {surviving: S, merged: T}}); this.sim.div.dispatchEvent(e); } computeForces() { const gravity = this.sim.getOption('param.gravity'); if (this.objects.length < 2) return; this.forEachObject(obj => { obj.forces = []; }); this.forEachObject((A, i) => { this.forEachObject(B => { const r = sub(B.position, A.position); const dSquared = square(r); const d = Math.sqrt(dSquared); // If the objects are overlapping, don't apply gravity; // This should give them a chance to fully rebound, and avoid accidental capture. // if (d <= A.radius + B.radius) { // return; // } const F = gravity * A.mass * B.mass / dSquared; const Fx = F * r.x / d; const Fy = F * r.y / d; // Equal and opposite forces A.forces.push({x: Fx, y: Fy}); B.forces.push({x: -Fx, y: -Fy}); if (A.currentPosition && B.currentPosition) { A.workDoneByForces += dot({x: Fx, y: Fy}, sub(A.position, A.currentPosition)); B.workDoneByForces += dot({x: Fx, y: Fy}, sub(B.position, B.currentPosition)); } }, {alive: true, startWith: i + 1}); }); // Additional computations this.forEachObject(obj => { // Acceleration obj.acceleration = obj.getAcceleration(); }); } computeSystemCenter(objects = []) { // Determine center of mass const {totalMass, count, totalMassLocation} = this.reduce((acc, obj) => { if (objects.length && !objects.includes(obj)) return acc; return { count: acc.count + 1, totalMass: acc.totalMass + obj.mass, totalMassLocation: add(acc.totalMassLocation, mult(obj.position, obj.mass)), }; }, { totalMassLocation: {x: 0, y: 0}, totalMass: 0, count: 0, }); const centerOfMass = count ? div(totalMassLocation, totalMass) : zero; // Determine average momentum const netMomentum = this.reduce((acc, obj) => { if (objects.length && !objects.includes(obj)) return acc; return add(acc, mult(obj.velocity, obj.mass)); }, zero, {alive: true}); return {totalMass, count, totalMassLocation, centerOfMass, netMomentum}; } computeSystemAngularMomentum(centerOfMass) { if (!centerOfMass) { const sys = this.computeSystemCenter(); centerOfMass = sys.centerOfMass; } return this.reduce((acc, obj) => { // Angular momentum for each object is m * s / d // where d is the distance of the object from the global center of mass // and s is the magnitude of the cross product of v and r const r = sub(obj.position, centerOfMass); const s = cross(obj.velocity, r); const d = Math.sqrt(r.x ** 2 + r.y ** 2); return acc + obj.mass * s / d; }, 0); } }