Replace tune2, tune3 and tune4 with a single search

documentation/BUILTIN.md

@@ -411,6 +411,9 @@ Transpose a matrix. For modal transposition see rotate().
 ### trunc(*interval*)
 Truncate value towards zero to the nearest integer.
 
+### tune(*vals*, *searchRadius*, *weights*)
+Attempt to combine the given vals into a more Tenney-Euclid optimal val. Weights are applied multiplicatively on top of Tenney weights of the subgroup basis.
+
 ### unshift(*interval*, *scale = $$*)
 Prepend an interval at the beginning of the current/given scale.
 
@@ -728,15 +731,6 @@ Obtain a copy of the current/given scale quantized to subharmonics of the given
 ### trap(*message*)
 Produce a function that fails with the given message when called.
 
-### tune2(*a*, *b*, *numIter = 1*, *weights = niente*)
-Find a combination of two vals that is closer to just intonation.
-
-### tune3(*a*, *b*, *c*, *numIter = 1*, *weights = niente*)
-Find a combination of three vals that is closer to just intonation.
-
-### tune4(*a*, *b*, *c*, *d*, *numIter = 1*, *weights = niente*)
-Find a combination of four vals that is closer to just intonation.
-
 ### u(*scale = ££*)
 Obtain a undertonal reflection of the popped/given overtonal scale.
 

documentation/tempering.md

