Skip to content
Snippets Groups Projects
Unverified Commit 3ce4ac05 authored by Jamie Temple's avatar Jamie Temple
Browse files

feat: adaptive (not-working)

parent 96784233
Branches
No related tags found
No related merge requests found
...@@ -3,10 +3,50 @@ import { PipelineObserver } from "../core/Pipeline"; ...@@ -3,10 +3,50 @@ import { PipelineObserver } from "../core/Pipeline";
import { Particle } from './physics/particle'; import { Particle } from './physics/particle';
import { SimulationData } from "./SimulationData"; import { SimulationData } from "./SimulationData";
function adaptive(x: Array<number>, h: number, eMax = 1e-10, n: number = 5): number {
const rk4 = SimulationData.INTEGRATORS[2];
// 1. calculate approximate
let xFull = rk4.step([...x], h, Simulation.derivation);
let xHalf = rk4.step([...x], h/2, Simulation.derivation);
xHalf = rk4.step(xHalf, h/2, Simulation.derivation);
// 2. calculate error
let eList = [];
for (let i = 0; i < x.length; i++) {
eList.push(Math.abs(xFull[i] - xHalf[i]));
}
let e = 0;
for (let i = 0; i < eList.length; i++) {
e += Math.pow(eList[i], 2);
}
if (e < 1e-8) {
return h;
}
e = Math.sqrt(e); // total error
// 3. e = k * h^n
// e.g. h = 0.001 mit n = 5 => e = k * 0.001^5
// => k = e / h^n
// => k * h^n = eMax
// <=> h = pow(eMax / k, 1/n)
// => h = pow(eMax / (e / h^n), 1/n)
let newH = Math.pow(eMax / (e / Math.pow(h, n)), 1 / n);
if (isNaN(newH) || !isFinite(newH)) {
newH = h;
}
return newH;
}
export class Simulation implements PipelineObserver { export class Simulation implements PipelineObserver {
public static readonly G = new THREE.Vector3(0, -9.807, 0); public static readonly G = new THREE.Vector3(0, -9.807, 0);
public static current: Particle; public static current: Particle;
public static h = .001;
constructor() { constructor() {
SimulationData.generateTextile(); SimulationData.generateTextile();
...@@ -14,20 +54,34 @@ export class Simulation implements PipelineObserver { ...@@ -14,20 +54,34 @@ export class Simulation implements PipelineObserver {
update(_deltaTime: number): void { update(_deltaTime: number): void {
const h = _deltaTime / 100;
for (let t = 0; t < _deltaTime; t += h) { if (SimulationData.paused) {
return;
}
for (let t = 0; t < _deltaTime; t += Simulation.h) {
SimulationData.pList?.forEach(pSub => { SimulationData.pList?.forEach(pSub => {
pSub.forEach(p => { pSub.forEach(p => {
if (p.usePhysics) { if (p.usePhysics) {
Simulation.current = p; Simulation.current = p;
let x = [p.position.x, p.position.y, p.position.z, let x = [p.position.x, p.position.y, p.position.z,
p.velocity.x, p.velocity.y, p.velocity.z]; p.velocity.x, p.velocity.y, p.velocity.z];
x = SimulationData.INTEGRATORS[SimulationData.integratorType].step(x, h, Simulation.derivation);
if (SimulationData.stepSizeMethod == 1) {
Simulation.h = adaptive(x, Simulation.h);
}
x = SimulationData.INTEGRATORS[SimulationData.integrationMethod]
.step(x, Simulation.h, Simulation.derivation);
p.setState(x[0], x[1], x[2], x[3], x[4], x[5]); p.setState(x[0], x[1], x[2], x[3], x[4], x[5]);
} }
}); });
}) })
} }
console.log(Simulation.h);
} }
......
...@@ -9,12 +9,15 @@ document.addEventListener('keydown', (event: KeyboardEvent) => { ...@@ -9,12 +9,15 @@ document.addEventListener('keydown', (event: KeyboardEvent) => {
if (event.key === 'e') { if (event.key === 'e') {
SimulationData.switchParticlePhysics(); SimulationData.switchParticlePhysics();
} else if (event.key === 'w') {
SimulationData.paused = !SimulationData.paused;
} }
}); });
class SimulationData implements PipelineGUI, PipelineRenderable { class SimulationData implements PipelineGUI, PipelineRenderable {
public static container: THREE.Group = new THREE.Group(); public static container: THREE.Group = new THREE.Group();
public static paused: boolean = false;
// cloth settings // cloth settings
public static readonly MAX_TOTAL_MASS = 10; public static readonly MAX_TOTAL_MASS = 10;
...@@ -52,19 +55,15 @@ class SimulationData implements PipelineGUI, PipelineRenderable { ...@@ -52,19 +55,15 @@ class SimulationData implements PipelineGUI, PipelineRenderable {
public static sColorBend: number = 0x0000ff; public static sColorBend: number = 0x0000ff;
// integrator // integrator
public static readonly STEP_SIZE = {Fixed:0, Adaptive:1};
public static readonly INTEGRATORS: Integrator[] = [ public static readonly INTEGRATORS: Integrator[] = [
new EulerIntegrator(), new EulerIntegrator(),
new MidpointIntegrator(), new MidpointIntegrator(),
new RungeKuttaIntegrator new RungeKuttaIntegrator()
] ]
public static integratorType: number = 2; public static integrationMethod: number = 2;
public static stepSize: number = 0; public static stepSizeMethod: number = 0;
public static delta: number = 1e-3;
public static setEuler = () => { SimulationData.integratorType = 0; }
public static setMidpoint = () => { SimulationData.integratorType = 1; }
public static setRungeKutta = () => { SimulationData.integratorType = 2; }
public static resetExternalForce(): void { public static resetExternalForce(): void {
SimulationData.externalForce.x = 0; SimulationData.externalForce.x = 0;
...@@ -347,6 +346,7 @@ class SimulationData implements PipelineGUI, PipelineRenderable { ...@@ -347,6 +346,7 @@ class SimulationData implements PipelineGUI, PipelineRenderable {
const general = gui.addFolder('General'); const general = gui.addFolder('General');
general.add(SimulationData, 'paused').listen().name("Pause Simulation");
general.add(SimulationData, 'totalMass', SimulationData.MIN_TOTAL_MASS, SimulationData.MAX_TOTAL_MASS) general.add(SimulationData, 'totalMass', SimulationData.MIN_TOTAL_MASS, SimulationData.MAX_TOTAL_MASS)
.step(0.1).name('Total Mass').onChange(SimulationData.changedParticleSize); .step(0.1).name('Total Mass').onChange(SimulationData.changedParticleSize);
general.add(SimulationData, 'width', SimulationData.MIN_WIDTH, SimulationData.MAX_WIDTH) general.add(SimulationData, 'width', SimulationData.MIN_WIDTH, SimulationData.MAX_WIDTH)
...@@ -359,9 +359,9 @@ class SimulationData implements PipelineGUI, PipelineRenderable { ...@@ -359,9 +359,9 @@ class SimulationData implements PipelineGUI, PipelineRenderable {
const integration = gui.addFolder('Integrator'); const integration = gui.addFolder('Integrator');
integration.add(SimulationData, 'setEuler').name('Euler').onChange(()=>{console.log(SimulationData.integratorType)}); integration.add(SimulationData, 'integrationMethod', { 'Euler': 0, 'Midpoint': 1, 'Runge-Kutta': 2 }).name('Method');
integration.add(SimulationData, 'setMidpoint').name('Midpoint').onChange(()=>{console.log(SimulationData.integratorType)}); integration.add(SimulationData, 'stepSizeMethod', { 'Fixed': 0, 'Adaptive': 1 }).name('Step Size').onChange((e) => { SimulationData.stepSizeMethod = e; console.log(e);});
integration.add(SimulationData, 'setRungeKutta').name('Runge Kutta (4th order)').onChange(()=>{console.log(SimulationData.integratorType)}); integration.add(SimulationData, 'delta', 1e-10, 1e-6, 1e-10).name('Error Max');
integration.open(); integration.open();
const extForce = general.addFolder('External Force'); const extForce = general.addFolder('External Force');
......
abstract class Integrator { abstract class Integrator {
abstract step(x: Array<number>, h: number, f: (x: Array<number>) => Array<number>): Array<number>; abstract step(x: Array<number>, h: number, f: (x: Array<number>) => Array<number>): Array<number>;
public static adaptiveStep(integrator: Integrator, x: Array<number>, h: number,
f: (x: Array<number>) => Array<number>, errorMax: number): number {
let xFull = integrator.step(x, h, f);
let xHalf = integrator.step(x, h / 2, f);
let error: Array<number> = [];
for (let i = 0; i < xFull.length; i++) {
error.push(Math.abs(xFull[i] - xHalf[i]));
}
let errorApproxMax = error.reduce((a, b) => Math.max(a, b));
return h * Math.sqrt(errorMax / errorApproxMax);
}
} }
class EulerIntegrator extends Integrator { class EulerIntegrator extends Integrator {
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment