# Math library

By: Justin Meiners

Vector math and linear algebra library. This code should follow ES5 for maximum browser compatibility. Keep performance in mind. Generally, prefer loops to avoid unnecessary array allocations.

## 1. Utilities

``````function lerp(a, b, t) {
return (1 - t) * a + t * b;
};

function clamp(low, high, t) {
if (t < low) {
return low;
} else if (t > high) {
return high;
} else {
return t;
}
};
``````

Used by 1

## 2. Vector

This is a 2D vector library.

Class definition

``````function Vec(x, y) {
this.x = x;
this.y = y;
};

Vec.fromAngle = function(angle) {
return new Vec(Math.cos(angle), Math.sin(angle));
};

Vec.zero = function() {
return new Vec(0, 0);
};

Vec.prototype.copy = function() {
return new Vec(this.x, this.y);
};

@{vec arithmetic}
@{vec algorithms}
``````

Used by 1

### Arithematic operations

Generally, the prototypes variants modify the vector in place.

``````Vec.prototype.dot = function(b) {
return this.x * b.x + this.y * b.y;
};

Vec.prototype.lenSqr = function() {
return this.dot(this);
};

Vec.prototype.len = function() {
return Math.sqrt(this.lenSqr());
};

Vec.prototype.orthogonal = function() {
return new Vec(-this.y, this.x);
};

this.x += b.x;
this.y += b.y;
return this;
};

Vec.prototype.negate = function() {
this.x = -this.x;
this.y = -this.y;
return this;
};

Vec.prototype.sub = function(b) {
this.x -= b.x;
this.y -= b.y;
return this;
};

Vec.prototype.normalized = function() {
return Vec.scale(this, 1 / this.len());
};

};

Vec.sub = function(a, b) {
return a.copy().sub(b);
};

Vec.scale = function(a, s) {
return new Vec(a.x * s, a.y * s);
};

Vec.distSqr = function(a, b) {
return Vec.sub(a, b).lenSqr();
};

Vec.dist = function(a, b) {
return Math.sqrt(Vec.distSqr(a, b));
};

Vec.negate = function(v) {
return v.copy().negate();
};

Vec.min = function(a, b) {
return new Vec(Math.min(a.x, b.x), Math.min(a.y, b.y));
};

Vec.max = function(a, b) {
return new Vec(Math.max(a.x, b.x), Math.max(a.y, b.y));
};
``````

Used by 1

### Algorithms

``````Vec.prototype.inBounds = function(min, max) {
return this.x >= min.x && this.y >= min.y &&
this.x <= max.x && this.y <= max.y;
};

Vec.prototype.inCircle = function(o, r) {
return Vec.distSqr(this, o) < r * r;
};

Vec.lerp = function(a, b, t) {
return new Vec(lerp(a.x, b.x, t), lerp(a.y, b.y, t));
};

Vec.centroid = function(points) {
var sum = points.reduce(function(acc, x) {
});

return Vec.scale(sum, 1 / points.length);
};

Vec.pathLen = function(path) {
var sum = 0;
for (var i = 1; i < path.length; ++i) {
sum += Vec.dist(path[i], path[i-1]);
}
return sum;
};

Vec.bounds = function(points) {
return [points.reduce(Vec.min), points.reduce(Vec.max)];
};
``````

Used by 1

## 3. 2x2 Matrix

``````function Matrix() {
this.m = new Float64Array(2 * 2);
};

Matrix.from = function(a, b, c, d) {
var A = new Matrix();
A.m[0] = a; A.m[1] = b;
A.m[2] = c; A.m[3] = d;
return A;
};

Matrix.fromScale = function(sx, sy) {
return Matrix.from(sx, 0, 0, sy);
};

Matrix.fromAngle = function(angle) {
var x = Math.cos(angle);
var y = Math.sin(angle);
return Matrix.from(x, -y, y, x);
};

Matrix.prototype.index = function(row, col) {
return row * 2 + col;
};

Matrix.prototype.get = function(row, col) {
return this.m[this.index(row, col)];
};

Matrix.prototype.set = function(row, col, x) {
this.m[this.index(row, col)] = x;
};

@{matrix operations}
``````

Used by 1

