Solution de type PAVA pour un simple QP
Laisser $l,u\in\mathbb{R}^n$, et considérez le QP:
$$\min_{l\le x\le u} {(\Delta x)^\top (\Delta x)}$$
où $\Delta x=[x_2-x_1,\,x_3-x_2,\,\dots,\,x_n-x_{n-1}]^\top$.
C'est-à-dire que nous voulons minimiser le changement au carré des éléments de $x$ sujet à $x$ être au-dessus de la limite inférieure $l$ et en dessous de la limite supérieure $u$.
Mon intuition est que c'est assez simple pour qu'il devrait avoir une solution exacte de style "algorithme de violation adjacente groupée (PAVA)". Est-ce correct? Ce problème a-t-il été étudié dans la littérature antérieure?
A part: Bien sûr, comme indiqué ici, le problème peut avoir plusieurs solutions. Je me fiche de ce qui est retourné.
Plus loin: voici le code de type MATLAB pour une procédure de solution inefficace. Je pense qu'il y a un algorithme beaucoup plus efficace!
while true
xo = x;
x( 1 ) = max( l( 1 ), min( u( 1 ), x( 2 ) ) );
for i = 2 : ( n - 1 )
x( i ) = max( l( i ), min( u( i ), 0.5 * ( x( i - 1 ) + x( i + 1 ) ) ) );
end
x( n ) = max( l( n ), min( u( n ), x( n - 1 ) ) );
if all( abs( x - xo ) < 1e-12 )
break
end
end
Réponses
J'ai trouvé le problème intéressant, j'ai donc conçu un algorithme pour le résoudre (pas de style PAVA). Vous pouvez trouver une implémentation Python ici .
Une formulation géométrique
Les conditions d'optimalité sont, pour $1 \lt i \lt n$:
- $x_i - x_{i-1} = x_{i+1} - x_i~$ si $l_i \lt x_i \lt u_i$
- $x_i - x_{i-1} \leq x_{i+1} - x_i~$ si $x_i = u_i$
- $x_i - x_{i-1} \geq x_{i+1} - x_i~$ si $x_i = l_i$
En utilisant ceux-ci, nous pouvons le remplacer par un problème géométrique avec les mêmes solutions optimales. Pour chaque$i$, il y a une ligne verticale avec $x=i$ et $l_i \leq y \leq u_i$. Le but est de trouver le chemin le plus court passant par tous ces segments de gauche à droite. Tout coude de ce chemin correspondra à$x_i = l_i$ ou $x_i = u_i$.
Un exemple:
Et sa solution:
L'algorithme
L'idée de l'algorithme est de trouver où sera le prochain virage et de recommencer à partir de là. Connaissant le virage précédent, nous maintenons la pente minimale et maximale pour la ligne suivante, et qui$i$ est à la limite, jusqu'à ce que nous soyons obligés d'introduire un nouveau virage.
La boucle principale est donc:
def find_next_bend(n, l, u, x, previous_bend):
min_slope = -float("inf")
max_slope = float("inf")
min_slope_ind = -1
max_slope_ind = -1
for i in range(last_bend+1, n):
new_min_slope = (l[i] - x[previous_bend]) / (i-previous_bend)
new_max_slope = (u[i] - x[previous_bend]) / (i-previous_bend)
if new_min_slope > max_slope:
# Bend upwards
x[max_slope_ind] = u[max_slope_ind]
return max_slope_ind
if new_max_slope < min_slope:
# Bend downwards
x[min_slope_ind] = l[min_slope_ind]
return min_slope_ind
if new_max_slope < max_slope:
# Restrict max possible slope
max_slope_ind = i
max_slope = new_max_slope
if new_min_slope > min_slope:
# Restrict min possible slope
min_slope_ind = i
min_slope = new_min_slope
Il y a des majuscules pour le premier et le dernier virage, que je ne décrirai pas ici: elles sont commentées dans le code complet . La preuve de l'algorithme est laissée comme exercice au lecteur :)
Résultats de calcul
L'algorithme s'adapte parfaitement aux données aléatoires et trouve la solution optimale pour $n > 10^6$en secondes. Théoriquement, c'est seulement$O(n^2)$pour autant que je sache, mais a une complexité linéaire dans la pratique. Je suis convaincu qu'il est possible de trouver quelque chose qui$O(n \log n)$, car il présente une certaine similitude avec les algorithmes de coque convexe 2D, mais c'est amusant pour un autre jour;)