@@ Line 24: / Line 24: @@
 This is the {{w|JavaScript}} code I use to find purely-consistent EDOs by odd limit.
+=== Original ===
 <syntaxhighlight lang="javascript">
 function et_error(interval, et)
@@ Line 53: / Line 54: @@
          i++;
      }
+}
+</syntaxhighlight>
+=== Current ===
+(Note: This code was written with the assistance of Claude Opus 5, an artificial intelligence tool.)
+<syntaxhighlight lang="javascript">
+/**
+ * Generates a list of equal divisions of the octave that approximate the harmonic series well.
+ *
+ * An EDO qualifies when every odd harmonic from 3 to `h` falls within `thresh` of one of its steps, measured as a fraction of
+ * a step. Even harmonics are skipped: they are octave transpositions of odd ones, and the octave is exact in any EDO.
+ *
+ * Searches upward without an upper bound, logging each hit with the elapsed time, and returns once `max_edos` have been
+ * found. Qualifying EDOs always exist, but they thin out fast, so high `h` can take a very long time to find.
+ *
+ * @param {Number} h - Highest harmonic to check. Internally this value is rounded up to the current or next highest odd
+ *                     number. This means inputting 32 and 33 will both check through 33. Anything below 3 accepts every EDO.
+ *                     Usually this is set to one less than a power of two (e.g. 7, 15, 31, or 63)
+ * @param {Number} max_edos - How many EDOs to find before stopping.
+ * @param {Number} [start=1] - Lowest EDO to consider; "1 equal division of the octave" is simply whole-octave steps.
+ * @param {Number} [thresh=0.25] - Tolerance as a fraction of a step. Defaults to 0.25
+ * @returns {Number[]} - The list of EDOs found, in ascending order.
+ */
+function min_et(h, max_edos, start = 1, thresh = 0.25)
+{
+    h = Math.round(Math.abs(h));
+    h = h + 1 - (h % 2);
+    start = Math.max(1, Math.round(start));
+    // Odd primes <=h
+    let comp = new Int16Array(h + 1), primes = [];
+    for (let p = 3; p <= h; p += 2)
+    {
+        if (comp[p])
+        {
+            continue;
+        }
+        primes.push(p);
+        for (let q = p * p; q <= h; q += 2 * p)
+        {
+            comp[q] = 1;
+        }
+    }
+    let m = primes.length, out = [];
+    if (m === 0)
+    {
+        for (let n = start; out.length < max_edos; n++)
+        {
+            out.push(n);
+        }
+        return out;
+    }
+    // Per-prime tolerance: p, p^2, p^3, ... p^t collapse to |e_p| < thresh/t
+    let logs = new Float64Array(m), tol = new Float64Array(m);
+    for (let i = 0; i < m; i++)
+    {
+        logs[i] = Math.log2(primes[i]);
+        let t = 0;
+        for (let q = primes[i]; q <= h; q *= primes[i])
+        {
+            t++;
+        }
+        tol[i] = thresh / t;
+    }
+    // Composites with >= 2 distinct factors, bucketed by largest factor
+    let mixed = Array.from({length: m}, () => []);
+    for (let k = 9; k <= h; k += 2)
+    {
+        let n = k, top = -1, d = 0;
+        let term = [];
+        for (let i = 0; i < m && n > 1; i++)
+        {
+            let p = primes[i];
+            if (n % p)
+            {
+                continue;
+            }
+            let e = 0;
+            while (n % p === 0)
+            {
+                n /= p;
+                e++;
+            }
+            term.push(i, e);
+            top = i;
+            d++;
+        }
+        if (d >= 2)
+        {
+            mixed[top].push(Int32Array.from(term));
+        }
+    }
+    // 3-gap data for the 3-limit filter
+    let err = new Float64Array(m);
+    let a = logs[0], b = tol[0], L = 2 * b;      // Target arc has length L
+    let jm = 0, dm = 0, jn = 0, dn = 0;
+    for (let j = 1; ; j++)
+    {
+        let t = j * a; t -= Math.floor(t);
+        if (!jm && t < L)
+        {
+            jm = j; dm = t;
+        }
+        if (!jn && t > 1 - L)
+        {
+            jn = j; dn = 1 - t;
+        }
+        if (jm && jn)
+        {
+            break;
+        }
+    }
+    let jmn = jm + jn, dmn = dm - dn, lim = L - dm;
+    let N = start, x = start * a + b; x -= Math.floor(x);
+    while (x >= L)
+    {
+        N++; x += a - 1;
+        if (x >= 1)
+        {
+            x -= 1;
+        }
+    }
+    // Search
+    let t0 = Date.now();
+    let steps = 0;
+    while (out.length < max_edos)
+    {
+        err[0] = x - b;
+        let ok = true;
+        for (let i = 1; i < m; i++)
+        {
+            let y = N * logs[i], e = y - Math.round(y);
+            if (e >= tol[i] || e <= -tol[i])
+            {
+                ok = false;
+                break;
+            }
+            err[i] = e;
+            let g = mixed[i];
+            for (let t = 0; t < g.length; t++)
+            {
+                let term = g[t];
+                let s = 0;
+                for (let q = 0; q < term.length; q += 2)
+                {
+                    s += term[q + 1] * err[term[q]];
+                }
+                s -= Math.round(s);
+                if (s >= thresh || s <= -thresh)
+                {
+                    ok = false;
+                    break;
+                }
+            }
+            if (!ok)
+            {
+                break;
+            }
+        }
+        if (ok)
+        {
+            out.push(N);
+            console.log(N + " | " + (Date.now() - t0) / 1000 + "s");
+        }
+        // Jump to the next N that passes the 3-limit test
+        if (x < lim)
+        {
+            N += jm;
+            x += dm;
+        }
+        else if (x >= dn)
+        {
+            N += jn;
+            x -= dn;
+        }
+        else
+        {
+            N += jmn; x += dmn;
+        }
+        if ((++steps & 0xFFFFF) === 0)
+        {
+            x = N * a + b;
+            x -= Math.floor(x);
+        }
+    }
+    return out;
 }
 </syntaxhighlight>

User:ArrowHead294/Purely consistent EDOs by odd limit/Appendix: Difference between revisions

ArrowHead294 (talk | contribs)