Constrained tuning: Difference between revisions

From Xenharmonic Wiki
Jump to navigation Jump to search
Expansion
Update the script and style of the tuning maps
Line 44: Line 44:
Otherwise, as a standard optimization problem, numerous algorithms exist to solve it, such as [[Wikipedia: Sequential quadratic programming|sequential quadratic programming]], to name one.  
Otherwise, as a standard optimization problem, numerous algorithms exist to solve it, such as [[Wikipedia: Sequential quadratic programming|sequential quadratic programming]], to name one.  


The following [https://www.python.org Python] code is an excerpt from [[Flora Canou]]'s [https://github.com/FloraCanou/te_temperament_measures/blob/54c9ec58acf0ad5adc7c1d96f2deaf1d2a503702/tuning_optimizer.py tuning optimizer]. Note: it depends on [https://scipy.org/ Scipy].  
The following [https://www.python.org Python] code is an excerpt from [[Flora Canou]]'s [https://github.com/FloraCanou/te_temperament_measures/blob/545cc301418253cf46785c6c59a0cc6754986b08/tuning_optimizer.py tuning optimizer]. Note: it depends on [https://scipy.org/ Scipy].  


{{Databox| Code |
{{Databox| Code |
<syntaxhighlight lang="python">
<syntaxhighlight lang="python">
# © 2020-2021 Flora Canou | Version 0.8
# © 2020-2021 Flora Canou | Version 0.10
# This work is licensed under the GNU General Public License version 3.
# This work is licensed under the GNU General Public License version 3.


Line 58: Line 58:
SCALAR = 1200 #could be in octave, but for precision reason
SCALAR = 1200 #could be in octave, but for precision reason


def weighted (matrix, subgroup, type = "tenney"):
def weighted (matrix, subgroup, wtype = "tenney"):
     if not type in {"tenney", "frobenius"}:
     if not wtype in {"tenney", "frobenius", "partch"}:
         type = "tenney"
         wtype = "tenney"
        raise Warning ("unknown weighter type, using default (\"tenney\")")


     if type == "tenney":
     if wtype == "tenney":
         weighter = np.diag (1/np.log2 (subgroup))
         weighter = np.diag (1/np.log2 (subgroup))
     elif type == "frobenius":
     elif wtype == "frobenius":
         weighter = np.eye (len (subgroup))
         weighter = np.eye (len (subgroup))
    elif wtype == "partch":
        weighter = np.diag (np.log2 (subgroup))
     return matrix @ weighter
     return matrix @ weighter


Line 71: Line 74:
     return linalg.norm (gen @ map - jip, ord = order)
     return linalg.norm (gen @ map - jip, ord = order)


def optimizer_main (map, subgroup = [], order = 2, weighter = "tenney", cons_monzo_list = np.array ([]), stretch_monzo = np.array ([]), show = True):
def optimizer_main (map, subgroup = None, wtype = "tenney", order = 2, cons_monzo_list = None, stretch_monzo = None, show = True):
     if len (subgroup) == 0:
     if subgroup is None:
         subgroup = PRIME_LIST[:map.shape[1]]
         subgroup = PRIME_LIST[:map.shape[1]]


     jip = np.log2 (subgroup)*SCALAR
     jip = np.log2 (subgroup)*SCALAR
     map_w = weighted (map, subgroup, type = weighter)
     map_w = weighted (map, subgroup, wtype = wtype)
     jip_w = weighted (jip, subgroup, type = weighter)
     jip_w = weighted (jip, subgroup, wtype = wtype)
     if order == 2 and not cons_monzo_list.size: #te with no constraints, simply use lstsq for better performance
     if order == 2 and cons_monzo_list is None: #te with no constraints, simply use lstsq for better performance
         res = linalg.lstsq (map_w.T, jip_w)
         res = linalg.lstsq (map_w.T, jip_w)
         gen = res[0]
         gen = res[0]
Line 84: Line 87:
     else:
     else:
         gen0 = [SCALAR]*map.shape[0] #initial guess
         gen0 = [SCALAR]*map.shape[0] #initial guess
         cons = {'type': 'eq', 'fun': lambda gen: (gen @ map - jip) @ cons_monzo_list} if cons_monzo_list.size else ()
         cons = () if cons_monzo_list is None else {'type': 'eq', 'fun': lambda gen: (gen @ map - jip) @ cons_monzo_list}
         res = optimize.minimize (error, gen0, args = (map_w, jip_w, order), method = "SLSQP", constraints = cons)
         res = optimize.minimize (error, gen0, args = (map_w, jip_w, order), method = "SLSQP", constraints = cons)
         print (res.message)
         print (res.message)
Line 90: Line 93:
             gen = res.x
             gen = res.x


     if stretch_monzo.size:
     if not stretch_monzo is None:
         gen *= (jip @ stretch_monzo)/(gen @ map @ stretch_monzo)
         gen *= (jip @ stretch_monzo)/(gen @ map @ stretch_monzo)


Line 97: Line 100:


     return gen
     return gen
optimiser_main = optimizer_main


</syntaxhighlight>
</syntaxhighlight>
Line 106: Line 111:
The pure-octave CTE tuning can be very different from [[POTE tuning]]. Take 7-limit meantone as an example, the POTE [[tuning map]]:  
The pure-octave CTE tuning can be very different from [[POTE tuning]]. Take 7-limit meantone as an example, the POTE [[tuning map]]:  


{{val| 1200.000 1896.495 2785.980 3364.949 }}
<math>\langle \begin{matrix} 1200.000 & 1896.495 & 2785.980 & 3364.949 \end{matrix} ]</math>


This is a little bit flatter than [[quarter-comma meantone]], with all the primes tuned flat.  
This is a little bit flatter than [[quarter-comma meantone]], with all the primes tuned flat.  
Line 112: Line 117:
The pure-octave CTE tuning map:  
The pure-octave CTE tuning map:  


{{val| 1200.000 1896.952 2787.809 3369.521}}
<math>\langle \begin{matrix} 1200.000 & 1896.952 & 2787.809 & 3369.521 \end{matrix} ]</math>


This is a little bit sharper than quarter-comma meantone, with prime 3 tuned flat and 5 and 7 sharp.  
This is a little bit sharper than quarter-comma meantone, with prime 3 tuned flat and 5 and 7 sharp.  

Revision as of 11:14, 1 February 2022

The CTE tuning (constrained Tenney-Euclidean tuning) is TE tuning under the constraints of some purely tuned intervals (i.e. eigenmonzos). While the TE tuning can be viewed as a least squares problem, the CTE tuning can be viewed as an equality-constrained least squares problem. For a rank-r temperament, specifying m eigenmonzos will yield r - m degrees of freedom to be optimized.

The most significant form of CTE tuning is pure-octave constrained. For higher-rank temperaments, it may make sense to add multiple constraints, such as the pure-{2, 3} CTE tuning.

Definition

Given a temperament mapping A and the JIP J0, denote the Tenney-weighted temperament mapping by V = AW, and the Tenney-weighted JIP by J = J0W. If the tuning is contrained by the eigenmonzo list B, the CTE tuning is equivalent to the following optimization problem:

Minimize

[math]\displaystyle{ \lVert GV - J \rVert }[/math]

subject to

[math]\displaystyle{ (GA - J_0)B = O }[/math]

where G is the generator list, and O the zero matrix.

The problem is feasible if

  1. rank (B) ≤ rank (A), and
  2. Each column in B and N (A) are linearly independent.

Computation

The tuning can be solved in the method of Lagrange multiplier. The solution is given by

[math]\displaystyle{ \begin{bmatrix} G^{\mathsf T} \\ \Lambda^{\mathsf T} \end{bmatrix} = \begin{bmatrix} VV^{\mathsf T} & AB \\ (AB)^{\mathsf T} & O \end{bmatrix}^{-1} \begin{bmatrix} VJ^{\mathsf T}\\ (J_0 B)^{\mathsf T} \end{bmatrix} }[/math]

which is almost an analytical solution. Notice we introduced the vector of lagrange multipliers Λ, with length equal to the number of constraints. The lagrange multipliers have no concrete meaning for the resulting tuning, so they can be discarded.

Otherwise, as a standard optimization problem, numerous algorithms exist to solve it, such as sequential quadratic programming, to name one.

The following Python code is an excerpt from Flora Canou's tuning optimizer. Note: it depends on Scipy.

Code
# © 2020-2021 Flora Canou | Version 0.10
# This work is licensed under the GNU General Public License version 3.

import numpy as np
from scipy import optimize, linalg
np.set_printoptions (suppress = True, linewidth = 256)

PRIME_LIST = [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61]
SCALAR = 1200 #could be in octave, but for precision reason

def weighted (matrix, subgroup, wtype = "tenney"):
    if not wtype in {"tenney", "frobenius", "partch"}:
        wtype = "tenney"
        raise Warning ("unknown weighter type, using default (\"tenney\")")

    if wtype == "tenney":
        weighter = np.diag (1/np.log2 (subgroup))
    elif wtype == "frobenius":
        weighter = np.eye (len (subgroup))
    elif wtype == "partch":
        weighter = np.diag (np.log2 (subgroup))
    return matrix @ weighter

def error (gen, map, jip, order = 2):
    return linalg.norm (gen @ map - jip, ord = order)

def optimizer_main (map, subgroup = None, wtype = "tenney", order = 2, cons_monzo_list = None, stretch_monzo = None, show = True):
    if subgroup is None:
        subgroup = PRIME_LIST[:map.shape[1]]

    jip = np.log2 (subgroup)*SCALAR
    map_w = weighted (map, subgroup, wtype = wtype)
    jip_w = weighted (jip, subgroup, wtype = wtype)
    if order == 2 and cons_monzo_list is None: #te with no constraints, simply use lstsq for better performance
        res = linalg.lstsq (map_w.T, jip_w)
        gen = res[0]
        print ("L2 tuning without constraints, solved using lstsq. ")
    else:
        gen0 = [SCALAR]*map.shape[0] #initial guess
        cons = () if cons_monzo_list is None else {'type': 'eq', 'fun': lambda gen: (gen @ map - jip) @ cons_monzo_list}
        res = optimize.minimize (error, gen0, args = (map_w, jip_w, order), method = "SLSQP", constraints = cons)
        print (res.message)
        if res.success:
            gen = res.x

    if not stretch_monzo is None:
        gen *= (jip @ stretch_monzo)/(gen @ map @ stretch_monzo)

    if show:
        print (f"Generators: {gen} (¢)", f"Tuning map: {gen @ map} (¢)", sep = "\n")

    return gen

optimiser_main = optimizer_main

Constraints can be added from the parameter cons_monzo_list of the optimizer_main function.

Versus POTE tuning

The pure-octave CTE tuning can be very different from POTE tuning. Take 7-limit meantone as an example, the POTE tuning map:

[math]\displaystyle{ \langle \begin{matrix} 1200.000 & 1896.495 & 2785.980 & 3364.949 \end{matrix} ] }[/math]

This is a little bit flatter than quarter-comma meantone, with all the primes tuned flat.

The pure-octave CTE tuning map:

[math]\displaystyle{ \langle \begin{matrix} 1200.000 & 1896.952 & 2787.809 & 3369.521 \end{matrix} ] }[/math]

This is a little bit sharper than quarter-comma meantone, with prime 3 tuned flat and 5 and 7 sharp.

It can be speculated that POTE tends to result in biased tunings whereas CTE less so.

Special constraint

The special eigenmonzo Wj removes the weighted tuning bias, where j is the all-ones monzo:

[math]\displaystyle{ W \vec j = [ \begin{matrix} 1 & 1/\log_2 (3) & \ldots & 1/\log_2 (p) \end{matrix} \rangle }[/math]

The monzo introduces weight to the constraint:

[math]\displaystyle{ (GV - J)\vec j = 0 }[/math]

It can be regarded as a distinct optimum, the TOCTE tuning (Tenney ones constrained Tenney-Euclidean tuning).

For equal temperaments

Since the one degree of freedom of equal temperaments is determined by the constraint, optimization is not involved. It is thus reduced to TOC tuning (Tenney ones constrained tuning). This constrained tuning demonstrates interesting properties.

The step size g can be found by

[math]\displaystyle{ g = 1/\operatorname {mean} (V) }[/math]

The edo number n can be found by

[math]\displaystyle{ n = 1/g = \operatorname {mean} (V) }[/math]

Unlike TE or TOP, the optimal edo number space in TOC is linear with respect to A. That is, if A = A1 + A2, then

[math]\displaystyle{ n = \frac {AW \vec j}{J \vec j} \\ = \frac {(A_1 + A_2) W \vec j}{J \vec j} \\ = \frac {A_1 W \vec j}{J \vec j} + \frac {A_2 W \vec j}{J \vec j} \\ = n_1 + n_2 }[/math]

As a result, the relative error space is also linear with respect to A.

For example, the relative errors of 12ettoc5 (12et in 5-limit TOC) is

[math]\displaystyle{ E_\text {r, 12} = \langle \begin{matrix} -1.55\% & -4.42\% & +10.08\% \end{matrix} ] }[/math]

That of 19ettoc5 is

[math]\displaystyle{ E_\text {r, 19} = \langle \begin{matrix} +4.08\% & -4.97\% & -2.19\% \end{matrix} ] }[/math]

As 31 = 12 + 19, the relative errors of 31ettoc5 is

[math]\displaystyle{ E_\text {r, 31} = \langle \begin{matrix} +2.52\% & -9.38\% & +7.88\% \end{matrix} ] }[/math]