|
const factors = [3, 5, 7, 11, 13, 17, 19, 23, 29, 31].map(x => ({ |
|
factor: x, |
|
// the following polynomial was justified more by experiment than anything else, |
|
// the goal is that a 1/n^3 progression should mean that we converge on the |
|
// basics quickly and only get weak influence from the 11:2 or 13:2 ratios. |
|
weight: 4 / ((x - 3) ** 3 + 24) // 1/3, 1/4, 1/11, 1/67, 1/128... |
|
})); |
|
|
|
/** |
|
* @type |
|
* a combo like [1, 1, -2, 0, 0] means that we want 3*5/(7**2), powers of the |
|
* digits. we try to assemble all combos where sum(map(abs)(seq)) <= 5. We use |
|
* a set of strings so that we do not have to worry about counting the same |
|
* combo twice. The choice of 5 is to get the syntonic comma 81/80 in there, |
|
* hopefully. |
|
*/ |
|
const comboStrings = new Set([]); |
|
/** |
|
* To do this we use "choosers" which choose an index and add +/- 1 to it. |
|
* @returns {Iterable<[number, number]>} |
|
*/ |
|
function* chooser() { |
|
for (let i = 0; i < factors.length; i++) { |
|
yield [i, 1]; |
|
yield [i, -1]; |
|
} |
|
} |
|
const blankCombo = factors.map(() => 0); |
|
for (const [n_a, k_a] of chooser()) { |
|
for (const [n_b, k_b] of chooser()) { |
|
for (const [n_c, k_c] of chooser()) { |
|
for (const [n_d, k_d] of chooser()) { |
|
const combo = blankCombo.slice(0); |
|
combo[n_a] += k_a; |
|
combo[n_b] += k_b; |
|
combo[n_c] += k_c; |
|
combo[n_d] += k_d; |
|
comboStrings.add(combo.join(" ")); |
|
} |
|
} |
|
} |
|
} |
|
for (const [n_a, k_a] of chooser()) { |
|
for (const [n_b, k_b] of chooser()) { |
|
for (const [n_c, k_c] of chooser()) { |
|
for (const [n_d, k_d] of chooser()) { |
|
for (const [n_e, k_e] of chooser()) { |
|
const combo = blankCombo.slice(0); |
|
combo[n_a] += k_a; |
|
combo[n_b] += k_b; |
|
combo[n_c] += k_c; |
|
combo[n_d] += k_d; |
|
combo[n_e] += k_e; |
|
comboStrings.add(combo.join(" ")); |
|
} |
|
} |
|
} |
|
} |
|
} |
|
// if this got in there who cares |
|
comboStrings.delete(blankCombo.join(" ")); |
|
|
|
// Now we want to bring these combinations to within one octave, and we want to |
|
// also not double-count our errors: log(a/b) = -log(b/a) and then when we add |
|
// to get it to live within the interval (0, 1) we will find that b/a is just as |
|
// far away from its best rational approximant as a/b was. So we can solve both |
|
// problems at once by just only considering the fractions that start out |
|
// positive. I would also like to formally bound this by the weight of the 3^5 |
|
// term, keeping only 4th-order terms with higher weight than that. |
|
|
|
/** |
|
* @type {Array<{factor: number, weight: number, frac: [number, number], prob: number, logFrac: number}>} |
|
*/ |
|
const combinations = []; |
|
const ignoreWeight = factors[0].weight ** 6; |
|
|
|
for (const comboString of comboStrings) { |
|
/** |
|
* @type {number[]} |
|
*/ |
|
const combos = comboString.split(" ").map(x => Number(x)); |
|
let numer = 1, |
|
denom = 1; |
|
let totalWeight = 1; |
|
for (let i = 0; i < factors.length; i++) { |
|
const { factor, weight } = factors[i]; |
|
const n = combos[i]; |
|
if (n > 0) { |
|
numer *= factor ** n; |
|
totalWeight *= weight ** n; |
|
} else if (n < 0) { |
|
const m = -n; |
|
denom *= factor ** m; |
|
totalWeight *= weight ** m; |
|
} |
|
} |
|
if (totalWeight > ignoreWeight && numer > denom) { |
|
let factor = numer / denom; |
|
while (factor > 2) { |
|
factor /= 2; |
|
denom *= 2; |
|
} |
|
combinations.push({ |
|
factor, |
|
weight: totalWeight, |
|
fraction: [numer, denom], |
|
prob: 0, |
|
note: Math.log(numer / denom) / Math.log(2) |
|
}); |
|
} |
|
} |
|
|
|
combinations.sort((x, y) => y.weight - x.weight); |
|
const sumWeights = combinations.reduce((acc, x) => acc + x.weight, 0); |
|
for (const combo of combinations) { |
|
combo.prob = combo.weight / sumWeights; |
|
} |
|
|
|
function centsOff(notesPerOctave, note) { |
|
const est = Math.round(note * notesPerOctave) / notesPerOctave; |
|
return 1200 * (est - note); |
|
} |
|
function error(notesPerOctave, noteProbs) { |
|
let squareSum = 0; |
|
for (const { note, prob } of noteProbs) { |
|
squareSum += prob * centsOff(notesPerOctave, note) ** 2; |
|
} |
|
return squareSum ** 0.5; |
|
} |
|
function roundedString(n, color) { |
|
let rounded = Math.round(n * 100) / 100; |
|
let out = String(rounded); |
|
if (out.indexOf(".") === -1) { |
|
out += ".0"; |
|
} |
|
const withLastDigit = /\.\d$/.exec(out) ? out + "0" : out; |
|
if (color && Math.abs(rounded) <= 6) { |
|
return `*${withLastDigit}*` |
|
} |
|
if (color && Math.abs(rounded) > 12) { |
|
return `**${withLastDigit}**` |
|
} |
|
return withLastDigit; |
|
} |
|
// n is the number of notes per octave in an equal temperament |
|
const output = [ |
|
["notes", "expect-err", "actual-err"].concat( |
|
combinations.slice(0, 12).map(x => x.fraction.join(":")) |
|
) |
|
]; |
|
output.push(output[0].map(() => "---")); |
|
const expect = [...Array(1e7)].map((_, i) => ({ |
|
note: (i + 0.5) / 1e7, |
|
prob: 1e-7 |
|
})); |
|
for (let notesPerOctave = 6; notesPerOctave < 50; notesPerOctave++) { |
|
const expectedError = error(notesPerOctave, expect); |
|
const actualError = error(notesPerOctave, combinations); |
|
if (actualError > expectedError) { |
|
continue; |
|
} |
|
output.push( |
|
[ |
|
String(notesPerOctave), |
|
roundedString(expectedError), |
|
roundedString(actualError) |
|
].concat( |
|
combinations |
|
.slice(0, 12) |
|
.map(x => roundedString(centsOff(notesPerOctave, x.note), true)) |
|
) |
|
); |
|
} |
|
const colLengths = output[0].map(() => 0); |
|
for (let i = 0; i < output.length; i++) { |
|
for (let c = 0; c < colLengths.length; c++) { |
|
colLengths[c] = Math.max( |
|
colLengths[c], |
|
output[i][c].length + (i > 0 ? 1 : 0) |
|
); |
|
} |
|
} |
|
for (let i = 0; i < output.length; i++) { |
|
const row = output[i]; |
|
console.log('|', row.join(" | "), '|'); |
|
} |