@@ -267,11 +267,10 @@ const S = @2.7.11
 ### errorTE
 In order to compare the quality of vals w.r.t. just intonation SonicWeave provides the `errorTE` helper that measures [RMS TE error](https://en.xen.wiki/w/Tenney-Euclidean_temperament_measures#TE_error) in cents. Providing explicit subgroups is highly recommended so that irrelevant higher primes do not interfere with the measure.
 
-### tune2
-The helper `tune2` takes two vals and tries to find a combination that's closer to just intonation, effectively performing constrained Tenney-Euclidean optimization (CTE). E.g. `tune2([email protected], [email protected])` finds 31p. Given more iterations `tune2([email protected], [email protected], 7)` finds an equal temperament in the thousands that's virtually indistinguishable from the true meantone CTE tuning.
+### tune
+The helper `tune` takes an array of vals and tries to find a combination that's closer to just intonation, effectively performing constrained Tenney-Euclidean optimization (CTE). E.g. `tune([[email protected], [email protected]])` finds 31p. Given a large search radius `tune2([email protected], [email protected], 200)` finds an equal temperament in the thousands that's virtually indistinguishable from the true meantone CTE tuning.
 
-### tune3 and tune4
-The helpers `tune3` and `tune4` do the same but try combinations of 3 or 4 vals instead. Assuming the vals are linearly independent the process corresponds to rank-3 and rank-4 CTE optimization.
+With more than two linearly independent vals the process corresponds to higher rank CTE optimization.
 
 ### Discovering vals
 Every just intonation subgroup has a generalized patent val sequence that dances around approximations to the just intonation point i.e. the perfectly pure tuning.

examples/keen12.sw

@@ -18,4 +18,4 @@ labelAbsoluteFJS
 (* See https://en.xen.wiki/w/Val#Sparse_Offset_Val_notation for the alternative notation for 12d. *)
 (* Doesn't actually matter here because we're stacking fifths, but communicates the intended interpretation. *)
 (* Same as 56@. *)
-tune(12[v7]@, 22@, 2)
+tune([12[v7]@, 22@], 2)

src/__tests__/temper.spec.ts

@@ -1,8 +1,13 @@
 import {describe, expect, it} from 'vitest';
-import {TuningMap, combineTuningMaps, vanishCommas} from '../temper';
+import {
+  TuningMap,
+  combineTuningMaps,
+  intCombineTuningMaps,
+  vanishCommas,
+} from '../temper';
 import {LOG_PRIMES, applyWeights, dot, unapplyWeights} from 'xen-dev-utils';
 
-describe('Val combination optimizer', () => {
+describe('Mapping combination optimizer', () => {
   it('computes TE meantone', () => {
     const p12: TuningMap = [12, 19, 28].map(c => (c / 12) * LOG_PRIMES[0]);
     const p19: TuningMap = [19, 30, 44].map(c => (c / 19) * LOG_PRIMES[0]);
@@ -47,6 +52,22 @@ describe('Val combination optimizer', () => {
   });
 });
 
+describe('Val combination search', () => {
+  it('combines 12p and 19p into 31p', () => {
+    const p12: TuningMap = [12, 19 / Math.log2(3), 28 / Math.log2(5)];
+    const p19: TuningMap = [19, 30 / Math.log2(3), 44 / Math.log2(5)];
+    const coeffs = intCombineTuningMaps([1, 1, 1], [p12, p19], 1);
+    expect(coeffs).toEqual([1, 1]);
+  });
+
+  it('combines 5p and 7p into 31p', () => {
+    const p5: TuningMap = [5, 8 / Math.log2(3), 12 / Math.log2(5)];
+    const p7: TuningMap = [7, 11 / Math.log2(3), 16 / Math.log2(5)];
+    const coeffs = intCombineTuningMaps([1, 1, 1], [p5, p7], 3);
+    expect(coeffs).toEqual([2, 3]);
+  });
+});
+
 describe('Comma vanisher', () => {
   it('computes TE meantone', () => {
     const syntonic = [-4, 4, -1];

src/parser/__tests__/expression.spec.ts

@@ -2703,4 +2703,41 @@ describe('SonicWeave expression evaluator', () => {
     );
     expect(fraction).toBe('11/9');
   });
+
+  it('can combine two vals to approach the JIP', () => {
+    const thirtyOne = evaluate('tune([[email protected], [email protected]])') as Val;
+    expect(thirtyOne.value.toIntegerMonzo()).toEqual([31, 49, 72]);
+  });
+
+  it('can combine two vals to approach the JIP (Wilson metric)', () => {
+    const eighty = evaluate(
+      'tune([[email protected], [email protected]], 5, [1 % 2, 3 /_ 2 % 3, 5 /_ 2 % 5])'
+    ) as Val;
+    expect(eighty.value.toIntegerMonzo()).toEqual([126, 200, 293]);
+  });
+
+  it('can combine three vals to approach the JIP', () => {
+    const fourtyOne = evaluate('tune([[email protected], [email protected], [email protected]])') as Val;
+    expect(fourtyOne.value.toIntegerMonzo()).toEqual([41, 65, 95, 115]);
+  });
+
+  it('can combine four vals to approach the JIP', () => {
+    const val = evaluate('tune([[email protected], [email protected], [email protected], [email protected]])') as Val;
+    expect(val.value.toIntegerMonzo()).toEqual([72, 114, 167, 202, 249]);
+  });
+
+  it("doesn't move from 31p when tuned with 5p", () => {
+    const p31 = evaluate('str(tune([[email protected], [email protected]]))');
+    expect(p31).toBe('<31 49 72]');
+  });
+
+  it('moves from 31p when tuned with 5p if given a large enough radius', () => {
+    const p31 = evaluateExpression('str(tune([[email protected], [email protected]], 6))');
+    expect(p31).toBe('<191 302 444]');
+  });
+
+  it('tunes close to CTE meantone if given a large enough radius', () => {
+    const cents = evaluate('str(cents(3/2 tmpr tune([[email protected], [email protected]], 200), 4))');
+    expect(cents).toEqual('697.2143');
+  });
 });

src/parser/__tests__/stdlib.spec.ts

@@ -6,7 +6,7 @@ import {
   getSourceVisitor,
   parseAST,
 } from '../../parser';
-import {Interval, Val} from '../../interval';
+import {Interval} from '../../interval';
 import {builtinNode, track} from '../../stdlib';
 import {Fraction} from 'xen-dev-utils';
 
@@ -366,30 +366,6 @@ describe('SonicWeave standard library', () => {
     );
   });
 
-  it('can combine two vals to approach the JIP', () => {
-    const thirtyOne = evaluateExpression('tune2([email protected], [email protected])') as Val;
-    expect(thirtyOne.value.toIntegerMonzo()).toEqual([31, 49, 72]);
-  });
-
-  it('can combine two vals to approach the JIP (Wilson metric)', () => {
-    const eighty = evaluateExpression(
-      'tune2([email protected], [email protected], 2, [log(2)/2, log(3)/3, log(5)/5])'
-    ) as Val;
-    expect(eighty.value.toIntegerMonzo()).toEqual([126, 200, 293]);
-  });
-
-  it('can combine three vals to approach the JIP', () => {
-    const fourtyOne = evaluateExpression('tune3([email protected], [email protected], [email protected])') as Val;
-    expect(fourtyOne.value.toIntegerMonzo()).toEqual([41, 65, 95, 115]);
-  });
-
-  it('can combine four vals to approach the JIP', () => {
-    const val = evaluateExpression(
-      'tune4([email protected], [email protected], [email protected], [email protected])'
-    ) as Val;
-    expect(val.value.toIntegerMonzo()).toEqual([94, 149, 218, 264, 325]);
-  });
-
   // Remember that unison (0.0 c) is implicit in SonicWeave
 
   it('has harmonic segments sounding downwards', () => {
@@ -1705,19 +1681,4 @@ describe('SonicWeave standard library', () => {
       '[[email protected], [email protected], [email protected], [email protected], [email protected]]'
     );
   });
-
-  it("doesn't move from 31p when tuned with 5p", () => {
-    const p31 = evaluateExpression('str(tune2([email protected], [email protected]))');
-    expect(p31).toBe('<31 49 72]');
-  });
-
-  it('moves from 31p when tuned with 5p if given enough time', () => {
-    const p31 = evaluateExpression('str(tune2([email protected], [email protected], 3))');
-    expect(p31).toBe('<191 302 444]');
-  });
-
-  it('tunes close to CTE meantone if given enough time', () => {
-    const scale = expand('cents(3/2 tmpr tune2([email protected], [email protected], 7), 4)');
-    expect(scale).toEqual(['697.2143']);
-  });
 });

src/stdlib/builtin.ts

@@ -80,7 +80,12 @@ import {
 } from './public';
 import {scaleMonzos} from '../diamond-mos';
 import {valToSparseOffset, valToWarts} from '../warts';
-import {TuningMap, combineTuningMaps, vanishCommas} from '../temper';
+import {
+  TuningMap,
+  combineTuningMaps,
+  intCombineTuningMaps,
+  vanishCommas,
+} from '../temper';
 const {version: VERSION} = require('../../package.json');
 
 // === Library ===
@@ -1770,6 +1775,55 @@ function nextGPV(
 nextGPV.__doc__ = 'Obtain the next generalized patent val in the sequence.';
 nextGPV.__node__ = builtinNode(nextGPV);
 
+function tune(
+  this: ExpressionVisitor,
+  vals: SonicWeaveValue,
+  searchRadius: SonicWeaveValue,
+  weights: SonicWeaveValue
+) {
+  if (!Array.isArray(vals)) {
+    throw new Error('An array of vals is required.');
+  }
+  if (!vals.length) {
+    throw new Error('At least one val is required.');
+  }
+  for (const val of vals) {
+    if (!(val instanceof Val)) {
+      throw new Error('An array of vals is required.');
+    }
+  }
+  const vs = vals as Val[];
+  const basis = vs[0].basis;
+  for (const val of vs.slice(1)) {
+    if (!val.basis.equals(basis)) {
+      throw new Error('Vals bases must agree in tuning.');
+    }
+  }
+  const radius =
+    searchRadius === undefined ? 1 : upcastBool(searchRadius).toInteger();
+  this.spendGas((2 * radius + 1) ** vals.length);
+  const jip = basis.value.map(m => m.totalCents());
+  const ws = valWeights(weights, basis.size);
+  const wvals = vs.map(val =>
+    applyWeights(
+      unapplyWeights(
+        basis.value.map(m => m.dot(val.value).valueOf()),
+        jip
+      ),
+      ws
+    )
+  );
+  const coeffs = intCombineTuningMaps(ws, wvals, radius);
+  let result = vs[0].mul(fromInteger(coeffs[0]));
+  for (let i = 1; i < coeffs.length; ++i) {
+    result = result.add(vs[i].mul(fromInteger(coeffs[i])));
+  }
+  return result.abs();
+}
+tune.__doc__ =
+  'Attempt to combine the given vals into a more Tenney-Euclid optimal val. Weights are applied multiplicatively on top of Tenney weights of the subgroup basis.';
+tune.__node__ = builtinNode(tune);
+
 function tenneyHeight(
   this: ExpressionVisitor,
   interval: SonicWeaveValue
@@ -3080,6 +3134,7 @@ export const BUILTIN_CONTEXT: Record<string, Interval | SonicWeaveFunction> = {
   TE,
   errorTE,
   nextGPV,
+  tune,
   tenneyHeight,
   wilsonHeight,
   respell,

src/stdlib/prelude.ts

@@ -333,133 +333,6 @@ riff enumerate(array = $$) {
   return [[i, array[i]] for i in array];
 }
 
-riff tune2(a, b, numIter = 1, weights = niente) {
-  "Find a combination of two vals that is closer to just intonation.";
-
-  let error = (v => errorTE(v, weights, true));
-  if (isFunction(weights)) {
-    error = weights;
-  }
-
-  numIter = real(numIter) + 1r;
-  while (--numIter) {
-    const combos = [
-      a,
-      b,
-      a + b,
-      2*a + b,
-      2*b + a,
-
-      abs (a - b),
-      abs (2*a - b),
-      abs (2*b - a),
-    ];
-
-    [a, b] = sort(combos, (u, v) => error(u) ~- error(v));
-  }
-  return abs a;
-}
-
-riff tune3(a, b, c, numIter = 1, weights = niente) {
-  "Find a combination of three vals that is closer to just intonation.";
-
-  let error = (v => errorTE(v, weights, true));
-  if (isFunction(weights)) {
-    error = weights;
-  }
-
-  numIter = real(numIter) + 1r;
-  while (--numIter) {
-    const combos = [
-      (* Corners *)
-      a,
-      b,
-      c,
-
-      (* Edge midpoints *)
-      a + b,
-      a + c,
-      b + c,
-
-      (* Midpoint *)
-      a + b + c,
-
-      (* Corner extensions *)
-      2 * a + b + c,
-      2 * b + a + c,
-      2 * c + a + b,
-
-      (* Edge extensions *)
-      2 * (a + b) + c,
-      2 * (a + c) + b,
-      2 * (b + c) + a,
-    ];
-
-    [a, b, c] = sort(combos, (u, v) => error(u) ~- error(v));
-  }
-  return abs a;
-}
-
-riff tune4(a, b, c, d, numIter = 1, weights = niente) {
-  "Find a combination of four vals that is closer to just intonation.";
-
-  let error = (v => errorTE(v, weights, true));
-  if (isFunction(weights)) {
-    error = weights;
-  }
-
-  numIter = real(numIter) + 1r;
-  while (--numIter) {
-    const combos = [
-      (* Corners *)
-      a,
-      b,
-      c,
-      d,
-
-      (* Edge midpoints *)
-      a + b,
-      a + c,
-      a + d,
-      b + c,
-      b + d,
-      c + d,
-
-      (* Face midpoints *)
-      a + b + c,
-      a + b + d,
-      a + c + d,
-      b + c + d,
-
-      (* Midpoint *)
-      a + b + c + d,
-
-      (* Corner extensions *)
-      2 * a + b + c + d,
-      2 * b + a + c + d,
-      2 * c + a + b + d,
-      2 * d + a + b + c,
-
-      (* Edge extensions *)
-      2 * (a + b) + c + d,
-      2 * (a + c) + b + d,
-      2 * (a + d) + b + c,
-      2 * (b + c) + a + b,
-      2 * (b + d) + a + c,
-      2 * (c + d) + a + b,
-
-      (* Face extensions *)
-      2 * (a + b + c) + d,
-      2 * (a + b + d) + c,
-      2 * (a + c + d) + b,
-      2 * (b + c + d) + a,
-    ];
-
-    [a, b, c, d] = sort(combos, (u, v) => error(u) ~- error(v));
-  }
-  return abs a;
-}
-
 riff supportingGPVs(initialVal, commas, count = 5, weights = niente, maxIter = 1000) {
   "Obtain generalized patent vals in the same sequence as the initial val that make the given commas vanish.";
   if (not isArray(commas)) {