Correction sur la trace de l'algo
[Projet_Recherche_Operationnelle.git] / rapport / ProjetOptimRO.tex
index e10b20c7f863aa6aaebe21c62c5df305948b0e7e..ccce4d2d1ae6ad6280595a3c715deff9bdb60f71 100644 (file)
@@ -21,6 +21,8 @@
 \usepackage{tocbibind}
 \usepackage{lmodern}
 \usepackage{enumitem}
+\usepackage{algorithm2e}
+\usepackage{algorithmic}
 
 
 %%%%%Marges & en-t\^etes
@@ -233,7 +235,11 @@ Définissons quelques notions supplémentaires de base nécessaires à la suite
  $ A \subset \mathbb{R}^n $ est un fermé $ \iff A = \overline{A} $.
 \end{Rmq}
 \begin{Def}
- Soient une fonction $ f: \mathbb{R}^n \longrightarrow \mathbb{R} $ et $ x^\ast \in \mathbb{R}^n $.
+ Soient $ f : \mathbb{R}^n \longrightarrow \mathbb{R} $ et $ S \subset \mathbb{R}^n $. On définit $ \mathrm{argmin} $ de $ f $ sur $ S $ par :
+ $$ \underset{x \in S}{\mathrm{argmin}} f(x) = \{ x \in \mathbb{R}^n \ | \ x \in S \land \forall y \in S \ f(y) \geq f(x) \} $$
+\end{Def}
+\begin{Def}
+ Soient une fonction $ f : \mathbb{R}^n \longrightarrow \mathbb{R} $ et $ x^\ast \in \mathbb{R}^n $.
  \newline
  On dit que $ f $ est continue en $ x^\ast $ si
  $$ \forall \varepsilon \in \mathbb{R}_{+}^{*} \ \exists \alpha \in \mathbb{R}_{+}^{*} \ \forall x \in \mathbb{R}^n \ \norme{x - x^\ast} \leq \alpha \implies |f(x) - f(x^\ast)| \leq \varepsilon $$
@@ -250,7 +256,7 @@ Définissons quelques notions supplémentaires de base nécessaires à la suite
  cette dérivée.
 \end{Def}
 \begin{Def}
- Soient une fonction $ f: \mathbb{R}^n \longrightarrow \mathbb{R} $
+ Soient une fonction $ f : \mathbb{R}^n \longrightarrow \mathbb{R} $
  et $ x^\ast, h \in \mathbb{R}^n $.
  \newline
  On dit que $ f $ est différentiable en $ x^\ast $ si il existe une application linéraire $ d_{x^\ast}f $ de $ \mathbb{R}^n $ dans $ \mathbb{R} $ telle que
@@ -276,10 +282,10 @@ Définissons quelques notions supplémentaires de base nécessaires à la suite
  \]
 \end{Def}
 \begin{Rmq}
- $ \forall h \in \mathbb{R}^n \ d_{x^\ast}f(h) = \langle \nabla f(x^\ast),h \rangle $
+ $ \forall h \in \mathbb{R}^n \ d_{x^\ast}f(h) = \langle \nabla f(x^\ast),h \rangle = \nabla f(x^\ast)^\top h $
 \end{Rmq}
 \begin{Def}
- Soit $ f: \mathbb{R}^n \longrightarrow \mathbb{R} $ un fonction de classe $ \mathcal{C}^2 $.
+ Soit $ f: \mathbb{R}^n \longrightarrow \mathbb{R} $ une fonction de classe $ \mathcal{C}^2 $.
  On définit la matrice hessienne de $ f $ en $ x^\ast $ par :
  $$ H[f](x^\ast) =
   \begin{pmatrix}
@@ -304,8 +310,8 @@ Définissons quelques notions supplémentaires de base nécessaires à la suite
 
 \subsection{Conditions d'existence d'un extremum}
 
-On peut démontrer que $ \mathcal{C }$ est un ensemble fermé de $ \mathbb{R}^n $ si $ g $ et $ h $ sont continues.
-On peut en déduire que si $ J $ est continue, $ \mathcal{C }$ est un ensemble fermé et borné de $ \mathbb{R}^n $.
+On peut démontrer que $ \mathcal{C}$ est un ensemble fermé de $ \mathbb{R}^n $ si $ g $ et $ h $ sont continues.
+On peut en déduire $ \mathcal{C} $ est un ensemble fermé et borné de $ \mathbb{R}^n $.
 \begin{Th}[Théorème de Weierstrass]
  Soient $ \mathcal{C} \neq \emptyset \subset \mathbb{R}^n $ un fermé borné et $ f : \mathcal{C} \longrightarrow \mathbb{R} $ une fonction continue.
  \newline
@@ -314,7 +320,7 @@ On peut en déduire que si $ J $ est continue, $ \mathcal{C }$ est un ensemble f
  \newline
  De la même façon, il existe un maximum global de $ J $ sur $ \mathcal{C} $.
 \end{Th}
