Découvrir Scilab/Optimisation d'une fonction

Optimisation d'une fonction

---

Soit une fonction ƒ de ${ℝ}^n\to ℝ$. L'optimisation consiste à trouver le vecteur x (vecteur à n composantes) donnant la valeur minimale de ƒ.

Optimisation non linéaire

L'optimisation non linéaire est le cas général. On utilise pour cela la fonction optim :

[fopt, xopt] = optim(costf, x0)

où

• costf est la « fonction de coût » de ƒ ; c'est une fonction qui renvoit la valeur de la fonction ƒ en x et le gradient de ƒ en x, défini sous la forme
  function [f, g, ind] = costf(x, ind)
  où f désigne ƒ(x), g est le gradient et ind est un index, un entier permettant de modifier le comportement de costf ;
• x0 est une estimation de la solution ;
• xopt est le vecteur trouvé ;
• fopt = ƒ(xopt), valeur estimée du minimum.

On peut indiquer des bornes inférieures et supérieure de x, sous al forme de vecteurs x_inf et x_sup :

[fopt, xopt] = optim(f, x0, "b", xinf, xsup)

La fonction optim() est une méthode de quasi-Newton utilisant les critères de Wolfe.

Par exemple :

Minimum de la fonction pseudo-convexe différentiable $y=x^2+5\sin x$ avec scilab.

function [f, g, ind] = cout(x, ind)
    f = x.*x + 5*sin(x);
    g = 2*x + 5*cos(x);
endfunction

n = 40;
x = linspace(-2*%pi, 2*%pi, n);
y = cout(x);

x0 = 0;

[fopt, xopt] = optim(cout, x0);

plot(x, y);
plot(xopt, fopt, "k+");

Si l'on ne connaît pas de forme analytique de la fonction dérivée, on peut utiliser la fonction numderivative() :

function [y] = fonction(x)
    y = x.*x + 5*sin(x);
endfunction


function [f, g, ind] = cout(x, ind)
    f = fonction(x);
    g = numderivative(fonction, x);
endfunction

n = 40;
x = linspace(-2*%pi, 2*%pi, n);
y = cout(x);

x0 = 0;

[fopt, xopt] = optim(cout, x0);

plot(x, y);
plot(xopt, fopt, "k+");

Optimisation linéaire

La fonction ƒ est une application linéaire ; elle peut donc s'écrire :

$f(𝐱)={𝐜}^{\mathrm{t}}𝐱$

où c est un vecteur de ${ℝ}^n$ et c^t est sa transposée. On peut par ailleurs restreindre le domaine de recherche à un polyèdre convexe décrit par m inéquations :

$\mathrm{A}𝐱⩽𝐛$

où

• A est une matrice m×n ;
• b est un vecteur de m dimensions.

Scilab dispose de la commande karmarkar() qui permet de résoudre le problème avec l'algorithme de Karmarkar.

La syntaxe la plus simple permet de résoudre le problème sur la frontière du polyèdre, donc avec des égalités partout :

xopt = karmarkar(Ae, be, c)

résout

