Comment produire des tracés graphiques vectoriels ultra-compacts avec un ajustement de distance orthogonal

1. Présentation

sont un moyen naturel de visualiser des fonctions mathématiques. Avec les graphiques vectoriels, une fonction est approchée par des segments de courbes de Bézier cubiques connectées qui sont rastérisées (c’est-à-dire converties en pixels) lorsqu’elles sont affichées. [1]. La rastérisation étant retardée, les images graphiques vectorielles sont naturellement plus portables que les images basées sur des pixels, car le rendu peut être adapté à son environnement d’affichage. Quel que soit le zoom que vous effectuez sur un tracé graphique vectoriel, vous verrez toujours des segments de lignes nets ; alors que si vous zoomez suffisamment sur une image pixellisée, vous finirez par voir les blocs granuleux qui représentent les pixels.

Bien que les graphiques vectoriels constituent un excellent format de traçage, produire de bons graphiques vectoriels peut s’avérer difficile. Considérons, par exemple, comment nous pourrions adopter cet exemple de tracé de matplotlib [2] pour la fonction f

import matplotlib.pyplot as plt
import numpy as np
def f
  return np.exp(-t) * np.cos(2*np.pi*t)
t = np.arange(0.0, 5.0, 0.02)
plt.plot(t, f
plt.savefig('fig.svg')

Matplotlib se rapproche du tracé lisse avec 250 segments de ligne connectés par morceaux (affichés avec des couleurs alternées). Bien que l’intrigue semble bonne, elle est beaucoup plus vaste qu’elle ne devrait l’être. Si, à la place, nous produisions une image tramée comme un PNG, les détails de la façon dont la courbe est construite n’auraient pas d’importance ; mais avec un graphique vectoriel, les segments de ligne individuels sont parcourus et conservés dans l’image de sortie. Avec quelques ajustements, nous pouvons cependant améliorer considérablement la taille sans sacrifier la qualité de l’image.

La primitive de base des graphiques vectoriels est la cubique paramétrique, représentée sous la forme d’une courbe de Bézier cubique. Avec les courbes de Bézier cubiques par morceaux, nous disposons de beaucoup plus de boutons que nous pouvons ajuster pour approximer les fonctions que si nous nous limitions à des segments cubiques linéaires par morceaux ou par morceaux (non paramétriques). En utilisant l’algorithme d’ajustement de distance orthogonal (ODF) de cet article, nous pouvons produire un tracé de la fonction qui ne nécessite que 8 segments de Bézier et est visuellement impossible à distinguer du graphique matplotlib. Ci-dessous, je montre la représentation de l’intrigue dans TikZ. Plus tard, je montrerai comment nous pouvons facilement transformer la commande TikZ en une image SVG pour le Web à l’aide de MetaPost.

\draw  (0.00000, 2.50000) ..  controls (0.39088, 2.10882) and (0.61233, -1.05415) ..
       (1.18750, 0.36510) ..  controls (1.60393, 1.31193) and (1.71252, 2.10800) ..
       (2.37500, 0.95127) ..  controls (2.88922, 0.15365) and (3.15477, 0.95162) ..
       (3.56250, 1.11921) ..  controls (3.98299, 1.31658) and (4.26123, 0.78781) ..
       (4.75000, 0.82415) ..  controls (5.02221, 0.81871) and (5.38203, 1.12388) ..
       (5.93750, 0.99939) ..  controls (5.96426, 1.01981) and (6.36975, 0.82102) ..
       (7.12500, 0.95127) ..  controls (8.08129, 1.04760) and (7.44859, 0.87986) ..
       (9.50000, 0.96171);

L’algorithme s’appuie sur le travail d’Alvin Penner [3]. Étant donné une fonction paramétrique f, elle ajuste d’abord une série de Chebyshev pour se rapprocher de f. Pour les fonctions analytiques, l’interpolation dans les nœuds Chebyshev permet une convergence géométrique rapide [4]. C’est bien mieux que la convergence typique du quatrième ordre que vous obtiendrez avec une interpolation en splines cubiques et bien mieux que la convergence quadratique que vous obtiendrez avec une interpolation linéaire par morceaux. [5]. En interpolant dans les nœuds Chebyshev, nous pouvons réduire le nombre de fois où nous devons évaluer f. À l’aide d’un optimiseur de région de confiance, l’algorithme recherche ensuite une courbe de Bézier cubique qui se rapproche de manière optimale de la fonction cible. Si la distance orthogonale maximale entre la courbe de Bézier ajustée et la série de Chebyshev pour f est inférieure à un seuil cible, alors super, nous avons terminé ; sinon, nous divisons le domaine en deux et répétons le processus. Je détaille les étapes ci-dessous.

Algorithme F (adapter un chemin de Bézier à une fonction). Étant donné une fonction paramétrique f

F1. En utilisant l’algorithme développé par Jared Aurentz et Lloyd Trefethen de [7]ajustez les séries Chebyshev f~_x , f~_y à f_x, f_y. Lorsque f est analytique ou différentiable à un degré élevé, l’ajustement est généralement proche de la précision machine, donc à l’avenir, nous supposons que nous pouvons utiliser f~ comme proxy de f et que toute perte de précision sera négligeable.

F2. Utilisation d’un optimiseur de région de confiance [8] et l’algorithme de Penner [3]ajustez une courbe de Bézier g pour minimiser

Plus de détails sur cette étape sont fournis au §3.

F3. Calculer

Si M est inférieur au seuil cible, terminer ; sinon, réglez

et répétez les étapes F1 à F3 pour f_l et f_r jusqu’à ce que le seuil soit atteint. Voir §4 pour plus de détails.

Dans la section suivante, je décris comment adapter des fonctions arbitraires à l’algorithme F à l’aide du package Python bbai (https://github.com/rnburn/bbai).

2. Visite des fonctionnalités

Le code ci-dessous montre la tâche de base consistant à adapter une fonction à une fenêtre spécifiée :

from bbai.graphics import BezierPath
import numpy as np
def f
  return np.exp(-t * t)
path = BezierPath(
     dst_xmin=0, dst_xmax=9.5,
     dst_ymin=0, dst_ymax=2)
path.fit(f, -2, 2)
print(path.tikz_)

sorties

\draw (0.000, 0.000)..controls (2.684, 0.092) and (3.273, 1.952)
..    (4.750, 2.000)..controls (6.229, 1.951) and (6.815, 0.092)
..    (9.500, 0.000);

Par défaut, la bibliothèque mettra à l’échelle le tracé afin que la fonction tienne simplement dans la fenêtre définie par dst_xmin, dst_xmax, dst_ymin et dst_ymax ; mais cela peut être modifié en spécifiant également une fenêtre source avec src_xmin, src_xmax, src_ymin et src_ymax. L’algorithme utilise 1,0 × 10^-2 comme seuil de distance orthogonale maximale par défaut qui, pour la perspective, est un centième de la largeur de ligne par défaut de TikZ.

Nous pouvons également adapter des fonctions paramétriques en fournissant à la fois une fonction x et une fonction ay.

from bbai.graphics import BezierPath
import numpy as np
R, r, d = 5, 3, 5
def fx
  return (R-r)*np.cos
def fy
  return (R-r)*np.sin
path = BezierPath(
    dst_xmin=0, dst_xmax=2.5,
    dst_ymin=0, dst_ymax=2.5)
path.fit(fx, fy, 0, 2*np.pi*r*d/R)
print(path.tikz_)

sorties

\draw (2.500, 1.250)..controls (2.533, 1.060) and (1.400, 0.745)
..    (0.850, 0.578)..controls (-0.049, 0.321) and (-0.145, 0.403)
..    (0.148, 0.876)..controls (0.197, 0.986) and (1.203, 2.351)
..    (1.405, 2.450)..controls (1.594, 2.564) and (1.722, 2.598)
..    (1.716, 1.250)..controls (1.722, -0.033) and (1.609, -0.085)
..    (1.405, 0.049)..controls (1.203, 0.149) and (0.197, 1.514)
..    (0.149, 1.624)..controls (-0.137, 2.086) and (-0.067, 2.185)
..    (0.851, 1.921)..controls (1.203, 1.809) and (2.534, 1.455)
..    (2.5000, 1.2500);

3. L’algorithme d’ajustement

Cette section décrit plus en détail l’étape F2 où nous ajustons une courbe de Bézier, g, pour approximer f~. L’algorithme est similaire à l’algorithme de Penner de [3] mais avec quelques modifications. Une courbe de Bézier cubique est paramétrée par quatre points. Dans l’algorithme d’ajustement, nous faisons correspondre les extrémités de f~, ce qui nous donne deux points ou quatre paramètres que nous pouvons ajuster. Notons θ les paramètres réglables ; soit B_θ la courbe de Bézier cubique avec des paramètres θ qui correspondent aux extrémités f~(a) et f~(b) ; et définir

Le calcul de base nous dit que s_t doit soit être un point final, soit satisfaire l’équation

Comme B_θ est une cubique paramétrique, les valeurs de s qui satisfont à l’équation sont les racines d’une quintique qui peut être facilement résolue en trouvant les valeurs propres de sa matrice collègue associée. [4].

Nous pouvons utiliser une quadrature de Clenshaw Curtis pour approximer la fonction objectif

Le gradient et la Hesse de h peuvent être calculés à l’aide du théorème des fonctions implicites. Voir [3] pour les équations. Nous pouvons maintenant mettre en place l’algorithme d’ajustement.

Algorithme B (adapter une courbe de Bézier à une fonction). Étant donné une fonction paramétrique f~

B1. Si possible, sélectionnez θ_0 pour que B_{θ_0} corresponde à la courbure de f~ en a et b ; sinon, sélectionnez θ_0 pour que B_{θ_0} soit le segment de droite qui passe par f(a) et f(b).

B2. À partir de θ_0 et en utilisant h, ∇h et ∇^2 h, appliquez un optimiseur de région de confiance [8] jusqu’à ce que nous soyons suffisamment proches d’un optimum ou que nous ayons dépassé un nombre maximum d’étapes prédéterminé

Le principal avantage de l’utilisation d’un optimiseur de région de confiance pour l’étape B2 est qu’il ne restera pas bloqué aux points de selle. En utilisant des informations de second ordre combinées à une « région de confiance » adaptative, un optimiseur de région de confiance peut toujours progresser même lorsque ∇h est proche de zéro et ∇^2h est indéfini.

4. L’algorithme de distance maximale

Cette section décrit l’étape F3 plus en détail où, étant donné une courbe de Bézier cubique, g, nous calculons la distance orthogonale maximale de g à f. Avec s_t défini comme au §3, mettez

et définissons l’algorithme suivant :

Algorithme M (trouver la distance orthogonale maximale). Étant donné une fonction paramétrique f~

M1. Définissez D_max ← 0.

M2. En utilisant la méthode de bissection, partez de l’intervalle [a, b] et trouver un maximum local t_0 de r. Définissez D_max ← max(D_max, r(t_0)).

M3. Mettez s_0 ← s_{t_0} et laissez s_l0′ , s_l0′′ , s_r0′ et s_r0′′ désigner les dérivées gauche et droite de

Notez que les valeurs de gauche et de droite peuvent différer. Définir les fonctions

et

M4. Trouver les plus petits h_l>0 et h_r>0 tels que

Si h_l et h_r appropriés existent, répétez les étapes M2 à M4 pour [a, t_0 −h_l] et [t_0 +h_r,b]. Sinon, renvoyez D_max.

Rappelez-vous de l’algorithme F que f~ est représenté comme une série de Chebyshev, c’est donc une étape facile pour calculer les représentations de la série de Chebyshev pour r~_l et r~_r. Ainsi, l’étape M4 peut être accomplie en appliquant l’algorithme de recherche de racine de Boyd [6].

La figure ci-dessous montre le résultat de l’algorithme de distance maximale en essayant d’ajuster la fonction f

\draw (0.000, 1.800)..controls (4.484, 4.660) and (8.448, -3.442)..(9.500, 1.800);

Nous pouvons voir que même s’il existe plusieurs optima locaux, l’algorithme M identifie correctement l’optimum global.