-On en déduit que $ \mathcal{P} $ admet au moins une solution dans le cas où $ J, g ,h $ sont continues \cite{LJK,RON}. L'étude de la convexité de $ J $ sur $ \mathcal{C} $ permet d'explorer l'unicité de la solution \cite{LJK,RON}.
+Si $ J $ est continue, on en déduit que $ \mathcal{P} $ admet au moins une solution dans le cas où $ J, g ,h $ sont continues \cite{LJK,RON}. L'étude de la convexité de $ J $ sur $ \mathcal{C} $ permet d'explorer l'unicité de la solution \cite{LJK,RON}.
 
 \subsection{Conditions de caractérisation d'un extremum}
 
@@ -322,7 +328,7 @@ Dans le cas où $ J, g, h $ sont continûment différentiable et ses dérivées
 \newline
 On peut en déduire que une condition nécessaire et suffisante pour que $ x^\ast \in \mathring{\mathcal{C}} $ soit un des extremums locaux de $ J $ est que $ \nabla J(x^\ast) = 0 $. Mais si $ x^\ast \in \overline{\mathcal{C}}\setminus\mathring{\mathcal{C}} $ (la frontière de $ \mathcal{C} $) alors $ \nabla J(x^\ast) $ n'est pas nécessairement nul. Il sera par conséquent nécessaire de trouver d'autres caratérisations d'un extremum local \cite{FEA,WAL}.
 
-\subsubsection{Conditions de Karuch-Kuhn-Tucker}\label{KKT}
+\subsubsection{Conditions nécessaires de Karuch-Kuhn-Tucker ou \textit{KKT}}\label{KKT}
 
 \begin{Th}
  Soient $ x^\ast \in \mathbb{R}^n $, $ I = \{ 1,\ldots,p \} $ et $ J = \{ 1,\ldots,q \} $.
@@ -336,12 +342,26 @@ On peut en déduire que une condition nécessaire et suffisante pour que $ x^\as
  et
  $$ \forall i \in I \ \exists \mu_i \in \mathbb{R}_{+} \land \forall j \in J \ \exists \lambda_j \in \mathbb{R} \ \nabla J(x^\ast) + \sum_{i \in I}\mu_i{\nabla g_i(x^\ast)} + \sum_{j \in J}\lambda_j{\nabla h_j(x^\ast)} = 0 \land \forall i \in I \ \mu_i \nabla g_i(x^\ast) = 0 $$
  On appelle $ (\mu_i)_{i \in I}$ les multiplicateurs de Kuhn-Tucker et $ (\lambda_j)_{j \in J}$ les multiplicateurs de Lagrange.
+ \newline
+ On nomme également les conditions \textit{KTT} conditions nécessaires d'optimalité du premier ordre.
 \end{Th}
 \begin{proof}
  Elle repose sur le lemme de Farkas \cite{FEA,RON}.
 \end{proof}
 Il est à noter que une condition d'égalité peut se répresenter par deux conditions d'inégalité : $ \forall x \in \mathbb{R}^n \ \forall i \in \{ 1,\ldots,q \} \ h_i(x) = 0 \iff h_i(x) \leq 0 \land h_i(x) \geq 0 $ \cite{FEA}, ce qui peut permettre de réécrire le problème $ \mathcal{P} $ en éliminant les contraintes d'égalités et change la forme des conditions \textit{KKT} à vérifier mais rajoute $ 2q $ conditions d'inégalités et donc $ 2q $ multiplicateurs de Kuhn-Tucker.
+\begin{Def}
+ On appelle un point admissible $ x^\ast \in \mathcal{C} $ un point critique de $ \mathcal{P} $ si il statisfait les conditions \textit{KKT}.
+\end{Def}
+
+\subsubsection{Conditions suffisantes du deuxième ordre}
 
+\begin{Th}
+ Les conditions suffisantes en plus de celles \textit{KKT} pour que $ x^\ast \in \mathcal{C} $ soit un minimum local de $ J $ sont :
+ \begin{enumerate}[label=(\roman*)]
+  \item relâchement complémentaire dual\footnote{La définition de cette notion ne sera pas donnée car elle n'est pas nécessaire pour l'étude de la méthode PQS.} strict en $ x^\ast $.
+  \item $ \forall v \in \mathbb{R}^n \land v \neq 0 \ \langle H_x[L](x^\ast,\lambda,\mu)v,v \rangle > 0 $.
+ \end{enumerate}
+\end{Th}
 
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
@@ -453,14 +473,13 @@ RECHERCHE LINÉAIRE : PAS OPTIMAL. $ s_k $ solution du problème $ \displaystyle
 
 \hrulefill
 \newline
-Illustrées par les méthodes de descente de gradient, aucune de ces deux stratégies ne
-s’est révélée réellement convaincante : si la première peut être “risquée” du point de vue de la convergence, la seconde est souvent loin d’être triviale à mettre en oeuvre (sauf dans le cas quadratique) et généralement inutilement coûteuse : en effet, à quoi bon calculer très précisément un pas optimal dans une direction qui n’est peut-être pas la bonne ? (comme c’est par exemple le cas pour la méthode de plus profonde descente). Les recherches linéaires modernes reposent sur l’idée qu’un pas de descente acceptable est un pas qui fait “suffisamment” décroître la fonction objectif. Reste alors à définir les pas qui sont acceptables et ceux qui ne le sont pas.
+Illustrées par les méthodes de descente de gradient, aucune de ces deux stratégies ne s’est révélée réellement convaincante : si la première peut être “risquée” du point de vue de la convergence, la seconde est souvent loin d’être triviale à mettre en oeuvre (sauf dans le cas quadratique) et généralement inutilement coûteuse : en effet, à quoi bon calculer très précisément un pas optimal dans une direction qui n’est peut-être pas la bonne ? (comme c’est par exemple le cas pour la méthode de plus profonde descente). Les recherches linéaires modernes reposent sur l’idée qu’un pas de descente acceptable est un pas qui fait “suffisamment” décroître la fonction objectif. Reste alors à définir les pas qui sont acceptables et ceux qui ne le sont pas.
 \begin{Def}
  On appelle $ \varphi : s \in \mathbb{R} \longmapsto J(x + sd)$ la fonction mérite associée à $ J $ en $ x $.
 \end{Def}
 \begin{Def}
  Dans le cas où $ J $ est différentiable sur $ \mathcal{C} $, on dit que un algorithme de descente converge ssi
- $$ \lim\limits_{k \rightarrow +\infty} \norme{\nabla J(x_k)} = 0 $$
+ $$ \forall x_0 \in \mathbb{R}^n \lim\limits_{k \rightarrow +\infty} \norme{\nabla J(x_k)} = 0 $$
 \end{Def}
 
 \subsubsection{Principe de démonstration de convergence}
@@ -473,7 +492,19 @@ En sommant ces inégalités pour $ k $ variant de $ 0 $ à $ N - 1 $, on obtient
 $$ \forall N \in \mathbb{N} \ J(x_0) - J(x_N) \geq c \sum_{i=0}^{N-1}\norme{\nabla J(x_i)}^2 $$
 Si $ J $ est bornée inférieurement, alors nécessairement $ J(x_0 ) - J(x_N) $ est majorée et donc la somme partielle est majorée, et donc la série $ (\sum\limits_{i=0}^{N-1}\norme{\nabla J(x_i)}^2)_{N \in \mathbb{N}} $ converge, ce qui implique :
 $$ \lim\limits_{k \rightarrow +\infty} \norme{\nabla J(x_k)} = 0 $$
-L'étude plus détaillée de différents algorithmes de descente qui utilisent différentes méthodes de recherche linéaire pour optimiser $ \varphi $ et le choix d'une direction ainsi que leurs convergences sort du cadre de ce projet.
+\begin{Def}
+ On considère $ (x_n)_{n \in \mathbb{N}} $, la suite des itérés donnés par un algorithme convergent. On note $ x^\ast $ la limite de la suite $ (x_n)_{n \in \mathbb{N}} $ et on suppose que $ x_k \neq x^\ast $, pour tout $ k \in \mathbb{N} $. La convergence de l’algorithme est alors dite :
+ \begin{itemize}
+  \item linéaire si l'erreur décroît linéairement i.e. :
+        $$ \exists \tau \in ]0,1[ \ \lim_{k \rightarrow +\infty} \frac{\norme{x_{k+1} - x^\ast}}{\norme{x_k - x^\ast}} = \tau $$
+  \item superlinéaire si :
+        $$ \lim_{k \rightarrow +\infty} \frac{\norme{x_{k+1} - x^\ast}}{\norme{x_k - x^\ast}} = 0 $$
+  \item d'ordre $ p $ si :
+        $$ \exists \tau \geq 0 \ \lim_{k \rightarrow +\infty} \frac{\norme{x_{k+1} - x^\ast}}{\norme{x_k - x^\ast}^p} = \tau $$
+        En particulier, si $ p = 2 $, la convergence est dite quadratique.
+ \end{itemize}
+\end{Def}
+L'étude plus détaillée de différents algorithmes de descente qui utilisent différentes méthodes de recherche linéaire pour optimiser $ \varphi $ ainsi que leurs convergences sort du cadre de ce projet.
 
 \section{Méthode Newtonienne}
 
@@ -485,8 +516,8 @@ En supposant $ J $ de classe $ \mathcal{C}^2 $ et la matrice hessienne $ H[J](x_
 $$ x_{k+1} = x_k - H[J](x_k)^{-1} \nabla J(x_k), $$
 où $ d_k = -H[J](x_k)^{-1} \nabla J(x_k) $ est appelée direction de Newton. La direction $ d_k $ est également l’unique solution du problème :
 $$ \underset{d \in \mathbb{R}^n}{\mathrm{argmin}} \ J(x_k) + \langle \nabla J(x_k),d \rangle + \frac{1}{2}\langle H[J](x_k)d,d \rangle $$
-Autrement dit, $ d_k $ est le point de minimum global de l’approximation de second ordre de $ J $ au voisinage du point courant $ x_k $.
-A condition que la matrice $ H[J](x_k) $ soit définie positive à chaque itération, la méthode de Newton est bien une méthode de descente à pas fixe égal à $ 1 $.
+Autrement dit, $ d_k $ est le point de minimum global de l’approximation quadratique de $ J $ au voisinage du point courant $ x_k $.
+À condition que la matrice $ H[J](x_k) $ soit définie positive à chaque itération, la méthode de Newton est bien une méthode de descente à pas fixe égal à $ 1 $.
 \newline
 Les propriétés remarquables de cet algorithme sont :
 
@@ -522,12 +553,34 @@ Nous ne répondrons pas à ces questions qui sont hors du cadre de ce projet. Ce
 
 \section{Méthode PQS (ou SQP)}
 
-Nous supposons les fonctions $ J,g,h $ à valeurs réelles et de classe $ \mathcal{C}^1 $.
-Trouver une solution d’un problème d’optimisation sous contraintes fonctionnelles consiste à déterminer un point optimal $ x^\ast $ et des multiplicateurs associés $ (\lambda^\ast,\mu^\ast) $. Deux grandes familles de méthodes peuvent être définies pour la résolution des problèmes d’optimisation sous contraintes : les méthodes primales et les méthodes duales. Les approches primales se concentrent sur la détermination du point $ x^\ast $, les multiplicateurs $ (\lambda,\mu) $ ne servant souvent qu’à vérifier l’optimalité de $ x^\ast $. Les méthodes duales quant à elles mettent l’accent sur la recherche d’un multiplicateur en travaillant sur un problème d’optimisation déduit du problème initial par \textit{dualité}.
+Nous supposons les fonctions $ J,g,h $ à valeurs réelles et de classe $ \mathcal{C}^1 $. Trouver une solution d’un problème d’optimisation sous contraintes fonctionnelles consiste à déterminer un point optimal $ x^\ast $ et des multiplicateurs associés $ (\lambda^\ast,\mu^\ast) $. Deux grandes familles de méthodes peuvent être définies pour la résolution des problèmes d’optimisation sous contraintes : les méthodes primales et les méthodes duales. Les approches primales se concentrent sur la détermination du point $ x^\ast $, les multiplicateurs $ (\lambda,\mu) $ ne servant souvent qu’à vérifier l’optimalité de $ x^\ast $. Les méthodes duales quant à elles mettent l’accent sur la recherche des multiplicateurs en travaillant sur un problème d’optimisation déduit du problème initial par \textit{dualité}.
 
-\subsection{Algorithmes newtoniens}
+\subsection{Problème quadratique sous contraintes linéaires}
 
-Les algorithmes newtoniens sont basés sur la linéarisation d’équations caractérisant les solutions que l’on cherche, fournies par les conditions d’optimalité d’ordre $ 1 $. Ces algorithmes sont \textit{primaux-duaux} dans le sens où ils génèrent à la fois une suite primale $ (x_k )_{k \in \mathbb{N}} $ convergeant vers une solution $ \overline{x} $ du problème considéré, et une suite duale $ (\lambda_k)_{k \in \mathbb{N}} $ (resp. $ ((\lambda_k, \mu_k))_{k \in \mathbb{N}} $) de multiplicateurs convergeant vers un multiplicateur optimal $ \overline{\lambda} $ (resp. $ (\overline{\lambda},\overline{\mu}) $) associé à $ \overline{x} $.
+Nous introduisons les différentes approches développées pour la résolution des problèmes de programmation quadratique avec contraintes d'égalités et d’inégalités linéaires.
+\newline
+Ce type de problème quadratique se pose sous la forme :
+$$
+ \mathcal{PQ} \left \{
+ \begin{array}{l}
+  \displaystyle\min_{x \in \mathbb{R}^n} c^\top x + \frac{1}{2} x^\top \mathcal{Q} x \\
+  A^\top x + b \leq 0                                                                \\
+  A^{\prime^\top} x + b^\prime = 0
+ \end{array}
+ \right .
+$$
+où $$ \mathcal{Q} \in \mathcal{M}_n(\mathbb{R}) \ symétrique, c \in \mathbb{R}^n, A \in  \mathcal{M}_{n,p}(\mathbb{R}), b \in \mathbb{R}^p, A^\prime \in \mathcal{M}_{n,q}(\mathbb{R}), b^\prime \in \mathbb{R}^q $$
+Or
+$$  A^{\prime^\top} x + b^\prime = 0 \iff A^{\prime^\top} x + b^\prime \leq 0 \land   -A^{\prime^\top} x - b^\prime \leq 0 $$
+Donc le problème se ramène à :
+
+\subsubsection{Algorithme 1}
+
+\subsubsection{Algorithme 2}
+
+\subsection{Algorithmes Newtoniens}
+
+Les algorithmes newtoniens sont basés sur la linéarisation d’équations caractérisant les solutions que l’on cherche, fournies par les conditions d’optimalité d’ordre $ 1 $. Ces algorithmes sont \textit{primaux-duaux} dans le sens où ils génèrent à la fois une suite primale $ (x_k )_{k \in \mathbb{N}} $ convergeant vers une solution $ \overline{x} $ du problème considéré, et une suite duale $ (\lambda_k)_{k \in \mathbb{N}} $ (resp. $ ((\lambda_k, \mu_k))_{k \in \mathbb{N}} $) de multiplicateurs convergeant vers un multiplicateur optimal $ \overline{\lambda} $ (resp. $ (\overline{\lambda},\overline{\mu})) $ associé à $ \overline{x} $.
 
 \subsection{Algorithme PQS}
 
@@ -568,7 +621,7 @@ $$ \begin{pmatrix}
   h(x_k)
  \end{pmatrix}  $$
 où $ D_h(x) $ désigne la matrice jacobienne de l’application $ h : \mathbb{R}^n \longrightarrow \mathbb{R}^q $ définie par :
-$$ D_h(x)^\top = \begin{bmatrix} \nabla h_1(x)\ldots\nabla h_q(x) \end{bmatrix} $$
+$$ D_h(x)^\top = \begin{bmatrix} \nabla h_1(x)^\top\ldots\nabla h_q(x)^\top \end{bmatrix} $$
 Posons : $ H_k = H_x[L](x_k,\lambda_k), \ d = x_{k+1} - x_k $ et $ \mu = \lambda_{k+1} $. L'itération s'écrit donc :
 $$ \begin{pmatrix}
   H_k      & D_h(x_k)^\top \\
@@ -582,7 +635,7 @@ $$ \begin{pmatrix}
  \end{pmatrix} $$
 et est bien définie à condition que la matrice $ H_x[L](x_k,\lambda_k) $ soit inversible. Ce sera le cas si :
 \begin{enumerate}[label=(\roman*)]
- \item Les colonnes $ \nabla h_1(x_k),\ldots,\nabla h_q(x_k) $ de $ D_h(x_k)^\top $ sont linéairement indépendants : c’est l’hypothèse de qualification des contraintes.
+ \item Les colonnes $ \nabla h_1(x_k)^\top,\ldots,\nabla h_q(x_k)^\top $ de $ D_h(x_k)^\top $ sont linéairement indépendants : c’est la condition première de \textit{KTT} ou condition de qualification des contraintes.
  \item Quel que soit $ d \neq 0 $ tel que $ D_h(x_k)d = 0, \ d^\top H_k d > 0 $ : c’est la condition suffisante d’optimalité du second ordre dans le cas de contraintes d’égalité.
 \end{enumerate}
 Revenons à l’itération. Elle s’écrit encore :
@@ -658,8 +711,14 @@ $$
 $$
 où $ J: \mathbb{R}^n \longrightarrow \mathbb{R} $, $g: \mathbb{R}^n \longrightarrow \mathbb{R}^p$ et $h: \mathbb{R}^n \longrightarrow \mathbb{R}^q$ sont supposées au moins différentiables.
 \newline
-Selon le même principe qu’avec contraintes d’égalité seules, on linéarise les contraintes et on utilise une approximation quadratique du Lagrangien :
+Selon le même principe qu’avec contraintes d’égalité seules, on linéarise les contraintes et on utilise une approximation quadratique du Lagrangien à l'aide de développements de Taylor-Young en $ x_k $ et $ (x_k,\lambda_k,\mu_k) $ respectivement :
 $$ L(x,\lambda,\mu) = J(x) + \lambda^\top g(x) + \mu^\top h(x), \ \lambda \in \mathbb{R}_+^p \land \mu \in \mathbb{R}^q $$
+Soit à l'ordre 2 pour le Lagrangien :
+$$ L(x,\lambda,\mu) \approx L(x_k,\lambda_k,\mu_k) + \nabla L(x_k,\lambda_k,\mu_k)^\top (x - x_k) + \frac{1}{2} (x - x_k)^\top H[L](x_k,\lambda_k,\mu_k) (x - x_k) $$
+et à l'ordre 1 pour les contraintes :
+$$ g(x) \approx g(x_k) + \nabla g(x_k)^\top(x - x_k) $$
+$$ h(x) \approx h(x_k) + \nabla h(x_k)^\top(x - x_k) $$
+En posant $ d = x - x_k $ et $ H_k = H[L](x_k,\lambda_k,\mu_k) $, on obtient le sous problème quadratique $ \mathcal{PQ}_k $ :
 
 \hrulefill
 \newline
@@ -677,7 +736,7 @@ ALGORITHME PQS AVEC CONSTRAINTES D'ÉGALITÉ ET D'INEGALITÉ.
                \mathcal{PQ}_k \left \{
                \begin{array}{l}
                 \displaystyle\min_{d \in \mathbb{R}^n} \nabla J(x_k)^\top d + \frac{1}{2}d^\top H_k d \\
-                g_j(x_k) + \nabla g_j(x_k)^\top d = 0, \ \forall j \in \{1,\ldots,p\}                 \\
+                g_j(x_k) + \nabla g_j(x_k)^\top d \leq 0 \\, \ \forall j \in \{1,\ldots,p\}                 \\
                 h_i(x_k) + \nabla h_i(x_k)^\top d = 0, \ \forall i \in \{1,\ldots,q\}
                \end{array}
                \right .
@@ -694,6 +753,299 @@ Afin que le sous-programme quadratique $ \mathcal{PQ}_k $ admette une unique sol
 \newline
 Etant une méthode newtonienne, l’algorithme PQS converge localement quadratiquement pourvu que les points initiaux  $ (x_0,\lambda_0 ) $ (resp. $ (x_0,\lambda_0,\mu_0) $) soient dans un voisinage d’un point stationnaire $ \overline{x} $ et de ses multiplicateurs associés $ \overline{\lambda} $ (resp. $ (\overline{\lambda},\overline{\mu}) $). Bien entendu, il est possible de globaliser l’algorithme en ajoutant une étape de recherche linéaire.
 
+\subsection{Stratégie d'approximation de la hessienne}
+
+\subsubsection{Équation de sécante et approximation}
+
+L'approximation $ H_k $  de la hessienne du Lagrangien peut être obtenu par la relation :
+$$ \nabla L(x_{k+1},\lambda_{k+1},\mu_{k+1}) - \nabla L(x_{k},\lambda_{k+1},\mu_{k+1}) \approx H[L](x_{k+1},\lambda_{k+1},\mu_{k+1})(x_{k+1} - x_k) $$
+On construit une approximation $ H_{k+1} $ de $ H[L](x_{k+1},\lambda_{k+1},\mu_{k+1}) $ comme solution de l’équation :
+$$ H_{k+1}(x_{k+1} - x_k) = \nabla L(x_{k+1},\lambda_{k+1},\mu_{k+1}) - \nabla L(x_{k},\lambda_{k+1},\mu_{k+1}) $$
+appelée équation de sécante ou équation de quasi-Newton.
+\newline
+De façon similaire, on peut construire une approximation $ B_{k+1} $ de $ H[L](x_{k+1},\lambda_{k+1},\mu_{k+1})^{-1} $ comme solution de l’équation :
+$$ B_{k+1}(\nabla L(x_{k+1},\lambda_{k+1},\mu_{k+1}) - \nabla L(x_{k},\lambda_{k+1},\mu_{k+1})) = x_{k+1} - x_k $$
+Dans les deux cas, les équations de quasi-Newton forment un système sous-déterminé à $ n $ équations et $ n^2 $ inconnues. Il existe donc une infinité de matrices $ H_{k+1} $ pouvant convenir.
+\newline
+Une stratégie commune est de calculer $ (x_{k+1},\lambda_{k+1},\mu_{k+1}) $ pour une matrice $ H_k $ donnée et faire une mise à jour de $ H_k $ de rang 1 ou 2 :
+$$ H_{k+1} = H_k + U_k $$
+
+\subsubsection{Mises à jour DFP et BFGS}
+
+\subsection{Exemple d'utilisation de PQS}
+
+Considérons le problème $ \mathcal{P} $ suivant :
+$$
+ \mathcal{P} \left \{
+ \begin{array}{l}
+  \displaystyle\min_{(x,y,z) \in \mathbb{R}^3} J(x,y,z) = x^2 + y^2 + z^2 -r^2       \\
+  g(x,y,z) = (g_1(x,y,z), g_2(x,y,z)) = (x^2 + y^2 - r_1^2, x^2 + z^2 -r_2^2) \leq 0 \\
+ \end{array}
+ \right .
+$$
+où $$ (r,r_1,r_2) \in \mathbb{R}_+^3. $$
+\textit{Entrées} : $ J $ et $ g $ de classe $ \mathcal{C}^2 $, $ \varepsilon = 0.01 $ la précision, $ (x_0,y_0,z_0) = $ point initial et $ (\lambda_{0_1},\lambda_{0_2}) = $ multiplicateur initial.
+\newline
+Le Lagrangien $ L $ de $ \mathcal{P} $ : $$ L((x,y,z),(\lambda_1,\lambda_2)) = x^2 + y^2 + z^2 -r^2 + \lambda_1(x^2 + y^2 - r_1^2) + \lambda_2(x^2 + z^2 -r_2^2). $$
+\newline
+Le gradient de $ J $ : $$ \nabla J(x,y,z) = (\frac{\partial J}{\partial x}(x,y,z),\frac{\partial J}{\partial y}(x,y,z),\frac{\partial J}{\partial z}(x,y,z)) = (2x,2y,2z). $$
+\newline
+Le gradient de $ g $ : $$ \nabla g(x,y,z) = (\nabla g_1(x,y,z),\nabla g_2(x,y,z)) $$
+$$ = ((\frac{\partial g_1}{\partial x}(x,y,z),\frac{\partial g_1}{\partial y}(x,y,z),\frac{\partial g_1}{\partial z}(x,y,z)),(\frac{\partial g_2}{\partial x}(x,y,z),\frac{\partial g_2}{\partial y}(x,y,z),\frac{\partial g_2}{\partial z}(x,y,z)) $$
+$$ = ((2x,2y,0),(2x,0,2z)). $$
+\newline
+Le gradient du Lagrangien $ L $ :
+$$ \nabla L((x,y,z),(\lambda_1,\lambda_2)) = \nabla J(x,y,z) + \lambda_1 \nabla g_1(x,y,z) + \lambda_2 \nabla g_2(x,y,z)) $$
+\newline
+La matrice hessienne de $ J $ : $$ H[J](x,y,z) =
+ \begin{pmatrix}
+  \frac{\partial^2 J}{\partial^2 x}(x,y,z)         & \frac{\partial^2 J}{\partial x\partial y}(x,y,z) & \frac{\partial^2 J}{\partial x\partial z}(x,y,z) \\
+  \frac{\partial^2 J}{\partial y\partial x}(x,y,z) & \frac{\partial^2 J}{\partial^2 y}(x,y,z)         & \frac{\partial^2 J}{\partial y\partial z}(x,y,z) \\
+  \frac{\partial^2 J}{\partial z\partial x}(x,y,z) & \frac{\partial^2 J}{\partial z\partial y}(x,y,z) & \frac{\partial^2 J}{\partial^2 z}(x,y,z)         \\
+ \end{pmatrix} =
+ \begin{pmatrix}
+  2 & 0 & 0 \\
+  0 & 2 & 0 \\
+  0 & 0 & 2 \\
+ \end{pmatrix} = 2Id_{\mathbb{R}^3} $$
+On en déduit que $ H[J](x,y,z) $ est inversible et que $ H[J](x,y,z)^{-1} = \frac{1}{2}Id_{\mathbb{R}^3} $.
+
+\hrulefill
+
+\subsection{Trace d'éxécution de PQS}
+
+Utilisons le problème $ \mathcal{P} $ précédent :
+
+$$
+ \mathcal{P} \left \{
+ \begin{array}{l}
+  \displaystyle\min_{(x,y,z) \in \mathbb{R}^3} J(x,y,z) = x^2 + y^2 + z^2 -r^2       \\
+  g(x,y,z) = (g_1(x,y,z), g_2(x,y,z)) = (x^2 + y^2 - r_1^2, x^2 + z^2 -r_2^2) \leq 0 \\
+ \end{array}
+ \right .
+$$
+où $$ (r,r_1,r_2) \in \mathbb{R}_+^3. $$
+\textit{Entrées} : $ J $ et $ g $ de classe $ \mathcal{C}^2 $, $ \varepsilon = 0.01 $, $ (x_0,y_0,z_0) = (100, 100 ,0)$  et $(\lambda_{0_1},\lambda_{0_2}) = (1 , 1)$, les rayons : $r= 100$  et $r1 = r2 = 10$.
+\newline
+Le Lagrangien $ L $ de $ \mathcal{P} $ : $$ L((x,y,z),(\lambda_1,\lambda_2)) = x^2 + y^2 + z^2 -r^2 + \lambda_1(x^2 + y^2 - r_1^2) + \lambda_2(x^2 + z^2 -r_2^2). $$
+\newline
+Le Lagrangien $ L $ de $ \mathcal{P} $ avec les valeurs :
+ $ L((100,100,0),(1,1)) = 100^2 + 100^2 + 0^2 -100^2 + 1 * (100^2 +100^2 - 10^2) + \lambda_2(100^2 + 100^2 -10^2). $
+ $ L((100,100,0),(1,1)) = 1000 + 1000 - 1000 + (1000 + 1000 - 100) + (1000 + 1000 -100). $
+ $ L((100,100,0),(1,1)) = 4800. $
+
+\newpage
+\begin{algorithmfloat}[#Algo 1]
+ \caption {Trace d'éxécution du PQS du problème $ \mathcal{P} $}
+ \begin{algorithmic}
+ \REQUIRE $g(x_0,y_0,z_0)\leq 0$, $(x_0,y_0,z_0) = (10, 10 ,10)$
+ \ENSURE $\min_{(x,y,z) \in \mathbb{R}^3} J(x,y,z) = x^2 + y^2 + z^2 -r^2$ and \newline $g(x,y,z) = (g_1(x,y,z), g_2(x,y,z)) = (x^2 + y^2 - r_1^2, x^2 + z^2 -r_2^2) \leq 0 $
+ \STATE \textbf{Data :}
+ \STATE $k \leftarrow 0, (x_k, y_k, z_k)  \leftarrow (100, 100, 0), r \leftarrow 100$
+ \STATE $r_1 = r2 \leftarrow 10, \varepsilon \leftarrow 0.01$
+ \STATE $\lambda_1 = \lambda_2 = 1$
+ \STATE $ H[J](x,y,z)^{-1}\leftarrow  \begin{pmatrix}
+  0.5 & 0 & 0 \\
+  0 & 0.5 & 0 \\
+  0 & 0 & 0.5 \\ \end{pmatrix} $
+\newline
+
+ \STATE{//Calcule du gradient de $ J $ :}
+ \STATE $ \nabla J(x,y,z) = (2x_k,2y_k,2z_k)$  \hfill $  // résultat : (100,100,0) $
+\newline
+ \STATE {//calcule des deux sous partie de du gradient de $ g $: }
+ \STATE $ \nabla g_1(x_a,y_a,z_a) = ((2x_a,2y_a,0)$  \hfill $ //résultat : (20, 20, 0)$
+ \STATE $ \nabla g_2(x_a,y_a,z_a) = (2x_a,0,2z_a))$  \hfill $ //résultat : (20, 0, 20)$
+ \STATE $ \nabla g(x_k,y_k,z_k) = (\nabla g_1(x_k,y_k,z_k), \nabla g_2(x_k,y_k,z_k))$
+\newline
+ \WHILE{$ (\norme{\nabla L(x_k,\lambda_k,\mu_k)} > \varepsilon $ or k $ \leq 10)$}
+
+ \STATE { //première itération :}
+
+\STATE {//Calcule du gradient de $ L $ : }
+\STATE $\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2)) = \nabla J(x_k,y_k,z_k) + \lambda_1 \nabla g_1(x_a,y_a,z_a) + \lambda_2 \nabla g_2(x_a,y_a,z_a)) $ \hfill $// résultat : (220, 220, 40)$
+ \STATE $  \nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2)) = (x_L , y_L, z_L) $
+\newline
+ \STATE {//Calcule de la direction de la pente dk (méthode de Newton) : }
+ \STATE $ d_k = -H[J](x,y,z)^{-1}* J(x,y,z)$ \hfill $ //résultat : (-(50,50,0))$
+ \newline
+ \STATE {//Calcul nouvelles valeurs des coordonnées}
+ \STATE $ (x_{k+1},y_{k+1},z_{k+1}) = (x_k,y_k,z_k)+ d_k $ \hfill $ //résultat : (50,50,0)$
+ \newline
+ \STATE {//Incrémentation de k}
+ \STATE $ k \leftarrow k+1$\hfill $ //k = 1$
+\newline
+
+ \STATE {//Deuxième itération :}
+ \STATE{//Calcule du gradient de $ J $ :}
+ \STATE $ \nabla J(x,y,z) = (2x_k,2y_k,2z_k)$  \hfill $  // résultat : (100,100,0) $
+\newline
+\STATE {//Calcule du gradient de $ L $ : }
+\STATE $\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2)) = \nabla J(x_k,y_k,z_k) + \lambda_1 \nabla g_1(x_a,y_a,z_a) + \lambda_2 \nabla g_2(x_a,y_a,z_a)) $ \hfill $// résultat : (120, 120, 0)$
+ \STATE $  \varepsilon _1 = \norme{\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2))}$
+\newline
+ \STATE {//Calcule de la direction de la pente dk (méthode de Newton) : }
+ \STATE $ d_k = -H[J](x,y,z)^{-1}* J(x,y,z)$ \hfill $ //résultat : (-(25,25,0))$
+ \STATE {//Calcul nouvelles valeurs des coordonnées}
+ \STATE $ (x_{k+1},y_{k+1},z_{k+1}) = (x_k,y_k,z_k)+ d_k $ \hfill $ //résultat : (25,25,0)$
+ \newline
+ \STATE {//Incrémentation de k}
+\STATE $ k \leftarrow k+1$\hfill $ //k = 2$
+\newline
+
+\STATE {//Troisième itération :}
+\STATE{//Calcule du gradient de $ J $ :}
+\STATE $ \nabla J(x,y,z) = (2x_k,2y_k,2z_k)$  \hfill $  // résultat : (50,50,0) $
+\newline
+\STATE {//Calcule du gradient de $ L $ : }
+\STATE $\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2)) = \nabla J(x_k,y_k,z_k) + \lambda_1 \nabla g_1(x_a,y_a,z_a) + \lambda_2 \nabla g_2(x_a,y_a,z_a)) $ \hfill $// résultat : (70, 70, 0)$
+\STATE $  \varepsilon _1 = \norme{\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2))}$
+\newline
+\STATE {//Calcule de la direction de la pente dk (méthode de Newton) : }
+\STATE $ d_k = -H[J](x,y,z)^{-1}* J(x,y,z)$ \hfill $ //résultat : (-(12.5,12.5,0))$
+\STATE {//Calcul nouvelles valeurs des coordonnées}
+\newline
+\STATE $ (x_{k+1},y_{k+1},z_{k+1}) = (x_k,y_k,z_k)+ d_k $ \hfill $ //résultat : (12.5,12.5,0)$
+\STATE {//Incrémentation de k}
+\STATE $ k \leftarrow k+1$\hfill $ //k = 3$
+\newline
+
+\STATE {//Quatrième itération :}
+\STATE{//Calcule du gradient de $ J $ :}
+\STATE $ \nabla J(x,y,z) = (2x_k,2y_k,2z_k)$  \hfill $  // résultat : (25,25,0) $
+\newline
+\STATE {//Calcule du gradient de $ L $ : }
+\STATE $\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2)) = \nabla J(x_k,y_k,z_k) + \lambda_1 \nabla g_1(x_a,y_a,z_a) + \lambda_2 \nabla g_2(x_a,y_a,z_a)) $ \hfill $// résultat : (45, 45, 0)$
+\STATE $  \varepsilon _1 = \norme{\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2))}$
+\newline
+\STATE {//Calcule de la direction de la pente dk (méthode de Newton) : }
+\STATE $ d_k = -H[J](x,y,z)^{-1}* J(x,y,z)$ \hfill $ //résultat : (-(6.25,6.25,0))$
+\newline
+\STATE {//Calcul nouvelles valeurs des coordonnées}
+\STATE $ (x_{k+1},y_{k+1},z_{k+1}) = (x_k,y_k,z_k)+ d_k $ \hfill $ //résultat : (6.25,6.25,0)$
+\STATE {//Incrémentation de k}
+\newline
+\STATE $ k \leftarrow k+1$\hfill $ //k = 4$
+\STATE $ $
+
+\STATE {//Cinquième itération :}
+\STATE{//Calcule du gradient de $ J $ :}
+\STATE $ \nabla J(x,y,z) = (2x_k,2y_k,2z_k)$  \hfill $  // résultat : (12.5,12.5,0) $
+\newline
+\STATE {//Calcule du gradient de $ L $ : }
+\STATE $\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2)) = \nabla J(x_k,y_k,z_k) + \lambda_1 \nabla g_1(x_a,y_a,z_a) + \lambda_2 \nabla g_2(x_a,y_a,z_a)) $ \hfill $// résultat : (32.5, 32.5, 0)$
+\STATE $  \varepsilon _1 = \norme{\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2))}$
+\newline
+\STATE {//Calcule de la direction de la pente dk (méthode de Newton) : }
+\STATE $ d_k = -H[J](x,y,z)^{-1}* J(x,y,z)$ \hfill $ //résultat : (-(3.125,3.125,0))$
+\newline
+\STATE {//Calcul nouvelles valeurs des coordonnées}
+\STATE $ (x_{k+1},y_{k+1},z_{k+1}) = (x_k,y_k,z_k)+ d_k $ \hfill $ //résultat : (3.125,3.125,0)$
+\newline
+\STATE {//Incrémentation de k}
+\STATE $ k \leftarrow k+1$\hfill $ //k = 5$
+\newline
+
+\STATE {//Sixième itération :}
+\STATE{//Calcule du gradient de $ J $ :}
+\STATE $ \nabla J(x,y,z) = (2x_k,2y_k,2z_k)$  \hfill $  // résultat : (6.25,6.25,0) $
+\newline
+\STATE {//Calcule du gradient de $ L $ : }
+\STATE $\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2)) = \nabla J(x_k,y_k,z_k) + \lambda_1 \nabla g_1(x_a,y_a,z_a) + \lambda_2 \nabla g_2(x_a,y_a,z_a)) $ \hfill $// résultat : (26.25, 26.25, 0)$
+\STATE $  \varepsilon _1 = \norme{\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2))}$
+\newline
+\STATE {//Calcule de la direction de la pente dk (méthode de Newton) : }
+\STATE $ d_k = -H[J](x,y,z)^{-1}* J(x,y,z)$ \hfill $ //résultat : (-(1.5625,1.5625,0))$
+\STATE {//Calcul nouvelles valeurs des coordonnées}
+\newline
+\STATE $ (x_{k+1},y_{k+1},z_{k+1}) = (x_k,y_k,z_k)+ d_k $ \hfill $ //résultat : (1.5625,1.5625,0)$
+\STATE {//Incrémentation de k}
+\newline
+\STATE $ k \leftarrow k+1$\hfill $ //k = 6$
+\newline
+
+\STATE {//Septième itération :}
+\STATE{//Calcule du gradient de $ J $ :}
+\STATE $ \nabla J(x,y,z) = (2x_k,2y_k,2z_k)$  \hfill $  // résultat : (3.125, 3.125, 0) $
+\newline
+\STATE {//Calcule du gradient de $ L $ : }
+\STATE $\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2)) = \nabla J(x_k,y_k,z_k) + \lambda_1 \nabla g_1(x_a,y_a,z_a) + \lambda_2 \nabla g_2(x_a,y_a,z_a)) $ \hfill $// résultat : (23.125, 23.125, 0)$
+\STATE $  \varepsilon _1 = \norme{\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2))}$
+\newline
+\STATE {//Calcule de la direction de la pente dk (méthode de Newton) : }
+\STATE $ d_k = -H[J](x,y,z)^{-1}* J(x,y,z)$ \hfill $ //résultat : (-(0.78125,0.78125,0))$
+\STATE {//Calcul nouvelles valeurs des coordonnées}
+\newline
+\STATE $ (x_{k+1},y_{k+1},z_{k+1}) = (x_k,y_k,z_k)+ d_k $ \hfill $ //résultat : (0.78125,0.78125,0)$
+\STATE {//Incrémentation de k}
+\newline
+\STATE $ k \leftarrow k+1$\hfill $ //k = 7$
+\newline
+
+\STATE {//Huitième itération :}
+\STATE{//Calcule du gradient de $ J $ :}
+\STATE $ \nabla J(x,y,z) = (2x_k,2y_k,2z_k)$  \hfill $  // résultat : (1.5625, 1.5625, 0) $
+\newline
+\STATE {//Calcule du gradient de $ L $ : }
+\STATE $\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2)) = \nabla J(x_k,y_k,z_k) + \lambda_1 \nabla g_1(x_a,y_a,z_a) + \lambda_2 \nabla g_2(x_a,y_a,z_a)) $ \hfill $// résultat : (21.5625, 21.5625, 0)$
+\STATE $  \varepsilon _1 = \norme{\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2))}$
+\newline
+\STATE {//Calcule de la direction de la pente dk (méthode de Newton) : }
+\STATE $ d_k = -H[J](x,y,z)^{-1}* J(x,y,z)$ \hfill $ //résultat : (-(0.390625,0.390625,0))$
+\newline
+\STATE {//Calcul nouvelles valeurs des coordonnées}
+\STATE $ (x_{k+1},y_{k+1},z_{k+1}) = (x_k,y_k,z_k)+ d_k $ \hfill $ //résultat : (0.390625,0.390625,0)$
+\newline
+\STATE {//Incrémentation de k}
+\STATE $ k \leftarrow k+1$\hfill $ //k = 8$
+\newline
+
+\STATE {//neuvième itération :}
+\STATE{//Calcule du gradient de $ J $ :}
+\STATE $ \nabla J(x,y,z) = (2x_k,2y_k,2z_k)$  \hfill $  // résultat : (0.78125, 0.78125, 0) $
+\newline
+\STATE {//Calcule du gradient de $ L $ : }
+\STATE $\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2)) = \nabla J(x_k,y_k,z_k) + \lambda_1 \nabla g_1(x_a,y_a,z_a) + \lambda_2 \nabla g_2(x_a,y_a,z_a)) $ \hfill $// résultat : (20.78125, 20.78125, 0)$
+\STATE $  \varepsilon _1 = \norme{\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2))}$
+\newline
+\STATE {//Calcule de la direction de la pente dk (méthode de Newton) : }
+\STATE $ d_k = -H[J](x,y,z)^{-1}* J(x,y,z)$ \hfill $ //résultat : (-(0.1953125,0.1953125,0))$
+\newline
+\STATE {//Calcul nouvelles valeurs des coordonnées}
+\STATE $ (x_{k+1},y_{k+1},z_{k+1}) = (x_k,y_k,z_k)+ d_k $ \hfill $ //résultat : (0.1953125,0.1953125,0)$
+\newline
+\STATE {//Incrémentation de k}
+\STATE $ k \leftarrow k+1 \hfill  //k = 9$
+\newline
+
+\STATE {//Dixième itération :}
+\STATE{//Calcule du gradient de $ J $ :}
+\STATE $ \nabla J(x,y,z) = (2x_k,2y_k,2z_k)$  \hfill $  // résultat : (0.390625, 0.390625, 0) $
+\newline
+\STATE {//Calcule du gradient de $ L $ : }
+\STATE $\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2)) = \nabla J(x_k,y_k,z_k) + \lambda_1 \nabla g_1(x_a,y_a,z_a) + \lambda_2 \nabla g_2(x_a,y_a,z_a)) $ \hfill $// résultat : (20.390625, 20.390625, 0)$
+\STATE $  \varepsilon _1 = \norme{\nabla L((x_k,y_k,z_k),(\lambda_1,\lambda_2))}$
+\newline
+\STATE {//Calcule de la direction de la pente dk (méthode de Newton) : }
+\STATE $ d_k = -H[J](x,y,z)^{-1}* J(x,y,z)$ \hfill $ //résultat : (-(0.097665625,0.097665625,0))$
+\newline
+\STATE {//Calcul nouvelles valeurs des coordonnées}
+\STATE $ (x_{k+1},y_{k+1},z_{k+1}) = (x_k,y_k,z_k)+ d_k $ \hfill $ //résultat : (0.097665625,0.097665625,0)$
+\newline
+\STATE {//Incrémentation de k}
+\STATE $ k \leftarrow k+1$\hfill $ //k = 10$
+\newline
+\STATE {// Fin de la boucle "while" car nous avons atteint k =10, condition mettant fin à la //boucle}
+\newline
+
+ \ENDWHILE
+
+\end{algorithmic}
+\end{algorithmfloat}
+
+
+\hrulefill
+
 \bibliographystyle{plain}
 \bibliography{stdlib_sbphilo}