``````Matrix.prototype.transform = function(v) {
return new Vec(
this.m[0] * v.x + this.m[1] * v.y,
this.m[2] * v.x + this.m[3] * v.y
);
};

for (var i = 0; i < 4; ++i) {
this.m[i] += other.m[i];
};
return this;
};

Matrix.prototype.scale = function(scale) {
for (var i = 0; i < 4; ++i) {
this.m[i] *= scale;
};
return this;
};

Matrix.mul = function(A, B) {
return Matrix.from(
A.m[0] * B.m[0] + A.m[1] * B.m[2],
A.m[0] * B.m[1] + A.m[1] * B.m[3],
A.m[2] * B.m[0] + A.m[3] * B.m[2],
A.m[2] * B.m[1] + A.m[3] * B.m[3]
);
};

Matrix.prototype.det = function() {
return this.m[0] * this.m[3] - this.m[1] * this.m[2];
};

Matrix.inverse = function(A) {
var s = 1.0 / A.det();
return Matrix.from(A.m[3], -A.m[1], -A.m[2], A.m[0]).scale(s);
};

Matrix.prototype.copy = function() {
var copy = new Matrix();
copy.m.set(this.m);
return copy;
};
``````

Used by 1

## 4. Affine transform

``````function AffineTransform() {
this.m = Matrix.from(1, 0, 0, 1);
this.t = new Vec(0, 0);
};

AffineTransform.prototype.transform = function(v) {
};

AffineTransform.inverse = function(a) {
var b = new AffineTransform();
b.m = Matrix.inverse(a.m);
b.t = b.m.transform( a.t.copy().negate());
return b;
};
``````

Used by 1

## 5. Optimizer

This is an implementation of the Nelder-mead algorithm.

``````function multivarOptimize(initial, costFunc, options) {
var maxIterations = options.maxIterations || 10000;
var tolerance = options.tolerance || 0.001;

function predicate(iterations, image) {
return iterations > maxIterations || image < tolerance;
};

var n = initial.length;

function guessSizes(initial) {
var sizes = new Float64Array(n);
for (var j = 0; j < n; ++j) {
sizes[j] = (0.05 * Math.abs(initial[j])) + 0.00025;
}
return sizes;
};

var simplex = new Simplex(n);
simplex.positionAround(initial, guessSizes(initial));
return simplex.optimize(costFunc, predicate, options);
};

@{simplex}
``````

Used by 1

``````function Simplex(n) {
this.dimension = n;
this.vertices = [];
for (var i = 0; i < (n + 1); ++i) {
this.vertices.push(new SimplexVertex(n));
}
};

function SimplexVertex(n) {
this.point = new Float64Array(n);
this.image = 0;
};

SimplexVertex.prototype.copyFrom = function(other) {
this.point.set(other.point);
this.image = other.image;
};

SimplexVertex.prototype.swapWith = function(other) {
var ptemp = this.point;
var itemp = this.image;

this.point = other.point;
this.image = other.image;
other.point = ptemp;
other.image = itemp;
};

Simplex.prototype.positionAround = function(point, sizes) {
var n = this.dimension;
for (var i = 0; i < (n + 1); ++i) {
for (var j = 0; j < n; ++j) {
this.vertices[i].point[j] = point[j] + ((i - 1 == j) ? sizes[j] : 0.0);
}
}
};

Simplex.prototype.sort = function() {
this.vertices.sort(function(a, b) {
return a.image - b.image;
});
};

Simplex.prototype.getCentroid = function(output) {
var n = this.dimension;
// last point is not included
for (var j = 0; j < n; ++j) {
var x = this.vertices[0].point[j];
for (var i = 1; i < n; ++i) {
x += this.vertices[i].point[j];
}
output.point[j] = x / n;
}
};

Simplex.prototype.getNewPoint = function(centroid, lambda, output) {
var n = this.dimension;
for (var j = 0; j < n; ++j) {
output.point[j] = (1.0 + lambda) * centroid.point[j] - lambda * this.vertices[n].point[j];
};
};

Simplex.prototype.shrink = function(sigma) {
var n = this.dimension;
for (var i = 1; i < (n + 1); ++i) {
for (var j = 0; j < n; ++j) {
var x = this.vertices[0].point[j];
this.vertices[i].point[j] = x + sigma * (this.vertices[i].point[j] - x);
}
}
};

