Implement methods for obtaining TE optimal tuning maps

ref #365

documentation/BUILTIN.md

@@ -384,6 +384,9 @@ Return the higher prime tail of an interval starting from the given index. Prime
 ### tan(*x*)
 Calculate tan x.
 
+### TE(*valsOrCommas*, *weights*)
+Calculate Tenney-Euclid optimal PrimeMapping by combining the given vals or tempering out the given commas. Weights are applied multiplicatively on top of Tenney weights. Use a single large value for CTE. For vals the weights apply to the subgroup basis. A minimal prime subgroup is inferred from the commas, but the weights are for the primes in order if given.
+
 ### templateArg(*index*)
 Access the nth template argument when using the `sw` tag inside JavaScript.
 

documentation/advanced-dsl.md

@@ -513,6 +513,116 @@ Intrinsic behavior may be evoked explicitly by simply calling a value e.g. `3/2(
 Some expressions like `440Hz` or `440 Hz` appear similar to intrinsic calls and would correspond to `440 × (1 Hz)` but `600.0 Hz` is actually `600.0e × (1 Hz)`.
 It's legal to declare `let Hz = 'whatever'`, but the grammar prevents the `Hz` variable from invoking intrinsic behavior of integer literals from the right.
 
+## Advanced tempering
+While vals are technically powerful enough to describe any pure-octaves tuning to sufficient precision, it can be challenging to discover them from first principles.
+
+### Tempering out a comma
+[Tenney-Euclidean tuning](https://en.xen.wiki/w/Tenney-Euclidean_tuning) is a commonly used scheme for adjusting intervals in such a way that some of them coincide. This can help tame the complexity of pure just intonation. It is desirable to make these tiny adjustments with minimal damage.
+
+To make two intervals such as `9/8` and `10/9` coincide we *temper out* their geometric difference `81/80` using the `TE` helper.
+```ocaml
+(* Generate Pythagorean major scale *)
+rank2(3/2, 5, 1)
+
+(* Temper out the syntonic comma 81/80 *)
+TE(9/8 ÷ 10/9)
+```
+
+Now the `81/64` major thirds sounds more like `5/4` and the major chord resembles `4:5:6`. Compared to just intonation the scale maintains the nice property that there are only two step sizes.
+
+The downside is that the scale is now expressed in terms of real cents and all structural information such as the stack of pure fifths has been lost.
+
+#### Maintaining pure octaves
+TE optimal tunings damage all primes including the octave which is often undesirable.
+
+To destretch the octave after the fact use the stretching operator `~^` alongside the size-comparison operator `~/_`. The expression `2 ~/_ $[-1]` tells you how much wider the octave is compared to the last interval in the scale (the current equave).
+
+```ocaml
+"POTE meantone[7]"
+rank2(3/2, 5, 1)
+TE(81/80)
+
+(* Destretch octaves *)
+£ ~^ (2 ~/_ £[-1])
+```
+
+##### CTE
+[Constrained tunings](https://en.xen.wiki/w/Constrained_tuning) have more finesse when it comes to maintaining the size of important intervals like the octave. While SonicWeave doesn't support CTE exactly, you can get close enough by assigning a large weight to the octave.
+
+```ocaml
+"Near-CTE meantone[7]"
+rank2(3/2, 5, 1)
+(* Temper out the syntonic comma while making octaves 1000 times more important than anything else. *)
+TE(81/80, 1000)
+```
+
+#### Tempering out multiple commas
+Multiple commas can be tempered out by passing in an array to `TE`. Multiple primes can be weighted by passing in an array as the second argument.
+
+```ocaml
+"Unimarv[19] with a focus on near-pure octaves and undecimal harmony"
+
+(* Unimarv[19] 5-limit transversal *)
+25/24
+16/15
+9/8
+75/64
+6/5
+5/4
+32/25
+4/3
+45/32
+64/45
+3/2
+25/16
+8/5
+5/3
+128/75
+16/9
+15/8
+48/25
+2/1
+
+(* Temper out the marvel comma and the keenanisma *)
+TE(
+  [225/224, 385/384],
+  (* Give 1000 times more weight to octaves and 10 times more weight to prime 11. *)
+  [1000, 1, 1, 1, 10],
+)
+```
+
+### Combining multiple vals
+Linear combinations of vals share the same temperament. In the 11-limit unimarv can also be expressed as the 19 & 22 temperament. In SonicWeave we must be explicit about the subgroup `@2.3.5.7.11` or `@.11` for short and pass the vals as arguments to `TE`.
+
+```ocaml
+"Marveldene but using unimarv for no reason"
+
+(* Generate duodene *)
+eulerGenus(675)
+
+(* TE unimarv 19 & 22 *)
+TE([[email protected], [email protected]])
+```
+
+It is strongly advised to use an explicit subgroup when combining vals in `TE`. Depending on the runtime the ambient subgroup might contain tens of primes, most of little interest, only wasting compute and hurting the low prime accuracy of the tuning.
+
+#### Subgroup weights
+Because vals always come with a subgroup basis the weights are associated with it instead of the actual primes when combining vals.
+```ocaml
+"Barbados[5]"
+15/13
+4/3
+3/2
+26/15
+2/1
+(** Importance weights:
+ * 2/1: 200%
+ * 3/1: 100%
+ * 13/5: 300%
+ *)
+TE([[email protected]/5, [email protected]/5], [2, 1, 3])
+```
+
 ### Obscure types
 | Type              | Literal    | Meaning                                                    |
 | ----------------- | ---------- | ---------------------------------------------------------- |

examples/jubilismic-gs.sw

@@ -1,4 +1,5 @@
 
 "Generator sequence 4:5:6 around the half-octave (Jubilismic TE)"
 gs(geodiff(4:5:6), 10, 1\2, 2);
-PrimeMapping(1199.35339, 1901.955, 2779.34146, 3379.01816);
+(* Temper out the jubilisma 50/49 using a Tenney-Euclid optimal tuning. *)
+TE([50/49]);

examples/marveldene.sw

@@ -1,3 +1,4 @@
 "The Marveldene, a 2.3.5.7 classic. Septimal Aeolian minor oriented."
 eulerGenus(675, 15)
-PrimeMapping(1200.5978, 1901.3543, 2785.0245, 3369.7682)
+(* Temper out 225/224 by combining 3 linearly independent vals that support the marvel temperament. *)
+TE([12@.7, 31@.7, 41@.7])

examples/semitonismic-gs.sw

@@ -1,12 +1,7 @@
-"(A)GS(7/6, 8/7) over 2^1/2 (Semitonismic TE)"
+"(A)GS(7/6, 8/7) over 2^1/2 (Semitonismic POTE)"
 
 csgs([7/6, 8/7], 3, 2^1/2, 2)
-PrimeMapping(
-  1200.29451,
-  1902.25094,
-  niente,
-  niente,
-  niente,
-  niente,
-  4902.98721
-)
+(* Temper out the semitonisma. TE optimal. *)
+TE([289/288])
+(* De-stretch the octave. (Pure Octaves TE) *)
+£ ~^ (2 ~/_ £[-1])

package.json

@@ -68,7 +68,7 @@
   },
   "dependencies": {
     "moment-of-symmetry": "^0.8.2",
-    "xen-dev-utils": "^0.10.0"
+    "xen-dev-utils": "^0.10.2"
   },
   "engines": {
     "node": ">=12.0.0"

src/__tests__/interval.spec.ts

@@ -3,6 +3,7 @@ import {TimeMonzo, TimeReal} from '../monzo';
 import {Interval, ValBasis, intervalValueAs} from '../interval';
 import {FractionLiteral, NedjiLiteral} from '../expression';
 import {sw} from '../parser';
+import {dot} from 'xen-dev-utils';
 
 describe('Idempontent formatting', () => {
   it('has stable ratios (common factor)', () => {
@@ -221,4 +222,17 @@ describe('(Val) subgroup basis', () => {
     expect(basis.ortho[0].dot(basis.ortho[2]).n).toBe(0);
     expect(basis.ortho[1].dot(basis.ortho[2]).n).toBe(0);
   });
+
+  it('can fix subgroup maps to the standard basis', () => {
+    const basis = new ValBasis([
+      TimeMonzo.fromFraction('3/2'),
+      TimeMonzo.fromFraction('10/9'),
+      TimeMonzo.fromFraction('7/5'),
+    ]);
+    const map = [700, 200, 600];
+    const fixed = basis.standardFix(map);
+    expect(dot(basis.value[0].toIntegerMonzo(), fixed)).toBeCloseTo(700);
+    expect(dot(basis.value[1].toIntegerMonzo(), fixed)).toBeCloseTo(200);
+    expect(dot(basis.value[2].toIntegerMonzo(), fixed)).toBeCloseTo(600);
+  });
 });

src/__tests__/temper.spec.ts

@@ -0,0 +1,90 @@
+import {describe, expect, it} from 'vitest';
+import {TuningMap, combineTuningMaps, vanishCommas} from '../temper';
+import {LOG_PRIMES, applyWeights, dot, unapplyWeights} from 'xen-dev-utils';
+
+describe('Val 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]);
+    const meantone = combineTuningMaps(
+      [1, 1, 1],
+      [unapplyWeights(p12, LOG_PRIMES), unapplyWeights(p19, LOG_PRIMES)]
+    );
+    const map = applyWeights(meantone, LOG_PRIMES);
+    expect(Math.abs(dot(map, [-4, 4, -1]))).toBeCloseTo(0, 9);
+    expect(dot(map, [1, 0, 0])).toBeCloseTo(LOG_PRIMES[0], 2);
+    expect(dot(map, [-1, 1, 0])).toBeCloseTo(LOG_PRIMES[1] - LOG_PRIMES[0], 2);
+    expect(dot(map, [0, 0, 1])).toBeCloseTo(LOG_PRIMES[2], 2);
+    expect((((map[1] - map[0]) / map[0]) * 1200).toFixed(3)).toBe('696.239');
+    expect(map.map(c => ((c / LOG_PRIMES[0]) * 1200).toFixed(3))).toEqual([
+      '1201.397',
+      '1898.446',
+      '2788.196',
+    ]);
+  });
+
+  it('almost computes CTE meantone', () => {
+    const p5 = [5, 8, 12];
+    const p7 = [7, 11, 16];
+    const jip = LOG_PRIMES.slice(0, 3);
+    const weights = LOG_PRIMES.slice(0, 3);
+    weights[0] *= 0.005;
+    const meantone = combineTuningMaps(unapplyWeights(jip, weights), [
+      unapplyWeights(p5, weights),
+      unapplyWeights(p7, weights),
+    ]);
+    const map = applyWeights(meantone, weights);
+    expect(Math.abs(dot(map, [-4, 4, -1]))).toBeCloseTo(0, 7);
+    expect(dot(map, [1, 0, 0])).toBeCloseTo(LOG_PRIMES[0], 7);
+    expect(dot(map, [-1, 1, 0])).toBeCloseTo(LOG_PRIMES[1] - LOG_PRIMES[0], 2);
+    expect(dot(map, [0, 0, 1])).toBeCloseTo(LOG_PRIMES[2], 2);
+    expect((((map[1] - map[0]) / map[0]) * 1200).toFixed(4)).toBe('697.2143');
+    expect(map.map(c => ((c / LOG_PRIMES[0]) * 1200).toFixed(4))).toEqual([
+      '1200.0000',
+      '1897.2144',
+      '2788.8572',
+    ]);
+  });
+});
+
+describe('Comma vanisher', () => {
+  it('computes TE meantone', () => {
+    const syntonic = [-4, 4, -1];
+    const meantone = vanishCommas(
+      [1, 1, 1],
+      [applyWeights(syntonic, LOG_PRIMES)]
+    );
+    const map = applyWeights(meantone, LOG_PRIMES);
+    expect(Math.abs(dot(map, [-4, 4, -1]))).toBeCloseTo(0, 9);
+    expect(dot(map, [1, 0, 0])).toBeCloseTo(LOG_PRIMES[0], 2);
+    expect(dot(map, [-1, 1, 0])).toBeCloseTo(LOG_PRIMES[1] - LOG_PRIMES[0], 2);
+    expect(dot(map, [0, 0, 1])).toBeCloseTo(LOG_PRIMES[2], 2);
+    expect((((map[1] - map[0]) / map[0]) * 1200).toFixed(3)).toBe('696.239');
+    expect(map.map(c => ((c / LOG_PRIMES[0]) * 1200).toFixed(3))).toEqual([
+      '1201.397',
+      '1898.446',
+      '2788.196',
+    ]);
+  });
+
+  it('almost computes CTE meantone', () => {
+    const syntonic = [-4, 4, -1];
+    const jip = LOG_PRIMES.slice(0, 3);
+    const weights = LOG_PRIMES.slice(0, 3);
+    weights[0] *= 0.005;
+    const meantone = vanishCommas(unapplyWeights(jip, weights), [
+      applyWeights(syntonic, weights),
+    ]);
+    const map = applyWeights(meantone, weights);
+    expect(Math.abs(dot(map, [-4, 4, -1]))).toBeCloseTo(0, 7);
+    expect(dot(map, [1, 0, 0])).toBeCloseTo(LOG_PRIMES[0], 7);
+    expect(dot(map, [-1, 1, 0])).toBeCloseTo(LOG_PRIMES[1] - LOG_PRIMES[0], 2);
+    expect(dot(map, [0, 0, 1])).toBeCloseTo(LOG_PRIMES[2], 2);
+    expect((((map[1] - map[0]) / map[0]) * 1200).toFixed(4)).toBe('697.2143');
+    expect(map.map(c => ((c / LOG_PRIMES[0]) * 1200).toFixed(4))).toEqual([
+      '1200.0000',
+      '1897.2144',
+      '2788.8572',
+    ]);
+  });
+});

src/interval.ts

@@ -48,12 +48,17 @@ import {
   FractionalMonzo,
   LOG_PRIMES,
   PRIMES,
+  PRIME_CENTS,
+  add,
   applyWeights,
+  dot,
   fractionalLenstraLenstraLovasz,
   lenstraLenstraLovasz,
   primeLimit,
+  scale,
   unapplyWeights,
 } from 'xen-dev-utils';
+import {TuningMap} from './temper';
 
 /**
  * Interval domain. The operator '+' means addition in the linear domain. In the logarithmic domain '+' correspond to multiplication of the underlying values instead.
@@ -1451,9 +1456,7 @@ export class ValBasis {
         /** Fall through */
       }
     }
-    let basis: number[][] = this.value.map(m =>
-      m.primeExponents.map(f => f.valueOf())
-    );
+    let basis: number[][] = this.value.map(m => m.toMonzo());
     if (weighting === 'tenney') {
       basis = basis.map(pe => applyWeights(pe, LOG_PRIMES));
     }
@@ -1525,6 +1528,25 @@ export class ValBasis {
     return true;
   }
 
+  /**
+   * Fix a map in this basis to the stadard basis.
+   * @param map Tuning map of this basis' elements to cents.
+   * @returns Tuning map of primes to cents.
+   */
+  standardFix(map: TuningMap): TuningMap {
+    let result = PRIME_CENTS.slice(0, this.value[0].numberOfComponents);
+    const basis = this.value.map(m => m.toMonzo());
+    const dual = this.dual.map(m => m.toMonzo());
+    for (let i = 0; i < basis.length; ++i) {
+      const cents = dot(basis[i], result);
+      result = add(
+        result,
+        scale(dual[i], (map[i] - cents) / dot(basis[i], dual[i]))
+      );
+    }
+    return result;
+  }
+
   /**
    * Convert this basis to a virtual AST fragment compatible with wart notation.
    * @returns Array of wart basis elements.