$\{\begin{array}{c}\underset{𝐱}{\inf }{𝐜}^{\mathrm{t}}𝐱\\ {\mathrm{A}}_{\mathrm{e}}𝐱={𝐛}_{𝐞}\end{array}$

Exemple

Si l'on veut minimiser la fonction ƒ(x₁, x₂, x₃) = –x₁ – x₂ avec pour contraintes x₁ – x₂ = 0 et x₁ + x₂ + x₂ = 2, on a :

$\left(\begin{array}{ccc}-1& -1& 0\end{array}\right)\cdot \left(\begin{array}{c}{x}_1\\ x_2\\ x_3\end{array}\right)=0⟹𝐜=\left(\begin{array}{c}-1\\ -1\\ 0\end{array}\right)$

$\left(\begin{array}{ccc}1& -1& 0\\ 1& 1& 1\end{array}\right)\cdot \left(\begin{array}{c}{x}_1\\ x_2\\ x_3\end{array}\right)=\left(\begin{array}{c}0\\ 2\end{array}\right)⟹{\mathrm{A}}_{\mathrm{e}}=\left(\begin{array}{ccc}1& -1& 0\\ 1& 1& 1\end{array}\right)\text{ ; }{𝐛}_{𝐞}=\left(\begin{array}{c}0\\ 2\end{array}\right)$

La solution s'obtient donc par

Aeq = [1, -1, 0 ; 1, 1,1];
beq = [0 ; 2];
c = [-1 ; -1 ; 0];
xopt = karmarkar(Aeq, beq, c)

et le résultat est :

xopt  = 

   0.9999949
   0.9999949
   0.0000102

On peut y ajouter un système d'inéquations :

xopt = karmarkar(Ae, be, c, [], [], [], [], [], Ai, bi)

résout

$\{\begin{array}{c}\underset{𝐱}{\inf }{𝐜}^{\mathrm{t}}𝐱\\ {\mathrm{A}}_{\mathrm{e}}𝐱={𝐛}_{𝐞}\\ {\mathrm{A}}_{\mathrm{i}}𝐱⩽{𝐛}_{𝐢}\end{array}$

et si l'on n'a que des inéquations :

xopt = karmarkar([], [], c, [], [], [], [], [], Ai, bi)

résout

$\{\begin{array}{c}\underset{𝐱}{\inf }{𝐜}^{\mathrm{t}}𝐱\\ {\mathrm{A}}_{\mathrm{i}}𝐱⩽{𝐛}_{𝐢}\end{array}$

On peut indiquer un vecteur de départ x₀ :

xopt = karmarkar(A, b, c, x0)
xopt = karmarkar(A, b, c, x0, [], [], [], [], Ai, bi)
xopt = karmarkar([], [], c, x0, [], [], [], [], Ai, bi)

et l'on peut demander de calculer la valeur de ƒ(x_opt) :

[xopt,fopt] = karmarkar(…)

Voir aussi

• Modèle:Lien web

Optimisation quadratique

Soit ƒ une fonction quadratique de n variables (x_i)_{1 ≤ i ≤ n}, c'est-à-dire une combinaison linéaire de x_ix_j. Cette fonction peut se décrire par deux matrices, Q, définie positive de dimension n×n, et p, matrice colonne de dimension n :

si l'on appelle x la matrice colonne [x₁ ; x₂ ; … ; x_n]

ƒ(x₁, …, x_n) = ½^txQx + ^tpx

où ^tM est la transposée de la matrice M.

L'optimisation quadratique consiste à trouver le vecteur x donnant la valeur minimum de ƒ, en imposant de plus que la solution se trouve dans un espace convexe, ce qui se traduit par m contraintes linéaires : m_e conditions d'égalité

∑_{1 ≤ j ≤ n}C_ijx_j = b_i, 1 ≤ i ≤ m_e

et m - m_e conditions d'inégalité

∑_{1 ≤ j ≤ n}D_ijx_j ≤ d_i, 1 ≤ i ≤ m - m_e

soit sous forme matricielle

C₁x = b₁

C₂x ≤ b₂

où

• C₁ est une matrice m_e×n ;
• b₁ est une matrice colonne de dimension m_e ;
• C₂ est une matrice (m - m_e)×n ;
• b₂ est une matrice colonne de dimension (m - m_e).

On rassemble les matrices :

• C = [C₁ ; C₂] ;
• b = [b₁ ; b₂].

Pour résoudre un tel système, Scilab propose deux commandes.

La commande qpsolve, sous la forme :

[x] = qpsolve(Q, p, C, b, ci, cs, me)

utilise la fonction Fortran qpgen1.f (également appelée QP.solve.f), developpée par Berwin A. Turlach selon l'algorithme de Goldfarb/Idnani.

Les paramètres c_i et c_s sont des vecteurs colonne de dimension n contenant les limites inférieures et supérieures aux valeurs des x_i

c_i ≤ x ≤ c_s.

Si l'on n'impose pas de limite, on utilise des matrices vides [].

La commande qld, sous la forme :

[x, lagr] = qld(Q, p, C, b, ci, cs, me)

utilise la méthode de Powell modifiée par Schittkowski.

La variable lagr est un vecteur de multiplicateurs lagrangiens.

Notes

---

Résolution d'équations < ↑ > Matrices creuses