@{simplex optimizer}
``````

Used by 1

Nelder-Mead optimization iteratively transforms a simplex, moving it in the direction where the cost function decreases and contracting it around a local minimum when it finds one.

At each step, one of four transformations is applied: “reflect”, “expand”, “contract”, or “shrink”. The magnitude of each transformation is controlled by a coefficient.

``````Simplex.prototype.computeVertexImages = function(costFunc) {
for (var i = 0; i < (this.dimension + 1); ++i) {
this.vertices[i].image = costFunc(this.vertices[i].point);
}
this.sort();
};

Simplex.prototype.optimize = function(costFunc, stopPredicate, options) {
var rho = options.rho || 1.0; // coefficient for "reflect"
var chi = options.chi || 2.0; // coefficient for "expand"
var gamma = options.gamma || 0.5; // coefficient for "contract"
var sigma = options.sigma || 0.5; // coefficient for "shrink"
var debug = options.debug == undefined ? false : options.debug;

var n = this.dimension;

var pointR = new SimplexVertex(n);
var pointE = new SimplexVertex(n);
var pointC = new SimplexVertex(n);
var centroid = new SimplexVertex(n);

this.computeVertexImages(costFunc);
this.getCentroid(centroid);

var iterations = 1;

while (!stopPredicate(iterations, this.vertices[0].image)) {
var shrink = false;
this.getNewPoint(centroid, rho, pointR);
pointR.image = costFunc(pointR.point);

if (pointR.image < this.vertices[0].image) {
this.getNewPoint(centroid, rho * chi, pointE);
pointE.image = costFunc(pointE.point);

if (pointE.image < pointR.image) {
if (debug) console.log("expand");
this.vertices[n].copyFrom(pointE);
} else {
if (debug) console.log("reflect");
this.vertices[n].copyFrom(pointR);
}
} else {
if (pointR.image < this.vertices[n - 1].image) {
if (debug) console.log("reflect");
this.vertices[n].copyFrom(pointR);
} else {
if (pointR.image < this.vertices[n].image) {
this.getNewPoint(centroid, rho * gamma, pointC);
pointC.image = costFunc(pointC.point);

if (pointC.image <= pointR.image) {
if (debug) console.log("contract out");
this.vertices[n].copyFrom(pointC);
} else {
if (debug) console.log("shrink");
shrink = true;
}
} else {
this.getNewPoint(centroid, -gamma, pointC);
pointC.image = costFunc(pointC.point);

if (pointC.image < this.vertices[n].image) {
if (debug) console.log("contract inside");
this.vertices[n].copyFrom(pointC);
} else {
if (debug) console.log("shrink");
shrink = true;
}
}
}
}

if (shrink) {
this.shrink(sigma);
this.computeVertexImages(costFunc);
} else {
for (var i = n - 1; i >= 0 && this.vertices[i + 1].image < this.vertices[i].image; i--) {
this.vertices[i + 1].swapWith(this.vertices[i]);
}

}
this.getCentroid(centroid);
++iterations;

if (debug) console.log("best: ", this.vertices[0].point.join(", "), "score: ", this.vertices[0].image);
};

return {
vars: this.vertices[0].point,
cost: this.vertices[0].image,
iterations: iterations
};
};
``````

Used by 1

The Ackley function is a commonly used for testing optimization algorithms.

``````function makeAckley(a, b, c) {
return function(vars) {
var sumSquares = vars.reduce(function(total, x) {
}, 0);

var sumCos = vars.reduce(function(total, x) {
}, 0);

var n = vars.length;
return -a * Math.exp(-b * Math.sqrt(sumSquares /  n)) - Math.exp(sumCos / n) + a + Math.exp(1);
};
};

function testOptimizer() {
var f = makeAckley(20.0, 0.2, 2.0 * Math.PI);

var options = {
debug: true,
maxIterations: 100,
tolerance: 0.001
};

console.log(multivarOptimize([-2.1, -3.04, 4.5], f, options));
};
``````

Used by 1

``````// Created by Justin Meiners (2022)