Otimização/Imprimir

Origem: Wikilivros, livros abertos por um mundo aberto.
Ir para: navegação, pesquisa


Índice

Situação inicial


Embora se possa trabalhar com mínimos e máximos, ao longo dos próximos capítulos só trabalharemos com mínimos, pois achar o máximo de uma função f é equivalente a achar o mínimo da função  -f.

[editar] Mínimo global

Sejam  D \subset \mathbb{R}^n e  f:D \rightarrow \mathbb{R}. Para encontrarmos o mínimo global, devemos encontrar o  \operatorname{min} \; f(x), \forall \; x \in D .

Definição

Dizemos que um ponto  \bar{x} \in D é mínimo global, se  f(\bar{x}) \le f(x), \forall \; x \in D.

[editar] Máximo global

Seja  D \subset \mathbb{R}^n e  f:D \rightarrow \mathbb{R}. Para encontrarmos o máximo global, devemos encontrar o  \operatorname{max} \; f(x), \forall \; x \in D .

Definição

Dizemos que um ponto  \bar{x} \in D é máximo global, se  f(\bar{x}) \ge f(x), \forall \; x \in D.

[editar] Mínimo local

Seja  D \subset \mathbb{R}^n e  f:D \rightarrow \mathbb{R}. Para encontrarmos o mínimo local, devemos encontrar o  \operatorname{min} \; f(x), \forall \; x \in D \cap B_\epsilon (\bar{x}).

Definição

Dizemos que um ponto  \bar{x} \in D é mínimo local, se

 f(\bar{x}) \le f(x), \forall \; x \in D \cap B_\epsilon (\bar{x}) onde  B_\epsilon (\bar{x}) = \{ x \in D ; \;  \| x - \bar{x} \| < \epsilon \}.

[editar] Máximo local

Seja  D \subset \mathbb{R}^n e  f:D \rightarrow \mathbb{R}. Para encontrarmos o máximo local, devemos encontrar o  \operatorname{max} \; f(x), \forall \; x \in D \cap B_\epsilon (\bar{x}).

Definição

Dizemos que um ponto  \bar{x} \in D é máximo local, se  f(\bar{x}) \ge f(x), \forall \; x \in D \cap B_\epsilon (\bar{x}) onde  B_\epsilon (\bar{x}) = \{ x \in D ; \;  \| x - \bar{x} \| < \epsilon \}.

[editar] Exemplos

Seja  f: \mathbb{R}^n \rightarrow \mathbb{R}; D_1, D_2 \subset \mathbb{R} , tais que  D_2 \subset D_1 .

[editar] Exemplo 1

Mostrar que  \inf_{x \in D_1}f(x) \le \inf_{x \in D_2}f(x) .

Afirmação:  f(y)=\inf_{x \in D_1}f(x) \Rightarrow  f(y) \le f(x), \forall \; x \in D_1 e  f(z)=\inf_{x \in D_2}f(x) \Rightarrow  f(z) \le f(x), \forall \; x \in D_2.

Prova1: Tome  t \in D_2 \subset D_1 \Rightarrow t \in D_1 \Rightarrow f(y) \le f(t), \forall \; t \in D_2 \Rightarrow f(y) \le f(z) .

Prova2: Suponha por contradição que  \inf_{x \in D_2}f(x) < \inf_{x \in D_1}f(x) \Rightarrow f(z) < f(y). Mas  z \in D_2 \subset D_1 \Rightarrow z \in D_1 \Rightarrow f(y) \le f(z) . Logo  f(z) < f(y) \le f(z) . Contradição! A contradição foi supor que  \inf_{x \in D_2}f(x) < \inf_{x \in D_1}f(x) .

Portanto,  \inf_{x \in D_1}f(x) \le \inf_{x \in D_2}f(x) .

[editar] Exemplo 2

Seja  f: \mathbb{R}^n \rightarrow \mathbb{R}; D_1, D_2 \subset \mathbb{R} , tais que  \bar{x} \in D_2 \subset D_1 . Seja  f(\bar{x})\le f(x), \forall \; x \in D_1 .

Mostrar que  f(\bar{x})\le f(x), \forall \; x \in D_2 .

Suponha por contradição que  \exists \; z \in D_2 tal que  f(z) < f(\bar{x}) . Por  z \in D_2 \subset D_1 \Rightarrow f(\bar{x})\le f(z). Logo  f(z) < f(\bar{x})\le f(z) . Contradição! A contradição foi supor que  \exists \; z \in D_2 tal que  f(z) < f(\bar{x}) . Portanto  f(\bar{x})\le f(x), \forall \; x \in D_2 .


Existência de soluções globais


Definição

 M(f,D) \; é o conjuntos dos minimizadores de f em D, locais e globais.

Definição

Dizer que  \# M(f,D)=1 \; significa que o conjuntos dos minimizadores de f em D possui um mínimo e ele é global.

Definição

Seja  \bar{v} \in \; ]-\infty, \infty [; \bar{v} = \inf_{x \in D}f(x) \Rightarrow \bar{v} \in M(f,D) , onde  \bar{v} \; é chamado valor ótimo do problema e é um mínimo global.

[editar] Teorema de Weierstrass

Seja  f : D \rightarrow \mathbb{R} contínua em D compacto.

[editar] Então  M(f,D) \not = \empty

Suponha que f é ilimitada inferiormente, então  \forall \; k \in \mathbb{N}, \exist \; x^k \in D; f(x^k)<-k \Rightarrow \lim_{k \rightarrow \infty}f(x^k)=-\infty . Por outro lado, D é compacto e  \{ x^k\} \subset D . Como D é limitado, logo a  \{ x^k\} \; é limitada. Toda sequência limitada possui uma subsequência convergente. Assim  \{ x^k\} \; possui uma subsequência convergênte  \{ x^{k_j}\} \;, tal que  \; \lim_{j \rightarrow \infty}x^{k_j} = \bar{x} . Assim  f(\bar{x}) = f(\lim_{j \rightarrow \infty}x^{k_j}) = - \infty . Absurdo.

 \bar{v} = \inf_{x\in D}f(x) > - \infty , pela definição de ínfimo, dado  k \in \mathbb{N}, \; \exists \; x^k \in D tal que  \bar{v} \le f(x^k) < \bar{v}+ {1 \over k} \Rightarrow \bar{v} \le \lim_{k \rightarrow \infty}f(x^k) \le \bar{v} \Rightarrow \lim_{k \rightarrow \infty}f(x^k) = \bar{v} .

[editar] Curva de nível  L_{f,D}(c)

Definição

Seja  f:D \rightarrow \mathbb{R}, c \in \mathbb{R}, L_{f,D}(c) = \{ x\in D;  f(x) \le c \}

[editar] Corolário da curva de nível compacta

Sejam  D \in \mathbb{R}^n, f:D \rightarrow \mathbb{R} contínua em D. Se  \exists \; c \in \mathbb{R}, \empty \not =  L_{f,D}(c)\; é compacto.

[editar] Então  M(f,D) \not = \empty

Prova: Pelo Teorema de Weierstrass  M(f,D) \not = \empty , isto é,  \exists \bar{x} \in D; f(\bar{x})\le f(x), \forall \; x \in  L_{f,D}(c)\; ).

Mas se  x \in D \setminus  L_{f,D}(c)\;  \Rightarrow f(x) > c \ge f(\bar{x}) . Assim  f(\bar{x}) \le f(x), \forall \; x \in D , isto é,  M(f,D) \not = \empty .

[editar] Projeção de y sobre D

Definição

 P_D(y) = \bar{x} \Leftrightarrow \inf_{x \in D } \| x - y \| = \| \bar{x} - y \|

[editar] Corolário da projeção de y sobre D

 \empty \not = D \in \mathbb{R}^n é fechado.

[editar] Então  P_D(y) \not = \empty, \forall \; y \in \mathbb{R}^n

Tome  y \in \mathbb{R}^n e f_y:D \rightarrow \mathbb{R}; x \mapsto f_y(x)=\|x-y\| . É facil ver que  |\|x\|-\|y\|| \le \|x-y\| . Agora dado  \epsilon > 0, \exists \; \delta = \epsilon, \|x-y\| < \delta \Rightarrow |\|x\|-\|y\|| < \epsilon. Assim  f_y \; é contínua.

Por outro lado,  L_{f_y, D}(c)= \{ x \in D / f_y(x)\le c \}= \{ x \in D / \|x-y\|\le c \} = D \cap B_c(y) . Visto que  D , B_c(y) \; são fechados, temos que  D \cap B_c(y) é também fechado. Além disso, sendo  B_c(y) limitado, segue que  D \cap B_c(y) \subset B_c(y) é também limitado e conseqüentemente compacto. Como  L_{f_y, D}(c) = D \cap B_c(y) \Rightarrow L_{f_y, D}(c) é compacto.

Vimos que  f_y \; é contínua e  L_{f_y, D}(c) é compacto.. Tomando-se  c \in \mathbb{R} suficientemente grande, de tal forma que  L_{f_y, D}(c) \not = \empty . Pelo corolário da curva de nível,  M(f_y,D)\;  \not = \empty .

Mas  M(f_y,D)\;  = \{ \bar{x} \in D / f_y(\bar{x}) \le f_y(x), \forall \; x \in D \} = \{ \bar{x} \in D / \|\bar{x}-y\| \le \|x-y\|, \forall \; x \in D  \} =

 = \{ \bar{x} \in D / \inf_{x \in D}\|x-y\| = \|\bar{x}-y\| \} = P_D(y) \not = \empty

[editar] Exemplo

Seja  D= \{ x \in \mathbb{R}^n | x_1 + x_2 > 0 \} e  f: D \rightarrow \mathbb{R}, f(x)=x_1^2+x_2+{1 \over x_1 + x_2}

[editar] Mostrar que  M(f,D)\not = \empty .

Suponhamos que  L_{f,D}(c) \; é ilimitado para um  c \in \mathbb{R} \Rightarrow \exists \{ x^k\} \subset L_{f,D}(c) tal que  \lim_{k \rightarrow \infty }\|x^k\| . Se  \lim_{k \rightarrow \infty }\|x_1^k\| , isto é, dado  \lim_{k \rightarrow \infty }\|x_1^k\| . (...)


Condições de otimalidade para problemas sem restrições


Rekopis chopin.jpg Esta página é somente um esboço.
Ampliando-a você ajudará a melhorar o Wikilivros.

[editar] Teorema

 f:\mathbb{R}^n \rightarrow \mathbb{R} é diferenciável em  \bar{x} \in \mathbb{R}^n e  \bar{x} \in  M(f,B_\epsilon(\bar{x}))

[editar] Mostrar que  f'(\bar{x})=0 \;

[editar] Teorema

 f:\mathbb{R}^n \rightarrow \mathbb{R} é duas vezes diferenciável em  \bar{x} \in \mathbb{R}^n e  \bar{x} \in \;  M(f,B_\epsilon(\bar{x}))

[editar] Mostrar que  f'(\bar{x})=0 \; e  \langle f''(\bar{x})d,d \rangle \ge 0, \forall \; d \in \mathbb{R}^n

[editar] Teorema

 f:\mathbb{R}^n \rightarrow \mathbb{R} é duas vezes diferenciável em  \bar{x} \in \mathbb{R}^n ,  f'(\bar{x})=0 \; e  \langle f''(\bar{x})d,d \rangle > 0, \forall \; d \in \mathbb{R}^n \setminus \{ 0 \}

[editar] Mostrar que  \bar{x} \in  \;  M(f,B_\epsilon(\bar{x})\setminus \{ \bar{x} \})

Elementos de análise convexa


Rekopis chopin.jpg Esta página é somente um esboço.
Ampliando-a você ajudará a melhorar o Wikilivros.

[editar] Convexo

Definição

Dizemos que um conjunto  D \in \mathbb{R}^n é convexo quando  z(\alpha)=( 1 - \alpha ) x + \alpha y \in D, \forall \; x,y \in D, \forall \; \alpha \in [0,1] , onde  z(\alpha) \; é a combinação convexa de  x,y \; .

[editar] Teorema

Sejam  D \subset \mathbb{R}^n, f:\mathbb{R}^n \rightarrow \mathbb{R} um conjunto convexo e uma função diferenciável em  \bar{x} \in D . Seja também  \bar{x} \in   M(f,D \cap B_\epsilon(\bar{x}))\; .

[editar]  \langle f'(\bar{x}),x-\bar{x}\rangle \ge 0, \forall \; x \in D

[editar] Função Convexa

Seja  f:D \rightarrow \mathbb{R}, D \subset \mathbb{R}^n

Definição

Dizemos que uma função f é convexa se  f(( 1 - \alpha ) x + \alpha y) \le ( 1 - \alpha ) f(x) + \alpha f(y), \forall \; x,y \in D, \forall \; \alpha \in [0,1].

Definição

Dizemos que uma função f é estritamente convexa se  f(( 1 - \alpha ) x + \alpha y) < ( 1 - \alpha ) f(x) + \alpha f(y), \forall \; x,y \in D, x\not=y, \forall \; \alpha \in \; ]0,1[.

Definição

Dizemos que uma função f é  \lambda \; - fortemente convexa se  f(( 1 - \alpha ) x + \alpha y) \le ( 1 - \alpha ) f(x) + \alpha f(y) - \lambda \alpha (1-\alpha)\| x-y \|^2, \forall \; x,y \in D, \forall \; \alpha \in [0,1].

Definição

Dizemos que o epígrafo da função f é  E_f = \{ (x,c) / x \in D, c \in \mathbb{R}, f(x) \le c \} .

[editar] Teorema

Seja  f:D \rightarrow \mathbb{R}, D \in \mathbb{R}^n um conjunto convexo.

[editar] Mostrar que f é convexo  \Leftrightarrow E_f \; é convexo

[editar] Teorema da minimização convexa

Seja  f:D \rightarrow \mathbb{R}, D \in \mathbb{R}^n ambos convexos.

[editar] Mostrar que se  \bar{x} \in  M(f,B_\epsilon(\bar{x}))\;  \Rightarrow \bar{x} \in  M(f,D) \;

[editar] Mostrar que  M(f,D)\; é convexo

[editar] Mostrar que se f é estritamente convexa, então  \# M(f,D)) \;  \le 1 é convexo

[editar] Função Concava

Definição

Uma função  f \; é chamada concava se  (-f) \; é convexa em  D \; convexa, onde  f:D \rightarrow \mathbb{R} \; \mbox{e} \; D \subset \mathbb{R}^n

Conjuntos convexos


Rekopis chopin.jpg Esta página é somente um esboço.
Ampliando-a você ajudará a melhorar o Wikilivros.

[editar] Intersecção de conjuntos convexos é convexo

Sejam  D_j \subset \mathbb{R}^n, j \in I_n , conjuntos convexos, onde  I_n = \{ k \in \mathbb{N} / 1 \le k \le n, n \in \mathbb{N} \}

[editar] Mostrar que  D = \cap_{j \in I_n}D_j é um conjunto convexo

[editar] Conjunto Poliedral

Definição

Um conjunto é poliedral se é a intersecção de hiperplanos e semi-espaços

[editar] Um conjunto poliedral em  \mathbb{R}^n é convexo

[editar] O fecho e o interior de um conjunto convexo são convexos

[editar] A soma de convexos fechados é convexo e fechado

Sejam  D_j \subset \mathbb{R}^n, j=1,2 , conjuntos convexos e fechados. Um deles é limitado.

[editar] Mostrar que  \sum_{i=1}^{2}D_j é um conjunto convexo e fechado

[editar] Combinação convexa de p pontos

Definição

Seja  x^i \in \mathbb{R}^n, \alpha_i \in [0,1], i=1,...,p; \sum_{i=1}^{p}\alpha_i = 1 . A combinação convexa dos  x^i \in \mathbb{R}^n é o ponto \sum_{i=1}^{p}\alpha_i \cdot x^i

[editar] Teorema da combinação convexa

Um conjunto  D \in \mathbb{R}^n é convexo se, e somente se, a combinação convexa  \sum_{i=1}^{p}\alpha_ix_i \in D , \forall p \in \mathbb{N}, x^i \in D\;  \and \alpha_i \in [0,1], i=1,...,p; \sum_{i=1}^{p}\alpha_i=1,

[editar] Desigualdade de Jensen

Sejam  D \subset \mathbb{R}^n um conjunto convexo e  f:D \rightarrow \mathbb{R} uma função convexa, \forall p \in \mathbb{N}, x^i \in D\;  \and \alpha_i \in \mathbb{R}_+, i=1,...,p; \sum_{i=1}^{p}\alpha_i=1

[editar] Mostrar que  f(\sum_{i=1}^{p}\alpha_i x_i)\le \sum_{i=1}^{p}\alpha_i f(x_i)

[editar] Teorema de Carathéodory

Seja  x \in \mathbb{R}^n uma combinação convexa de pontos do conjunto  D \subset \mathbb{R}^n .

[editar] Mostrar que  \exists x_i \in D \and \alpha_i \in \mathbb{R}_+, i=1,...,n+1; x = \sum_{i=1}^{n+1}\alpha_i x_i, \sum_{i=1}^{n+1}\alpha_i =1

[editar] Fecho convexo

Definição

O fecho convexa de um conjunto qualquer D é o menor conjunto convexo que contem D e simbolizado por conv D.

Definição

O conjunto de todas as combinações convexas de pontos de D, simbolizaremos por  comb D = \{ x \in \mathbb{R}^n; x=\sum_{i=1}^{n+1}\alpha_i x_i, \sum_{i=1}^{n+1}\alpha_i =1, x_i \in D \and \alpha_i \in \mathbb{R}_+, i=1,...,p; p \in \mathbb{N} \} .

[editar] Corolário de um fecho convexo

Se  D \subset \mathbb{R}^n

[editar] Mostrar que conv D = comb D

[editar] Corolário da compacidade do conv D

Seja  D \subset \mathbb{R}^n compacto

[editar] Mostrar que conv D é compacto

Operador de projeção


[editar] Teorema da projeção

Seja  D \subset \mathbb{R}^n um conjunto convexo e fechado.

[editar] Mostrar que  \forall x \in \mathbb{R}^n, \exists \; P_D(x) e é única

[editar] Mostrar que  \bar{x} = P_D(x) \Leftrightarrow \bar{x} \in D, \langle x-\bar{x}, y-\bar{x} \rangle \le 0, \forall \;y \in D

[editar] Teorema do ponto minimizador projetado

Sejam  D \subset \mathbb{R}^n um conjunto convexo e fechado, e  f: \mathbb{R}^n \rightarrow \mathbb{R} um função diferenciável no ponto  \bar{x} \in D .  \bar{x} \in M(f,D \cap B_\epsilon(\bar{x}))

[editar] Mostrar que  P_D(\bar{x}-\alpha f'(\bar{x}))=\bar{x}, \forall \alpha \in \mathbb{R}_+

Funções convexas


[editar] Convexidade da soma de funções convexas

Sejam  D \subset \mathbb{R}^n um conjunto convexo e  f_i: D \rightarrow \mathbb{R}, i=1,...,p, funções convexas em D.  \forall \mu_i \in \mathbb{R}_+, i=1,...,p,

[editar] Mostrar que a função  f: D \rightarrow \mathbb{R}, f(x)=\sum_{i=1}^{p}\mu_i f_i(x) é convexa em D

[editar] Corolário de convexidade do supremo de funções convexas

Sejam  D \subset \mathbb{R}^n um conjunto convexo e  f_i: D \rightarrow \mathbb{R}, i \in I_p, funções convexas em D.  \forall \beta \in \mathbb{R}, f_i(x)\le \beta, \forall x \in D \and i \in I_n

[editar] Mostrar que a função  f: D \rightarrow \mathbb{R}, f(x)=\sup_{i\in I_n} f_i(x) é convexa em D

[editar] Corolário: Função composta de duas convexas é convexa

Sejam  g: \mathbb{R}^n \rightarrow \mathbb{R} uma função convexa e  \psi: \mathbb{R} \rightarrow \mathbb{R} um função convexa e nãodecrescente.

[editar] Mostrar que  f(x)=\psi(g(x)) \; é convexa

[editar] Corolário: Convexidade de conjunto de nível de funções convexas

Suponhamos que o conjunto  D \subset \mathbb{R}^n seja convexo e a função  f: D \rightarrow \mathbb{R} seja convexa em D.

[editar] Mostrar que  L_{f,D}(c) \; é convexo para todo  c \in \mathbb{R}

[editar] Uma função é convexa se os vetores coordenadas são funções convexas

Definição

Seja  D \subset \mathbb{R}^n um conjunto convexo. Se todas as funções  g_i:D \subset \mathbb{R}, i \in I_m são convexas em D, então  g:D \subset \mathbb{R}^m é convexa em D


Introdução


Este material está sendo elaborado com base nas notas de aula da disciplina Otimização II, do curso de Matemática Industrial oferecido pela UFPR, ministrada pelo professor Wilfredo, no segundo semestre letivo do ano de 2008.

O conteúdo do livro não precisa (nem deve) se limitar àquele que consta atualmente no índice. Sendo assim, a qualquer momento o livro pode ser revisto e ampliado.

Sinta-se a vontade para ler este ou quaisquer outros livros do projeto, melhorando-os conforme lhe for possível. Com isso estará ajudando a aumentar a quantidade e a qualidade dos textos didáticos disponíveis em língua portuguesa, ao mesmo tempo em que colaborará com o crescimento projeto Wikilivros como um todo.

Se tiver dúvidas ou sugestões sobre páginas específicas, utilize as páginas de discussão correspondentes para deixar um comentário a respeito.

Ainda há muito por fazer, mas cada um daqueles que contribuem acredita estar fazendo o possível para oferecer o melhor a todos.

Método de gradientes conjugados


Wikipedia
A Wikipédia tem mais sobre este assunto:
Método do gradiente conjugado

[editar] Algumas considerações históricas

  • Este método foi originalmente proposto por Hestenes e Stiefel, em 1952.
  • Seu objetivo inicial foi a resolução de problemas quadráticos sem restrições, mas logo o mesmo foi estendido para casos mais gerais.

[editar] O método

Este método pode ser considerado sob dois pontos de vista:

  • Como um método de descida, com busca linear exata;
  • Como um método de resolução de sistema linear, baseado em um processo de ortogonalização.
Definição

Um conjunto não vazio D \subset \mathbb{R}^n é dito convexo quando \forall x,y \in D e t \in \left[0,1\right] vale

tx+(1-t)y \in D

[editar] Exemplos de conjuntos convexos e côncavos


Definição

Uma função f: D \mapsto \mathbb{R} é dita convexa quando D \subset \mathbb{R}^n é convexo e \forall x,y \in D e \forall t \in \left[0,1\right] vale

f(tx+(1-t)y) \le tf(x) + (1-t) f(y)
Definição

Dado um conjunto convexo D \subset \mathbb{R}^n, uma função f: D \mapsto \mathbb{R} é dita fortemente convexa quando existe uma constante a>0 tal que f(x)-a\|x\| ^2 é convexa.

Exercício

Verifique que uma função quadrática f: \mathbb{R}^n \mapsto \mathbb{R} é fortemente convexa se existe uma matriz simétrica definida positiva A, um vetor a \in \mathbb{R}^n e um escalar \alpha \in \mathbb{R}\ de modo que f(x)=\frac{1}{2}x^\top A x + a^\top x + \alpha.

Resolução
Sendo uma função quadrática, tem-se f(x)=\frac{1}{2}x^\top A x + a^\top x + \alpha. A matriz A pode ser suposta simétrica, pois caso não seja, toma-se B = \frac{A+A^\top}{2} (simétrica), e segue x^\top B x = x^\top A x (verifique).

Além disso, se f é uma função fortemente convexa, então é estritamente convexa. Como f é duas vezes diferenciável (por ser uma função quadrática), a convexidade estrita implica que \nabla^2 f = A é definida positiva.


Nota: Uma matriz é definida positiva se, e somente se, todos os seus autovalores são positivos.

Tem-se:

\nabla f  :  \mathbb{R}^n \mapsto \mathbb{R}
\nabla^2 f  :  \mathbb{R}^n \mapsto (\mathbb{R}^n)^2

Sendo f(x)=\frac{1}{2}x^\top A x + a^\top x + \alpha, segue em particular que \nabla f = A x + a e \nabla^2 f = A = P^\top\Lambda P, onde \Lambda é uma matriz diagonal cujos elementos da diagonal são os autovalores de A e P é uma matriz onde as colunas são os autovetores correspondentes aos autovalores.

Note que A é uma matriz simétrica, pois é a matriz Hessiana de uma função com segundas derivadas parciais contínuas, e consequentemente vale \frac{\partial^2 f}{\partial x \partial y} = \frac{\partial^2 f}{\partial y \partial x} .

Para introduzir o método de direções conjugadas, serão consideradas somente funções quadráticas.

Uma condição necessária de primeira ordem para que x seja um ponto de mínimo para a função f é que \nabla f(x) = 0. Para o presente caso, a função f é convexa, então, a condição necessária \nabla f(x) = 0 também é suficiente.

Exercício

Prove que se A é uma matriz simétrica definida positiva, então f: \mathbb{R}^n \mapsto \mathbb{R} dada por f(x)=x^\top A x + a^\top x + \alpha possui um único ponto de mínimo.

Resolução
Uma vez que A é simétrica definida positiva, a função f é fortemente convexa. Mas toda função fortemente convexa, definida em um conjunto fechado não vazio possui um único minimizador, pois:
  • Os conjuntos de nível de uma função fortemente convexa são compactos;
  • Toda função contínua definida em um compacto tem algum minimizador (pelo teorema de Weierstrass);
  • Os minimizadores de uma função convexa são globais;
  • Funções fortemente convexas não podem ter mais de um minimizador.

Em particular, f possui um único ponto de mínimo.

No caso de uma função quadrática, tem-se \nabla f(x) = 0 \Leftrightarrow A x + a = 0, ou seja, x é solução do sistema linear Ax = -a.

A resolução de um sistema linear nem sempre pode ser feita numericamente de forma eficiente. Por exemplo, se a matriz do sistema é:

A = \begin{bmatrix} 10^{-20} & 1 \\ 1 & 10^{20}+1 \end{bmatrix}

A solução do sistema linear corresponde à interseção entre duas retas quase paralelas, e os erros de truncamento podem causar imprecisão na solução obtida computacionalmente.

Analiticamente, o sistema Ax = -a tem x = -A^{-1}a como solução. Então alguém poderia se perguntar: qual o problema em resolver esse sistema linear, se basta calcular a inversa da matriz A e multiplicar pelo vetor -a? A resposta é que o calculo da inversa de uma matriz em geral é impraticável computacionalmente, por ter custo muito alto. Por isso, nas situações práticas, onde as matrizes tem ordem bem maior do que 2 (digamos 1000), o cálculo de matrizes inversas não é uma opção.

Assim, com o intuito de desenvolver um método computacional para o cálculo de minimizadores, é preciso utilizar outras técnicas. Considere o seguinte:

Em um método de descida tem-se sempre uma sequencia \{ x_k, t_k, d_k \}\in \mathbb{N}, com x_{k+1} = x_k + t_kd_k e t_k é um minimizador de {f(x_k + td_k): t\in \mathbb{R}}

\nabla f(x) = Ax + a

e

0 = \nabla f(x_{k+1}) = A x_{k+1} + a = A(x_k + t_k d_k) + a = Ax_k + t_k A d_k + a

Logo, t_k A d_k = -(A x_k + a) e multiplicando por d_k^\top obtem-se t_k d_k^\top A d_k = -(d_k^\top A x_k + d_k^\top a). Consequentemente, o valor de t_k é dado por

t_k= \frac{-(d_k^\top A x_k + a^\top d_k)}{ d_k^\top A d_k }

Deste modo, o método consistirá de escolher em cada etapa k uma direção d_k, e calcular o coeficiente t_k pela fórmula anterior, para gerar o próximo ponto x_{k+1}. Mas como escolher a direção d_k?

Dado x_k e escolhido d_k, defina \theta: R \mapsto \mathbb{R} como \theta (t) = f(x_k + t d_k), ou seja, \theta é a restrição da função f à reta que passa pelo ponto x_k e que tem direção d_k. Logo, derivando a expressão de \theta em relação a t, obtem-se

\theta^\prime(t) = \nabla f(x_k + t d_k)^\top d_k

Então, no ponto de mínimo, x_{k+1}, tem-se

0 = \nabla f(x_{k+1})^\top d_k

Ou seja, a direção d_k a ser seguida a partir do ponto x_k é ortogonal ao gradiente da função f, no ponto x_{k+1}.

[editar] Esquema do método de descida

x_{k+1} = x_k + t_k d_k = (x_{k-1} + t_{k-1} d_{k-1}) + t_k d_k = \ldots = x_1 + t_1 d_1 + \ldots + t_k d_k = x_1 + \sum_{i = 1}^{k}t_i d_i

Seja \bar{x} o minimizador da função f. Tem-se

x_{k+1} - \bar{x} = x_1 - \bar{x} + \sum_{i = 1}^{k}t_i d_i

Mas 0 = \nabla f(\bar{x}) = A\bar{x} + a implica que a = -A\bar{x}, logo

\nabla f(x) = Ax + a = Ax - A\bar{x} = A(x - \bar{x})

e consequentemente

A(x_{k+1} - \bar{x}) = A(x_1 - \bar{x}) + \sum_{i = 1}^{k}t_i A d_i

Donde \nabla f(x_{k+1}) = \nabla f(x_{1}) + \sum_{i = 1}^{k}t_i A d_i. Portanto 0 = \nabla f(x_{k+1})^\top d_k = \nabla f(x_{1})^\top d_k + \sum_{i = 1}^{k}t_i d_i^\top A d_k.

Exercício

Provar que se A é uma matriz simétrica, definida positiva, então existe uma matriz simétrica B, de modo que A = B^2

Resolução
Sendo uma matriz simétrica, tem-se A = P^\top \Lambda P, com P unitária e
\Lambda =\begin{bmatrix}
\lambda_1 &        & 0      \\
          & \ddots &        \\
0         &        & \lambda_n
\end{bmatrix} = 
\begin{bmatrix}
\sqrt{\lambda_1} &        & 0      \\
          & \ddots &        \\
0         &        & \sqrt{\lambda_n}
\end{bmatrix}^2 = \Lambda_1^2

Logo A = P^\top\Lambda_1\Lambda_1 P = (P^\top\Lambda_1 P) (P^\top \Lambda_1 P) = B^2

Usando o resultado desse exercício, tem-se ainda que 0 = \nabla f(x_{k+1})^\top d_k = \nabla f(x_{1})^\top d_k + \sum_{i = 1}^{k}t_i (B d_i)^\top (B d_k)

Fazendo \delta = Bd, o método do gradiente conjugado escolhe as direções de descida tais que \delta_i^\top d_j = 0, \forall i \not = j. Mas quando i \not = j, tem-se na expressão apresentada anteriormente apenas 0 = \nabla f(x_1)^\top d_k + t_k (B d_k)^\top (B d_k) = \nabla f(x_1)^\top d_k + t_k d_k A d_k

Finalmente, tem-se o algoritmo para este método.

[editar] Algoritmo de Hestenes-Stiefel

Uma comparação da convergência do método de descida do gradiente com tamanho de passo ótimo (em verde) e o método do gradiente conjugado (em vermelho) para a minimização da forma quadrática com um sistema linear dado. O gradiente conjugado, assumindo aritmética exata, converge em no máximo n passos onde n é o tamanho da matriz do sistema (no exemplo, n=2).
Primeiro passo: Escolha x_0 \in \mathbb{R}^n
  Se \nabla f(x_0) = 0, então pare: \bar{x} = x_0
  Senão: d_0 = -\nabla f(x_0) = -A x_0 - a
  Calcular t_0 = \frac{\|\nabla f(x_0) \|^2}{d_0^\top A d_0}
  x_1 = x_0 + t_0 d_0


Passo iterativo k: Dado x_k \in \mathbb{R}^n
  Se \nabla f(x_k) = 0, então pare: \bar{x} = x_k
  Senão: d_k = -\nabla f(x_k) + \frac{\nabla f(x_k)^\top A d_k}{d_k^\top A d_k} d_k
  t_k = \frac{\|\nabla^2 f(x_k) \|^2}{d_k^\top A d_k}
  x_{k+1} = x_k + t_k d_k

Pode-se verificar facilmente que d_{k+1} \perp d_k. De fato, como d_{k+1} = -\nabla f(x_{k+1}) + \frac{\nabla f(x_{k+1})^\top A d_k}{d_k^\top Ad_k}d_k, tem-se Ad_{k+1} = -A\nabla f(x_{k+1}) + \frac{\nabla f(x_{k+1})^\top A d_k}{d_k^\top Ad_k}Ad_k. Logo, d_k^\top Ad_{k+1} = -\nabla f(x_{k+1})^\top A d_k + \frac{\nabla f(x_{k+1})^\top A d_k}{d_k^\top Ad_k} d_k^\top Ad_k = -\nabla f(x_{k+1})^\top A d_k + \nabla f(x_{k+1})^\top A d_k = 0.

Exercício

Provar que se y = Bx então \|y_1 - \bar{y}\|^2 = \|y_{k+1} - \bar{y}\|^2 + \sum_{i=1}^{k}t_i^2 \|\delta_i\|^2.

Resolução
Tem-se
x_{k+1} = x_1 + \sum_{i = 1}^{k}t_i d_i

Multiplicando ambos os membros por B, e trocando x_1 de lugar com x_{k+1} resulta:

-Bx_1 = -Bx_{k+1} + \sum_{i = 1}^{k}t_i Bd_i,

ou seja,

-y_1 = -y_{k+1} + \sum_{i = 1}^{k}t_i \delta_i,

somando \bar{y} em ambos os lados, segue que

\bar{y}-y_1 = \bar{y}-y_{k+1} + \sum_{i = 1}^{k}t_i \delta_i,

Então

\begin{array}{lcl}
\|\bar{y}-y_1 \|^2 & = & (\bar{y}-y_1)^\top(\bar{y}-y_1) \\
                 \ & = & \left(\bar{y}-y_{k+1} + \sum_{i = 1}^{k}t_i \delta_i\right)^\top\left(\bar{y}-y_{k+1} + \sum_{i = 1}^{k}t_i \delta_i\right)\\
                 \ & = & \|\bar{y}-y_{k+1}\|^2 + \sum_{i = 1}^{k}t_i^2 \|\delta_i\|^2 + 2\sum_{i \not= j}^{k}t_it_j \delta_i^\top\delta_j \\
                 \ & = & \|\bar{y}-y_{k+1}\|^2 + \sum_{i = 1}^{k}t_i^2 \|\delta_i\|^2
\end{array}

Sendo a última igualdade devida ao fato de \delta_i^\top\delta_j=0 para i \not = j.



[editar] Exemplos

Considere f: \mathbb{R}^2 \mapsto \mathbb{R} definida por f(x,y) = \frac{1}{2}(x^2 + y^2) = \frac{1}{2} \begin{bmatrix} x & y\end{bmatrix}\begin{bmatrix} 1 & 0 \\ 0 & 1\end{bmatrix}\begin{bmatrix} x\\ y\end{bmatrix}. Em outros termos, tomando u = \begin{bmatrix} x \\ y\end{bmatrix}, tem-se f(u)=\frac{1}{2}u^\top A u, onde A = \begin{bmatrix} 1 & 0 \\ 0 & 1\end{bmatrix} = I_{2\times 2}.

Pode-se aplicar o método de direções conjugadas ao seguinte problema

(P) \left\{\begin{matrix}
min f(u)\\
u \in \mathbb{R}^2
\end{matrix}\right.

Note, desde já, que o conjunto solução é S = \{ \begin{bmatrix} 0 \\ 0\end{bmatrix}\}.

Inicio
  • Toma-se x_0 arbitrário, por exemplo, x_0 = \begin{bmatrix} 2 \\ 1\end{bmatrix}.
  • Avalia-se o gradiente da função f neste ponto inicial: \nabla f (x_0) = A x_0 = I_{2\times 2} x_0 = x_0
Iteração 1
  • d_0 = -\nabla f(x_0) = \begin{bmatrix} -2 \\ -1\end{bmatrix}
  • t_0 = \frac{\| \nabla f(x_0) \|^2}{d_0^\top A d_0} = \frac{5}{5} = 1
  • x_1 = x_0 + t_0 d_0 = \begin{bmatrix} 2 \\ 1\end{bmatrix} + 1 \begin{bmatrix} -2 \\ -1\end{bmatrix} = \begin{bmatrix} 0 \\ 0\end{bmatrix}

A seguir, verifica-se se o gradiente se anula no novo ponto x_1:

  • \nabla f(x_1) = A \begin{bmatrix} 0 \\ 0\end{bmatrix} = \begin{bmatrix} 0 \\ 0\end{bmatrix}

Como o gradiente já é nulo, não é preciso fazer a segunda iteração, e o ponto x_1 é o (único) minimizador global de f.


Em um caso mais geral, considerando f: \mathbb{R}^2 \mapsto \mathbb{R} definida por f(x,y) = \frac{\lambda}{2}(x^2 + y^2) = \frac{1}{2} \begin{bmatrix} x & y\end{bmatrix}\begin{bmatrix} \lambda & 0 \\ 0 & \lambda\end{bmatrix}\begin{bmatrix} x\\ y\end{bmatrix}, tem-se cálculos muito parecidos em cada passo.

O conjunto solução continua sendo S = \{ \begin{bmatrix} 0 \\ 0\end{bmatrix}\}.

Inicio
  • Considere x_0 como no primeiro exemplo, ou seja, x_0 = \begin{bmatrix} 2 \\ 1\end{bmatrix}.
  • Avalia-se o gradiente da função f neste ponto inicial: \nabla f (x_0) = A x_0 = \lambda x_0
Iteração 1
  • d_0 = -\nabla f(x_0) = \lambda \begin{bmatrix} -2 \\ -1\end{bmatrix}
  • t_0 = \frac{\| \nabla f(x_0) \|^2}{d_0^\top A d_0} = \frac{5 \lambda^2}{5 \lambda^3} = \frac{1}{\lambda}
  • x_1 = x_0 + \frac{1}{\lambda} \lambda d_0 = \begin{bmatrix} 2\lambda \\ \lambda\end{bmatrix} + 1 \begin{bmatrix} -2\lambda \\ -\lambda\end{bmatrix} = \begin{bmatrix} 0 \\ 0\end{bmatrix}

A seguir, verifica-se se o gradiente se anula no novo ponto x_1:

  • \nabla f(x_1) = A \begin{bmatrix} 0 \\ 0\end{bmatrix} = \lambda\begin{bmatrix} 0 \\ 0\end{bmatrix} = \begin{bmatrix} 0 \\ 0\end{bmatrix}

Novamente, o gradiente se anula já na primeira iteração, de modo que x_1 é o minimizador global de f.


Exercício

Seja A \in \mathbb{R}^{n\times n} uma matriz simétrica definida positiva, cujos autovalores são todos iguais. Então começando de qualquer ponto x_0 \not = 0, o método fornece x_{n-1} como solução.

Um terceiro exemplo pode ser dado tomando A = \begin{bmatrix} 2 & -1 \\ -1 & 2\end{bmatrix} e f: \mathbb{R}^2 \mapsto \mathbb{R} definida por f(u)=\frac{1}{2}u^\top A u. Observe que tal matriz é simétrica e definida positiva:

\det (A - \lambda I) = (2-\lambda)(3-\lambda) - 1 = \lambda^2 - 4\lambda - 3 = (\lambda-3)(\lambda-1)

Logo, os autovalores de A são \lambda = 1 e \lambda = 3. Isso também implica que a função é fortemente convexa.

Aplicando o método:

Início
  • Toma-se um ponto arbitrário no plano, por exemplo x_0 = \begin{bmatrix} 10 \\ 20\end{bmatrix};
  • Verifica-se se tal ponto é o minimizador global, avaliando nele o gradiente da função:
\nabla f(x_0) = \begin{bmatrix} 2 & -1 \\ -1 & 2\end{bmatrix} \begin{bmatrix} 10 \\ 20\end{bmatrix} = \begin{bmatrix} 0 \\ 30\end{bmatrix} \not = \begin{bmatrix} 0 \\ 0\end{bmatrix}.
  • Já que o gradiente não se anulou no chute inicial, é preciso escolher uma direção e um comprimento de passo para determinar a próxima aproximação:
Iteração 1
d_0 = - \nabla f(x_0) = \begin{bmatrix} 0 \\ -30\end{bmatrix}
t_0 = \frac{\| \begin{bmatrix} 0 & 30\end{bmatrix} \|^2}{\begin{bmatrix} 0 & -30\end{bmatrix} \begin{bmatrix} 2 & -1 \\ -1 & 2\end{bmatrix} \begin{bmatrix} 0 \\ -30\end{bmatrix}} = \frac{900}{\begin{bmatrix} 0 & -30\end{bmatrix} \begin{bmatrix} -30 \\ -60\end{bmatrix}} = \frac{900}{1800} = \frac{1}{2}

Feitos esses cálculos, o próximo ponto é dado por

x_1 = x_0 + t_0 d_0 = \begin{bmatrix} 10 \\ 20\end{bmatrix} + \frac{1}{2}\begin{bmatrix} 0 \\ -30\end{bmatrix} = \begin{bmatrix} 10 \\ 20\end{bmatrix} + \begin{bmatrix} 0 \\ -15\end{bmatrix} = \begin{bmatrix} 10 \\ 5\end{bmatrix}

Para saber se será necessária uma nova iteração, ou se o minimizador foi encontrado, calcula-se o gradiente da função no ponto:

\nabla f(x_1) = \begin{bmatrix} 2 & -1 \\ -1 & 2\end{bmatrix} \begin{bmatrix} 10 \\ 5\end{bmatrix} = \begin{bmatrix} 15 \\ 0\end{bmatrix} \not = \begin{bmatrix} 0 \\ 0\end{bmatrix}.

Novamente, será preciso calcular uma nova direção e um novo comprimento de passo:

Iteração 2
d_0 = \begin{bmatrix} -15 \\ 0\end{bmatrix} + \beta \begin{bmatrix} 0 \\ -30\end{bmatrix} = \begin{bmatrix} -15 \\ -30\beta\end{bmatrix}

onde \beta, no algoritmo de Hestenes é dado por:

\beta
= \frac{\begin{bmatrix} 15 & 0\end{bmatrix} \begin{bmatrix} 2 & -1 \\ -1 & 2\end{bmatrix} \begin{bmatrix} 0 \\ -30\end{bmatrix}}{\begin{bmatrix} 0 & -30\end{bmatrix} \begin{bmatrix} 2 & -1 \\ -1 & 2\end{bmatrix} \begin{bmatrix} 0 \\ -30\end{bmatrix}}
= \frac{\begin{bmatrix} 15 & 0\end{bmatrix} \begin{bmatrix} 30 \\ -60\end{bmatrix}}{\begin{bmatrix} 0 & -30\end{bmatrix} \begin{bmatrix} 30 \\ -60\end{bmatrix}} = \frac{15\times 30}{(-30) \times (-60)} = \frac{1}{4}

Portanto

d_0 = -15 \begin{bmatrix} 1 \\ 1/2\end{bmatrix}

Além disso, o tamanho do passo é dado por

t_1
= \frac{\| \nabla f(x_1) \|^2}{d_0^\top A d_0}
= \frac{15^2}{15^2 \begin{bmatrix} 1 & 1/2\end{bmatrix} \begin{bmatrix} 2 & -1 \\ -1 & 2\end{bmatrix} \begin{bmatrix} 1 \\ 1/2\end{bmatrix}}
= \frac{1}{\begin{bmatrix} 1 & 1/2\end{bmatrix} \begin{bmatrix} 3/2 \\ 0\end{bmatrix}}
= \frac{1}{3/2}
= \frac{2}{3}

Portanto

x_2 = x_1 + t_1 d_1
= \begin{bmatrix} 10 \\ 5\end{bmatrix} - 15 \frac{2}{3}\begin{bmatrix} 1 \\ 1/2\end{bmatrix}
= \begin{bmatrix} 10 \\ 5\end{bmatrix} - 10 \begin{bmatrix} 1 \\ 1/2\end{bmatrix}
= \begin{bmatrix} 0 \\ 0\end{bmatrix}

Obviamente, este é o minimizador procurado (pois o método tem a propriedade de convergência quadrática, ou seja utiliza no máximo n iterações para chegar a solução quando aplicado a funções quadráticas definidas em \mathbb{R}^n)

Exercício

Implementar o algoritmo de Hestenes-Stiefel em alguma linguagem de programação, por exemplo em Scilab, ou Matlab.

Exercício

Seja  f um função quadrática fortemente convexa. Verifique as seguintes igualdades:

  •  \frac{-\nabla f(x_k)^\top d_k}{d_k^\top A d_k} = \frac{\|\nabla f(x_k) \|^2}{d_k^\top A d_k}
  •  \frac{\nabla f(x_k)^\top A d_{k-1}}{d_{k-1}^\top A d_{k-1}} = \frac{\|\nabla f(x_k) \|^2}{\|\nabla f(x_{k-1}) \|^2}
  •  \frac{\nabla f(x_k)^\top A d_k}{d_{k-1}^\top A d_{k-1}} = \frac{\nabla f(x_k)^\top (\nabla f(x_k) - \nabla f(x_{k-1}))}{\|\nabla f(x_{k-1})\|^2}

[editar] Implementação em Scilab

Abaixo é apresentada uma implementação deste algoritmo na linguagem de programação utilizada pelo Scilab.

A = [2 1; 1 2];
 
function [x] = min_gradiente_conjugado(xk)
  //Entrada: xk em R^n, qualquer "chute inicial"
  //  Saída: x, o ponto em que f assume o valor mínimo
 
  k        = 0        //Indica quantos loops já foram feitos
  epsilon  = 0.01
  n        = size(xk,1)
  g        = df(xk)
  seq      = zeros(n,n+1)
  seq(:,1) = xk
  while (norm(g) > epsilon) & (k <= n)
    if (k == 0)
      d = -g
    else
      d = Hestenes(g,d,A)
    end
    t  = busca_linear(g,d,A)
    xk = xk + t*d
    k  = k+1
    seq(:,k+1) = xk
    g  = df(xk)
  end
  plot(seq(1,:),seq(2,:))
  x = xk  
endfunction
 
 
function [fu] = f(u)
  fu=(1/2)*(u'*A*u)
endfunction
 
 
function [grf] = df(u)
  grf = A*u
endfunction
 
 
function [d] = Hestenes(g,d,A)
  d=-g + ((g'*A*d)/(d'*A*d))*d
endfunction
 
 
function [t] = busca_linear(g,d,A)
  t=(g'*g)/(d'*A*d)
endfunction

Para facilitar a compreensão do método, pode ser útil exibir as curvas de nível da função. Uma forma de implementar uma função com esse propósito é a seguinte:

function plotar_curvas
  qtd=101
  tam=max(seq)
  x=linspace(-tam,tam,qtd)
  y=x
  z=zeros(qtd,qtd)
  for i=1:qtd
    for j=1:qtd
      z(i,j)=f([x(i);y(j)])
    end
  end
  contour2d(x,y,z,10)
  a=gca();
  a.x_location = "middle"; 
  a.y_location = "middle"; 
endfunction

[editar] Algoritmo de Fletcher-Reeves

Crystal Clear app kaddressbook.png Um dos autores deste material sugeriu a adição de uma imagem para ilustrar o método de Fletcher-Reeves.


Esta versão é na verdade uma extensão do algoritmo anterior, permitindo a aplicação no caso de funções que não são quadráticas.

Primeiro passo: Escolha  x_0 \in \mathbb{R}^n
 Se  \nabla f(x_0) = 0, então pare:  \bar{x} = x_0
 Senão:  d_0 = -\nabla f(x_0) (como em todo método de descida)
 Calcular  t_0, através de uma busca linear
  x_1 = x_0 + t_0 d_0
Passo iterativo:
 Se  \nabla f(x_k) = 0, então pare:  \bar{x} = x_k
 Senão:  d_k = -\nabla f(x_k) + \frac{\|\nabla f(x_k)\|^2}{\|\nabla f(x_{k-1})\|^2}d_{k-1}
 Calcular  t_k, através de uma busca linear
  x_{k+1} = x_k + t_k d_k
  k = k+1


Crystal Clear app kaddressbook.png Um dos autores deste material sugeriu a adição de uma implementação do algoritmo acima em SciLab.

[editar] Algoritmo de Polak-Ribière

Crystal Clear app kaddressbook.png Um dos autores deste material sugeriu a adição de uma imagem para ilustrar o método de Fletcher-Reeves.

Uma outra versão é a seguinte:

Primeiro passo: Tomar  x_0 \in \mathbb{R}^n
 Se  \nabla f(x_0) = 0, então pare:  \bar{x} = x_0
 Senão:  d_0 = -\nabla f(x_0) (como em todo método de descida)
 Calcular  t_0, através de uma busca linear
  x_1 = x_0 + t_0 d_0
  k = 1
Passo iterativo:
 Se  \nabla f(x_k) = 0, então pare:  \bar{x} = x_k
 Senão:  d_k = -\nabla f(x_k) + \frac{\nabla f(x_k)^\top (\nabla f(x_k) - \nabla f(x_{k-1}))}{\|\nabla f(x_{k-1})\|^2}d_{k-1}
 Calcular  t_k, através de uma busca linear
  x_{k+1} = x_k + t_k d_k
  k = k+1


Crystal Clear app kaddressbook.png Um dos autores deste material sugeriu a adição de uma implementação do algoritmo acima em SciLab.
Exercício

Verificar que, no caso de uma função  f quadrática e fortemente convexa, os algoritmos de Hestenes-Stiefel, de Fletcher-Reeves e de Polak-Ribière são os mesmos.

Exercício

Seja f(x) = e^{-x} + e^x. Implemente o método de gradientes conjugados, e utilize o algoritmo para determinar o ponto de mínimo da função f. Note que o espaço é unidimensional, então o método de gradientes conjugados reduz-se ao método dos gradientes, com primeira direção -\nabla f(x_0). Observe ainda que f é uma função coerciva fortemente convexa.

[editar] Algoritmo auxiliar

Para o caso de funções não quadráticas, é preciso usar algum método de busca linear para a implementação do método dos gradientes conjugados, seja a versão de Fletcher-Reeves ou a de Polak-Ribière. Uma possibilidade é a busca de linear de Armijo (ver Izmailov & Solodov (2007), vol 2, pag. 65), cujo algoritmo é esboçado a seguir:

function busca_linear_Armijo (f, theta, alpha, delta, t0)
  while (alpha * pred > ared)
    t = d * t
  end
endfunction

com:

  • pred = -t\theta
  • \theta(t) = f(x+td)
  • \theta'(t) = \nabla f(x+td)^\top d


Crystal Clear app kaddressbook.png Implementar a regra de Armijo e corrigir o esboço acima.

Métodos de penalidades

Wikipedia
A Wikipédia tem mais sobre este assunto:
Métodos de penalidades

Os métodos que recebem este nome fazem parte de uma família de métodos baseados em:

  • Simplicidade conceitual;
  • Eficiência prática;

Algumas das primeiras funções de penalidade foram desenvolvidas por:

  • Courant (1943)
  • Ablate Brighham (1955) (qual o artigo?)
  • AV Fiacco, GP McCormick (1968) (qual o artigo?)

Para compreender este tipo de método, será considerado o seguinte problema:

(P) \left\{\begin{matrix}
min f(x)\\
g_i(x)\le 0 \forall i = 1, \ldots, p\\
h_j(x)  = 0 \forall j = 1, \ldots, q
\end{matrix}\right.

Adicionalmente, se for considerado o conjunto de pontos C = \{ x\in \mathbb{R}^n: g_i(x)\le 0, \forall i = 1, \ldots, p; h_j(x)  = 0, \forall j = 1, \ldots, q\}, o problema (P) se escreve ainda como

(P) \left\{\begin{matrix}
min f(x)\\
x \in C
\end{matrix}\right.

Uma primeira idéia para a resolução desse tipo de problema é fazer uso de funções indicadoras, conforme é definido a seguir:

Definição

Dado um conjunto não vazio C \subset \mathbb{R}^n, define-se a função indicadora I_C:\mathbb{R}^n \mapsto \mathbb{R} \cup \{ \infty \} por:

I_C(x) = \begin{cases}
     0, & \mbox{se }x    \in C\\
\infty, & \mbox{se }x\not\in C 
\end{cases}

Notas: A função indicadora I_C é às vezes denotada por \chi_{C}.

Para transformar um conjunto (P) em um problema sem restrições, pode-se proceder da seguinte maneira: Considerar o problema (PD), dado por:

(PD) \left\{\begin{matrix}
min f(x) + I_C(x)\\
x\in \mathbb{R}^n
\end{matrix}\right.

Considerando i_C:\mathbb{R} \mapsto \mathbb{R} \cup \{ \infty \} definida por:

i_C(t) = \begin{cases}
     0, & \mbox{se } t\le 0\\
\infty, & \mbox{se } t  > 0 
\end{cases}

Tem-se ainda o problema (PP), dado por:

(PP) \left\{\begin{matrix}
min f(x) + \sum_{i = 1}^{p} i(g_i(x)) + \sum_{j = 1}^{p} i(h_j(x))\\
x\in \mathbb{R}^n
\end{matrix}\right.

Comentários:

  • A grande vantagem desta idéia é que ela transforma um problema com restrições em um problema irrestrito.
  • A principal desvantagem é que a função \phi definida a seguir é descontínua em cada ponto x \in \partial C:
\phi:\mathbb{R}^n \mapsto \mathbb{R} \cup \{ \infty \}
\phi(x) = f(x) + \sum_{i = 1}^{p} i(g_i(x)) + \sum_{j = 1}^{p} i(h_j(x))

[editar] Método de penalidade exterior

Para contornar a desvantagem da descontinuidade da função apresentada anteriormente, surgem outras funções, como por exemplo:

Gráfico de uma função de penalidade
H:\mathbb{R} \mapsto \mathbb{R}
H(t) = \begin{cases}
  0, & \mbox{se } t\le 0\\
t^2, & \mbox{se } t  > 0 
\end{cases}
Definição

Dada uma função g:A \mapsto \mathbb{R}, definida em um conjunto arbitrário A, define-se a parte positiva de g, como:

g^+(x) = max \{0, g(x) \}


Crystal Clear app kaddressbook.png Um dos autores deste material sugeriu a adição de uma imagem para comparar uma função g(x) e sua correspondente g^+(x) = max \{0, g(x) \}.
Exercício

Verifique que para cada x \in \mathbb{R}^n tem-se, para cada i, a igualdade H(g_i(x)) = \left(g_i^+(x) \right)^2, onde g_i^+ é a parte positiva de g_i e H é a função definida no exemplo anterior.


Resolução
De fato, dada qualquer função g:A \mapsto \mathbb{R}, segue da própria definição de H que
H(g(x)) =
\begin{cases}
     0, & \mbox{se } g(x)\le 0\\
g(x)^2, & \mbox{se } g(x)  > 0
\end{cases}

Mas

g^+(x) = max \{0, g(x) \} =
\begin{cases}
   0, & \mbox{se } g(x)\le 0\\
g(x), & \mbox{se } g(x)  > 0
\end{cases}

implica que

(g^+(x))^2 = max \{0, g(x) \} =
\begin{cases}
       0, & \mbox{se } g(x)\le 0\\
(g(x))^2, & \mbox{se } g(x)  > 0
\end{cases}

Portanto, H(g(x)) = \left(g^+(x) \right)^2. Como g pode ser tomada arbitrariamente, tem-se em particular que H(g_i(x)) = \left(g_i^+(x) \right)^2.

Exercício

Verifique que para cada x \in \mathbb{R}^n tem-se, para cada j, a igualdade H\left(h_j(x)^2\right) = \left(h_j(x) \right)^2.


Crystal Clear app kaddressbook.png Um dos autores deste material sugeriu conferir o enunciado do exercício anterior. A afirmação parece ser falsa nos casos em que h_j(x)<0.

A idéia de aplicar penalizações aos pontos que não pertencem ao conjunto viável é formalizada na seguinte definição:

Definição

Seja H:\mathbb{R}^n \mapsto \mathbb{R}. A função H é chamada de função de penalidade exterior se possui as seguintes propriedades:

  • H(x) \ge 0, \forall x \in \mathbb{R}^n
  • H(x)  =  0, \Leftrightarrow x \in C
  • H é contínua

Nota: Lembre-se que C = \{ x\in \mathbb{R}^n: g_i(x)\le 0, \forall i = 1, \ldots, p; h_j(x)  = 0, \forall j = 1, \ldots, q\} é o conjunto viável do problema (P).

Em particular, as funções H \circ g_i são funções de penalidade exterior.

Definição

Seja \theta:\mathbb{R}^n \mapsto \mathbb{R}. A função \theta é dita coerciva se \lim_{\| x \| \to \infty}\theta(x) = +\infty


Crystal Clear app kaddressbook.png Um dos autores deste material sugeriu a adição de exemplos de funções de penalidade, juntamente com algumas imagens ilustrando os seus gráficos.

Nota: Conforme o Wikcionário, o termo coercivo significa: que coage; que reprime; que impõe pena; coercitivo. Nesse sentido, esse é um termo adequado ao tratar do conceito anterior, no contexto dos métodos de penalidade.

Exercício

Verifique que \theta é uma função coerciva se, e somente se, L_\theta (\lambda) = \{ x \in \mathbb{R}^n : \theta (x) \le \lambda \} é limitado para todo \lambda \in \mathbb{R}.

Resolução
Pela definição de limite, a afirmação
\lim_{\| x \| \to \infty}\theta(x) = +\infty

é equivalente a dizer que

\forall K > 0, \exists M >0 tal que \forall x tem-se \|x\|>M \Rightarrow \theta(x)>K.

Esta última implicação, é equivalente à

\theta(x)\le K \Rightarrow \|x\| \le M (sua contrapositiva).

Pela definição de conjunto de nível, isso equivale à x \in L_\theta(K) \Rightarrow \|x\| \le M. A existência de um número M com tal propriedade significa que L_\theta(K) é um conjunto limitado, donde conclui-se a equivalência entre coercividade e limitação dos conjuntos de níveis


Crystal Clear app kaddressbook.png Um dos autores deste material sugeriu a adição de uma imagem que ilustre geometricamente a relação entre coercividade e conjuntos de nível.
Exercício

Verifique que H:\mathbb{R}^n \mapsto \mathbb{R} definida por H(x) = \sum_{i = 1}^{p}\left( g_i^+(x)\right)^2 + \sum_{j = 1}^{q}\left( h_j(x)\right)^2 é uma função de penalidade exterior, onde conforme anteriormente, g_i^+ é a parte positiva de g_i.


Resolução
Primeiramente, H(x) é uma função não negativa, pois o quadrado de qualquer número real é não negativo, assim como a soma de números reais não negativos.

Em segundo lugar, tem-se H(x) = 0 se, e somente se, para cada índice i e cada índice j vale g_i^+(x) = 0 e h_j(x) = 0. Como g_i^+(x) = 0 \Leftrightarrow g_i \le 0, tem-se H(x)=0 \Leftrightarrow x \in C

Finalmente, como a soma e o produto de funções contínuas resulta em uma função contínua, segue que H(x) é contínua, pois g_i, h_j são contínuas.

Exercício

Verifique que se g: \mathbb{R}^n \mapsto \mathbb{R} é uma função contínua coerciva, então existe \bar{x} \in \mathbb{R}^n tal que g(\bar{x}) \le g(x), \forall x \in \mathbb{R}^n.


Resolução
Considere um ponto arbitrário \hat{x} \in \mathbb{R}^n. Tome \lambda = g(\hat{x}). Nestas condições, L_g(\lambda) é limitado (conforme um dos exercícios anteriores) e fechado (pois é pré-imagem de um conjunto fechado por uma função contínua), portanto compacto. Neste caso, conforme o teorema de Weierstrass, a função g possui algum ponto de mínimo no conjunto L_g(\lambda), ou seja, existe \bar{x} \in L_g(\lambda) tal que g(\bar{x}) \le g(x), \forall x \in L_g(\lambda).

Ne verdade, tal ponto \bar{x} é também um minimizador global da função g, pois se x \not\in L_g(\lambda) então g(x) > \lambda (pela definição do conjunto L_g(\lambda)).


Crystal Clear app kaddressbook.png Um dos autores deste material sugeriu a adição de uma imagem ilustrando a existência de minimizadores globais para funções coercivas.

[editar] Alguns conceitos da topologia

Em algumas situações, é interessante ter em mente que certos conceitos definidos no contexto da Otimização são, na verdade, instanciações de conceitos mais gerais, muitos deles provenientes da topologia. Alguns exemplos são apresentados a seguir.

Definição

Dado um conjunto X, uma coleção \Gamma \subset \mathcal{P} de subconjuntos de X é chamada de topologia se:

  • \emptyset, X \in \Gamma
  • \left\{ A_i \right\}_{i \in I} \subset \Gamma \Rightarrow \cup_{i \in I} A_i \in \Gamma
  • \left\{ B_i \right\}_{i =1}^p \subset \Gamma \Rightarrow \cap_{i =1}^p B_i \in \Gamma
Exemplos
  • Topologia euclidiana:
X = R
\Gamma = \Gamma_1 = \{ \emptyset \} \cup \{X\} \cup \left\{ A \subset \mathbb{R}: \forall x \in A, \exists \epsilon_x > 0, \left( x-\epsilon_x, x+\epsilon_x \right) \subset A \right\}

Em outras palavras, a topologia euclideana é a coleção de todos os conjuntos abertos contidos em \mathbb{R}. Pode-se verificar com facilidade que de fato são satizfeitas as três propriedades que definem uma topologia.

Outro exemplo muito comum é o seguinte:

  • Topologia euclideana estendida:
X = R \cup \{ \infty \}
\Gamma = \Gamma_1 \cup \left\{ A \subset \mathbb{R} \cup \{ \infty \}: \exists a \in \mathbb{R}, A = \left] a, \infty \right]\right\}

Em geral, a noção de limite seria caracterizada topologicamente da seguinte forma:

Definição

\lim_{n \to \infty} x_n = a \Leftrightarrow \forall I \ni x, \exists n_0 \in \mathbb{N}, \forall n \ge n_0 \Rightarrow x_n \in I

[editar] De volta ao método

Uma vez apresentados os conceitos iniciais, pode-se provar o seguinte teorema:

Teorema

Considere:

  • H: \mathbb{R}^n \mapsto \mathbb{R} uma função de penalidade exterior;
  • f: \mathbb{R} \mapsto \mathbb{R} uma função contínua;
  • C = \{ x\in \mathbb{R}^n: g_i(x)\le 0, \forall i = 1, \ldots, p; h_j(x)  = 0, \forall j = 1, \ldots, q\} um conjunto fechado;
  • \left\{ r_k \right\} uma sequência de termos positivos tal que \lim_{k \to \infty} r_k = 0.
Suponha que é válida uma das seguintes propriedades:
  1. f é coerciva.
  2. C é limitado e H é coerciva.
Se, para cada k, for escolhido \bar{x}(r_k) \in \arg \min \{ f(x) + \frac{1}{r_k}H(x)\},então:
  1. \left\{ \bar{x}(r_k) \right\} possui algum ponto de acumulação;
  2. Todo ponto de acumulação de \left\{ \bar{x}(r_k) \right\} é solução do problema (P);
  3. \lim_{k \to \infty} H(\bar{x}(r_k)) = 0.


Crystal Clear app kaddressbook.png Um dos autores deste material sugeriu a adição de uma imagem que ilustre geometricamente o significado do teorema acima.


Demonstração
Se (1) acontece, então f + \frac{1}{r}H é coerciva, pois
\lim_{\|x\| \to \infty} f(x) + \frac{1}{r}H(x) \ge \lim_{\|x\| \to \infty} f(x) = \infty

A desigualdade é válida pois H é uma função de penalidade, portanto não negativa, e a igualdade se deve à hipotese sobre f. Além de coerciva, tal função é contínua (pois é combinação linear de funções contínuas). Logo, a função possui um ponto de mínimo global, ou seja, existe \bar{x} \in \mathbb{R}^n tal que

f(\bar{x}) + \frac{1}{r}H(\bar{x}) \le f(x) + \frac{1}{r}H(x)

, para qualquer x \in \mathbb{R}^n.

Mas f é contínua e coerciva, então também existe x^* \in \mathbb{R}^n tal que

f(x^*) \le f(x), \forall x \in C

Por outro lado, se vale (2), então f + \frac{1}{r}H é contínua e C compacto. Analogamente, f é contínua e C compacto, donde tem-se algum x^* \in C tal que

f(x^*) \le f(x), \forall x \in C

Em ambos os casos, dada uma sequência \left\{r_k\right\} tal que 0 < r_{k+1} < r_k e \lim_{k \to \infty} r_k = \infty, defina-se \varphi: \mathbb{R}^n \times \mathbb{R} \mapsto \mathbb{R} por \varphi (x,r) = f(x) + \frac{1}{r}H(x). Logo,

\begin{array}{lcl}
\varphi(\bar{x}(r_k), r_k) &   = & f(\bar{x}(r_k)) + \frac{1}{r_k}H(\bar{x}(r_k)) \\
                         \ & \le & f(\bar{x}(r_{k+1})) + \frac{1}{r_k}H(\bar{x}(r_{k+1}) \\
                         \ & \le & f(\bar{x}(r_{k+1})) + \frac{1}{r_{k+1}}H(\bar{x}(r_{k+1}) \\
                         \ &   = & \varphi(\bar{x}(r_{k+1}), r_{k+1})
\end{array}

Sendo que a primeira desigualdade se deve ao fato de \bar{x}(r_k) ser um minimizador, por construção, e a segunda segue por que \{r_k\} é decrescente. Portanto, \varphi(\bar{x}(r_k), r_k) \le \varphi(\bar{x}(r_{k+1}), r_{k+1}), \forall k.

Além disso, tem-se as seguintes desigualdades:

f(\bar{x}(r_{ k })) + \frac{1}{r_{ k }}H(\bar{x}(r_{ k })) \le f(\bar{x}(r_{k+1})) + \frac{1}{r_{ k }}H(\bar{x}(r_{k+1}))
f(\bar{x}(r_{k+1})) + \frac{1}{r_{k+1}}H(\bar{x}(r_{k+1})) \le f(\bar{x}(r_{ k })) + \frac{1}{r_{k+1}}H(\bar{x}(r_{ k }))

Logo, somando os membros correspondentes, obtem-se:

\frac{1}{r_k}H(\bar{x}(r_k)) + \frac{1}{r_{k+1}}H(\bar{x}(r_{k+1})) \le \frac{1}{r_k}H(\bar{x}(r_{k+1})) + \frac{1}{r_{k+1}}H(\bar{x}(r_{k}))

ou seja,

\left( \frac{1}{r_{k+1}} - \frac{1}{r_{k}}\right) H(\bar{x}(r_{k+1})) \le \left( \frac{1}{r_{k+1}} - \frac{1}{r_{k}}\right) H(\bar{x}(r_{k}))

Portanto,

H(\bar{x}(r_{k+1})) \le H(\bar{x}(r_{k})), \forall k

Por outro lado, f(\bar{x}(r_k)) \le f(\bar{x}(r_k)) + \frac{1}{r_k}H(\bar{x}(r_k)) \le f(x^*) + \frac{1}{r_k}H(x^*) = f(x^*), \forall k, sendo primeira desigualdade válida por H ser não negativa e r_k positivo, e a segunda devida à própria definição de \bar{x}(r_k). Logo, f(\bar{x}(r_k)) \le f(x^*).

Se a primeira das hipóteses acontece, segue da coercividade e dessa última desigualdade que \{\bar{x}(r_k)\} \subset L_f(f(x^*)). Então \{\bar{x}(r_k)\}_{i=1}^{\infty} é uma sequência limitada.

Se ocorre a segunda, então H é coerciva, mas foi mostrado que H(\bar{x}(r_{k+1})) \le H(\bar{x}(r_{k})), consequentemente \{\bar{x}(r_{k+1})\} \subset L_H(H(\bar{x}(r_1))). Sendo H coerciva, conclui-se novamente que \{\bar{x}(r_k)\}_{i=1}^{\infty} é uma sequência limitada.

Portanto, em qualquer caso, \{\bar{x}(r_k)\}_{i=1}^{\infty} possui algum ponto de acumulação.

Seja \bar{x} um ponto de acumulação de \{\bar{x}(r_k)\}_{i=1}^{\infty}. Então existe \{\bar{x}(r_{k_n})\} \subset \{\bar{x}(r_k)\} tal que \lim_{n \to \infty}\bar{x}(r_{k_n}) = \bar{x}. Logicamente, \lim_{n \to \infty}r_{k_n} = 0. Pela continuidade de f e sabendo que f(\bar{x}(r_k)) \le f(x^*), se deduz que f(\bar{x}) = \lim_{n \to \infty}f(\bar{x}(r_{k_n})) \le f(x^*). Mas já havia sido verificado que f(x^*) \le f(x), \forall x \in C, então segue a igualdade f(\bar{x}) = f(x^*).

A sequência \{\varphi(\bar{x}(r_k), r_k)\} é crescente. Seja \bar{\varphi} = \lim_{n \to \infty}\varphi(\bar{x}(r_k),r_k) (por que razão ele existe?). Como \varphi(\bar{x}(r_k),r_k) = f(\bar{x}(r_k)) + \frac{1}{r_k}H(\bar{x}(r_k)), tem-se \bar{\varphi} = \lim_{n \to \infty}\frac{1}{r_k}H(\bar{x}(r_k)) = \lim_{n \to \infty}\varphi(\bar{x}(r_k),r_k) - \lim_{n \to \infty}f(\bar{x}(r_k)) = \bar{\varphi}-f(\bar{x}). Portanto, \lim_{n \to \infty}H(\bar{x}(r_{k_n})) = 0. Logo \lim_{k \to \infty}H(\bar{x}(r_k)), ou seja, \lim_{r \to 0}H(\bar{x}(r)).

Como H é contínua, \lim_{n \to \infty}H(\bar{x}(r_{k_n})) = H(\bar{x}), e este valor é nulo se, e somente se, \bar{x} \in C.

Logo, f(\bar{x}) \le f(x^*) \le f(\bar{x}), donde f(\bar{x}) = f(x^*). Assim, tem-se f(\bar{x}) = f(x^*) \le f(x), \forall x \in C, ou seja, \bar{x} é solução de (P).

Exercício

Dado o problema (P), considere m = \inf \left\{ f(x) : g_i(x) \le 0, \forall i \ \text{e}\ h_j = 0, \forall j\right\} (isso não quer dizer que o problema tenha solução). suponha-se que f, g e h são funções contínuas e que C = \left\{x: g_i(x) \le 0, \forall i ; h_j(x) = 0, \forall j\right\} seja não vazio, ou seja, que C é factível. Tome H: \mathbb{R}^n \mapsto \mathbb{R} como H(x) = \sum_{i=1}^{p}[g_i^+(x)]^2 + \sum_{j=1}^{q}[h_j(x)]^2, onde g_i^+ denota a parte positiva de g_i, como de costume. Considere ainda \varphi: \mathbb{R}^n \times \mathbb{R}^+ \mapsto \mathbb{R} dada por \varphi (x,r) = f(x) + \frac{1}{r}H(x) e, para cada r>0, seja m(r) = \inf_{x \in \mathbb{R}^n} \varphi(x,r) Nessas condições, provar que:

  1. H é uma função de penalidade exterior
  2. m>-\infty
  3. \varphi(x,r) = f(x), \forall x \in C, \forall r>0
  4. Se 0 < r < s então m(s) < m(r)
  5. Se f, g_i, h_j são covexas, então \varphi é convexa.
  6. Se f, g_i, h_j são diferenciáveis, então \varphi é diferenciável em x, e \nabla_x\varphi(x,r) = \nabla f(x) + \sum_{i=1}^{p}g_i^+(x)\nabla g_i(x) + \sum_{j=1}^{q}h_j(x)\nabla h_j(x)
  7. Se 0 < r_{k+1} < r_k e \lim_{k\to\infty} r_k = 0 e se \{x_k\} é uma sequência tal que \varphi(x_k,r_k) = m(r_k) e \lim_{k\to\infty} x_k = \bar{x} então \bar{x} é solução de (P).

[editar] Algoritmo de penalidade exterior

Primeiro passo: Escolha x_0 \in \mathbb{R}^n, r>0 e k=1

Passo iterativo k: Calcular x_k = \bar{x} (r_k), solução global de 
(P) \left\{\begin{matrix}
\min \varphi(x,r_k) & = & f(x) + \frac{1}{2r}H(x)\\
x \in \mathbb{R}^n &   &
\end{matrix}\right.
Escolha r_k tal que 0 < r_k< r_{k-1}

Nota: Neste algoritmo, H é uma função de penalidade exterior.

Exercício

Sejam f(x,y) = x^2 + y^2 e C = \{ (x,y)\in \mathbb{R}^2 : -x \le 0; -y \le 0; x + y - 1 \le 0 \}. Considere o problema: (P) \left\{\begin{matrix}
\min f(x,y)\\
(x,y) \in \mathbb{R}^n
\end{matrix}\right. Utilize o método de penalidade exterior para determinar o ponto de mínimo da função f sobre o conjunto C.

Resolução
Esboço: Pode-se tomar H(x,y) = \sum_{i=1}^{3} (g_i^+ (x,y))^2, onde:
g_1 (x,y) = -x
g_2 (x,y) = -y
g_3 (x,y) = x+y-1

Tem-se então o problema irrestrito: (P_r) \left\{\begin{matrix}
\min f(x,y) + \frac{1}{r}H(x,y)\\
(x,y) \in \mathbb{R}^2
\end{matrix}\right., onde pode ser aplicado o método de gradientes conjugados.

[editar] Implementação do algoritmo de penalidade exterior em SciLab

Considerando o problema

(P) \left\{\begin{matrix}
min f(x)\\
g_i(x)\le 0 \forall i = 1, \ldots, p\\
h_j(x)  = 0 \forall j = 1, \ldots, q
\end{matrix}\right.

com f uma função quadrática, A simétrica positiva definida, g_i, h_j lineares, e a função de penalidade, H, dada por:

H(x) = \sum_{i=1}^{p} (g_i^+ (x,y))^2 + \sum_{j=1}^{q} (h_j (x,y))^2

Pode-se usar o método dos gradientes conjugados para resolver o seguinte problema:

(P) \left\{\begin{matrix}
min f(x) + H(x)\\
x \in \mathbb{R}^n
\end{matrix}\right.

onde f + H terá sua matriz Hessiana definida positiva.


Crystal Clear app kaddressbook.png Incluir o código para uma implementação do método em SciLab.

[editar] Método de penalidade interior

Este método também é conhecido como método de barreira. Ele consiste em trabalhar com funções de penalidade tais que \forall x \in \partial C e qualquer que seja a sequência \{x_n\} \subset C = \{x: g_i(x) < 0\} para a qual x_n \to x, se tem que a função de penalidade tende a +\infty.


Crystal Clear app kaddressbook.png Um dos autores deste material sugeriu a adição de uma imagem para ilustrar uma sequência de pontos que tende a um ponto da fronteira de um conjunto C, de preferência junto com o gráfico de uma função de penalidade deste novo tipo.

Mas como conseguir esse tipo de função?

[editar] Exemplos

Considere o caso em que g_i(x) = a_i^\top x - b_i. Lembrando que a função logarítmica tem a propriedade:

\lim_{x\to 0^+}-\ln(x) = +\infty

pode-se tomar H_1(x) = -\sum_{i=1}^{p}\ln(b_i-a_i^\top x), pois

\lim_{a_i^\top x\to b_i}-\ln(b_i-a_i^\top x) = +\infty

implica que

\lim_{x\to \partial C} H_1(x) = +\infty

Analogamente, tem-se

\lim_{x\to 0^+} \frac{1}{x} = +\infty

então, uma outra função com a propriedade desejada é H_2(x) = \sum_{i=1}^{p}\frac{1}{b_i-a_i^\top x}.

O problema a ser resolvido quando se quer aplicar o método de barreira é o seguinte:

(P) \left\{\begin{matrix}
\min f(x)\\
g_i(x) \le 0; i = 1, \ldots, p
\end{matrix}\right.

onde C = \{ x: g_i(x) \le 0, \forall i\} é tal que C^0 = \{ x: g_i(x) < 0, \forall i\} \not = \emptyset

Observação: C^0 não é necessariamente igual ao interior do conjunto C. Por exemplo:


Crystal Clear app kaddressbook.png Um dos autores deste material sugeriu a adição de uma imagem para ilustrar um conjunto C cujo correspondente C0 seja diferente do interior de C.
Definição

Se diz que uma função B: C^0 \mapsto \mathbb{R} é uma função de barreira se ela possui as seguintes propriedades:

  1. B é limitada inferiormente em C.
  2. Para qualquer sequência \{ x_n\} \subset C^0 tal que \lim_{n \to \infty} x_n = x \in \partial C vale \lim_{n \to \infty} B(x_n) = +\infty
  3. B é contínua

[editar] Algoritmo de barreira

Primeiro passo: Escolha x_0 \in C^0, r>0 e t_0>0

Passo iterativo k: Calcular x_{k+1} = \bar{x} (t_k), solução global de 
(P) \left\{\begin{matrix}
\min f(x) + t_k B(x)\\
x \in C^0
\end{matrix}\right.
Escolha t_{k+1} tal que 0 < t_{k+1}< t_k

Observação: Neste método os termos da sequência \{ x_n\} estão sempre em C^0.

[editar] Implementação do algoritmo de barreira em SciLab

Considerando o problema

\left\{\begin{matrix}
min f(x)\\
g_i(x)\le 0 \forall i = 1, \ldots, p\\
\end{matrix}\right.

pode-se usar o método de barreira no conjunto:

C = \{ x \in \mathbb{R}^n: g_i(x) \le 0 \}

Pode-se usar qualquer das funções de penalidade interior apresentadas, por exemplo:

B(x) = \sum_{i=1}^{p}\frac{1}{g_i(x)}

ou ainda

B(x) = \sum_{i=1}^{p}-\ln(-g_i(x))

Como C^0 = \{ x \in \mathbb{R}^n: g_i(x) < 0 \}, o problema passa a ser:

\left\{\begin{matrix}
min f(x) + B(x)\\
x \in C^0
\end{matrix}\right.

Observação: Note que é necessário conhecer um ponto inicial x_0 \in C^0, para servir de primeira aproximação, ou ponto de partida para o algoritmo.

Métodos de região de confiança


Considere o seguinte problema de programação diferenciável não linear sem restrições:

(P) \left\{\begin{matrix}
\min f(x)\\
x \in \mathbb{R}^n
\end{matrix}\right.

onde f: \mathbb{R}^n \mapsto \mathbb{R} é de classe C^2 (\mathbb{R}^n).

Observação: Como f não é necessariamente convexa, a matriz \nabla^2f(x) pode não ser definida positiva, apesar de ser simétrica. Neste caso, o método de Newton ou suas variantes (direções conjugadas, quase Newton, etc) não servem.

O primeiro método de região de confiança (em inglês, Trust region method), foi introduzido por Powel em 1970 (qual artigo?) mas oficialmente introduzido por Dennis em 1978 (artigo?). Ele consiste no seguinte:

A cada iteração, se constrói um modelo quadrático e uma região de confiança

este princípio pode ser considerado como uma extensão da busca de Armijo unidimensional.

[editar] Breve revisão da busca de Armijo

Para entender a geometria do método de região de confiança, é bom lembrar a geometria da busca de Armijo unidimensional.

Seja d uma direção de descida de uma função f a partir do ponto x. Então \nabla f(x)^\top d < 0. Agora, considerando \theta : \mathbb{R} \mapsto \mathbb{R}, definida por \theta (t) = f(x + td), tem-se:

\theta'(t) = \nabla f(x+td)^\top d. Assim, \theta'(0) = \nabla f(x)^\top d < 0. Se \alpha \in (0,1), então \theta'(0) < \alpha \theta'(0) < 0.
Exercício

Mostrar que sempre existe t>0 tal que \frac{f(x+td) - f(x)}{t} = \frac{\theta(t) - \theta(0)}{t - 0} \le \alpha \theta'(0).

Resolução
A resolução deste exercício é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a na versão online deste material.

A busca de Armijo consiste em tomar t>0.

Se \frac{f(x + td) - f(x)}{t} > \alpha \theta'(0) > \theta'(0)

Mas como \theta'(0) = \lim_{t\to0} \frac{f(x+td)-f(x)}{t}, então existe algum t>0 tal que vale \frac{f(x + td) - f(x)}{t} \le \alpha \theta'(0) . A tal ponto, chama-se ponto de Armijo.

Busca de Armijo.svg

Introduzindo as seguinte notação

m(x+h) = f(x) + \nabla f(x)^\top h (modelo linear para a função f)
pred = m(x) - m(x+td) = -t\theta'(0) (redução predita por este modelo)
ared = f(x) - f(x+td) (redução real no valor da função)

a pergunta é:

Quando um ponto t vai ser ponto de Armijo com essa notação?

Tem-se:

\alpha \cdot pred = -\alpha t \theta'(0) \le f(x) - f(x + td) = ared

Logo, se \alpha \cdot pred \le ared segue que t é um ponto de Armijo.

Observações: Note que a essência da busca linear de Armijo é construir um modelo linear e um intervalo compacto [0,t], sendo t>0 e o ponto inicial da busca e logo procurar o ponto de Armijo em [0,t].

[editar] De volta ao método

O método de região de confiança será uma generalização da busca de Armijo, consistindo da construção de um modelo quadrático e uma região R, chamada de região de confiança, e nessa região calcular o novo iterando.

[editar] Algoritmo da região de confiança

Primeiro passo: Escolha x_0 \in \mathbb{R}^n, r>0, \alpha \in (0,1), \beta \in (0,1) e k=0.

Passo iterativo k: Enquanto \nabla f(x_k) \not = 0, construa o modelo quadrático:
 m_k(h) = \nabla f (x_k)^\top h + \frac{1}{2}h^\top H_k h
Calcule d, solução de
 (P) \left\{\begin{matrix}
\min m_k(h)\\
\|h\| \le r
\end{matrix}\right.
Tome ared = f (x_k) - f(x_k + d) e pred = -m_k(d)
Se \frac{ared}{pred} > \alpha, fazer x_{k+1} = x_k + d
Senão x_{k+1} = x_k, r = \beta r e k = k+1

Comentários: No algoritmo anterior, quando se tem um passo falho, a região de confiança sempre diminui. Seria bom incluir casos bons, onde a região deve crescer.

[editar] Algoritmo da região de confiança melhorado

Primeiro passo: Escolha x_0 \in \mathbb{R}^n, r>0, 0 < \alpha < \gamma < 1, \beta \in (0,1), \beta_2 > 1, \epsilon \in (0,1) e k=0.

Passo iterativo k: Enquanto \|\nabla f(x_k)\| > \epsilon , construa o modelo quadrático:
 m_k(h) = \nabla f (x_k)^\top h + \frac{1}{2}h^\top H_k h
Continue = 1
Enquanto (Continue = 1)
  Calcule d, solução de
   (P) \left\{\begin{matrix}
\min m_k(h)\\
\|h\| \le r
\end{matrix}\right.
  Tome ared = f (x_k) - f(x_k + d) e pred = -m_k(d)
  Se \frac{ared}{pred} > \alpha, fazer x_{k+1} = x_k + d
    Se r > \gamma, r = \beta_2 r; Continue = 0
    k = k+1
  Senão r = \beta_1 r


[editar] O subproblema quadrático

Conforme foi explicado, o método das regiões de confiança constrói um modelo quadrático da forma:

m(h) = g^\top h + \frac{1}{2}h^\top H h

onde, no método, g = \nabla f (x_k) e H = H_k. Em tal modelo, tem-se

  • g um vetor não nulo;
  • h uma matriz simétrica (que pode não ser definida positiva)

O problema quadrático é

(PQ) \left\{\begin{matrix}
\min m(h)\\
\|h\| \le r
\end{matrix}\right.

com r>0.

Este é o problema que será tratado a seguir.

Exercício

Provar que (PQ) sempre tem solução.

Resolução
Esboço: Como h é uma função contínua e o conjunto onde se quer minimizar tal função é uma bola, em dimensão finita, trata-se de minimizar uma função contínua em um compacto. Esse é um problema que tem solução, pelo teorema de Weierstrass.

O exercício anterior garante que o método está bem definido, quer dizer, todas as etapas podem ser realizadas.

O problema de mínimos quadrados


Este tipo de problema é muito frequente em ciências experimentais. Para ter um exemplo em mente durante a discussão que será feita mais adiante, considere as seguintes informações:

Segundo dados disponibilizados pelo IBGE, o Produto Interno Bruto per capita na cidade de São Paulo, no período de 2002 a 2005, foi (em reais): 17734, 19669, 20943 e 24083, respectivamente.

Com base nessas informações, como poderia ser feita uma previsão do valor correspondente ao ano seguinte (2006)?

Uma escolha possível seria supor que a cada ano o PIB aumenta uma aproximadamente constante, ou seja, usar um modelo linear para obter tal estimativa (que possivelmente será bem grosseira). Intuitivamente, bastaria analisar os dados disponíveis e a partir deles deduzir qual é o aumento que ocorre a cada ano. Depois, a previsão para 2006 seria aproximadamente igual à de 2005 somada com aquele aumento anual.

Esta idéia poderia até funcionar para o caso deste exemplo, mas o que fazer se a quantidade de dados disponíveis sobre algum fenômeno (ou alguma situação) for significativamente maior?

A melhor escolha, sem dúvida, é fazer uso de um computador para obter o modelo que melhor descreve o "comportamento" dos dados experimentais.

Em geral, os problemas de mínimos quadrados consistem em identificar os valores de determinados parâmetros, de modo que se satisfaçam certas equações h_i(x) = 0,\, i = 1, \ldots, p. No contexto do exemplo anterior, se procura um modelo linear para os dados, ou seja, uma função y = mx + n que os descreva da melhor forma possível. Assim, os parâmetros a considerar são m e n. Os valores ideais para essas variáveis seriam aqueles que verificassem as seguintes equações:

\left\{\begin{matrix}
2002m + n & = & 17734\\
2003m + n & = & 19669\\
2004m + n & = & 20943\\
2005m + n & = & 24083
\end{matrix}\right.

Note que a partir dessas equações poderiam ser definidas as funções h_i como:

h_1(m,n) = 2002m + n - 17734
h_2(m,n) = 2003m + n - 19669
h_3(m,n) = 2004m + n - 20943
h_4(m,n) = 2005m + n - 24083

É de se esperar que o sistema de equações obtido a pouco não admitirá uma solução exata, pois se tem mais equações do que variáveis.

Isso geralmente acontece, pois é comum haver uma quantidade p de equações bem maior que o número n de parâmetros a identificar. Em particular, quando todas as equações h_i são lineares em suas n variáveis, dificilmente existirá uma solução exata para o sistema linear resultante, pois este terá mais equações do que incógnitas (como no exemplo). Em geral, não é possível encontrar parâmetros que satisfaçam exatamente todas as equações. Por isso, costuma-se tentar identificar os parâmetros que "melhor se aproximam" de uma solução exata, em algum sentido.

Uma forma de obter uma solução aproximada (uma "quase-solução") resulta da seguinte observação: o valor de cada função h_i em uma solução exata x deveria ser zero. Se tal exigência é restritiva demais, e com ela não é possível encontrar qualquer solução, uma possibilidade seria exigir um pouco menos. Por exemplo, poderia ser exigido apenas que o valor de h_i, para i = 1, \ldots, p seja, em geral, pequeno. Uma das formas de capturar essa idéia em termos mais precisos é dizer que se pretende minimizar a soma dos quadrados dos valores de cada h_i. Em símbolos, o problema passaria a ser:

\min_{x \in \mathbb{R}^n} \sum_{i=1}^{p}\left[h_i(x)\right]^2

[editar] O caso linear

Neste caso, para cada índice i = 1, \ldots, p, a função h_i é afim linear, ou seja:

h_i:\mathbb{R}^n \mapsto \mathbb{R}
h_i(x) = a_i^\top x + b_i

onde a_i \in \mathbb{R}^n e b_i \in \mathbb{R} para cada i = 1, \ldots, p. Deste modo, pode-se definir uma função H como

H:\mathbb{R}^n \mapsto \mathbb{R}^p de modo que
H(x)
= \begin{bmatrix} h_1(x)\\ \vdots \\ h_p(x)\end{bmatrix}
= \begin{bmatrix} a_1^\top x + b_1\\ \vdots \\ a_p^\top x + b_p\end{bmatrix}
= \begin{bmatrix} a_1^\top\\ \vdots \\ a_p^\top\end{bmatrix}x + \begin{bmatrix} b_1\\ \vdots \\ b_p\end{bmatrix}

Motivando a introdução da seguinte notação:

A = \begin{bmatrix} a_1^\top\\ \vdots \\ a_p^\top\end{bmatrix},\quad b = \begin{bmatrix} b_1\\ \vdots \\ b_p\end{bmatrix}

Assim,

H(x) = Ax + b

Logo, buscar uma solução exata é o mesmo que procurar uma solução para o seguinte sistema linear:

Ax = -b

E uma solução aproximada poderia ser buscada a partir do seguinte problema de minimização:

\min_{x \in \mathbb{R}^n} \frac{1}{2}\|Ax + b\|^2
Exercício

Verifique que se houver uma solução exata, ela também será um minimizador de \sum_{i=1}^{p}\left[h_i(x)\right]^2. A recíproca é verdadeira?

Resolução
Primeiramente, observe que para qualquer x tem-se
\sum_{i=1}^{p}\left[h_i(\bar{x})\right]^2 \ge 0

Isso significa que 0 é uma cota inferior para os valores assumidos por esta soma de quadrados. Além disso, se \bar{x} é uma solução exata, então para cada índice i tem-se:

h_i(\bar{x}) = 0

e consequentemente

\sum_{i=1}^{p}\left[h_i(\bar{x})\right]^2 = \sum_{i=1}^{p} 0^2 = 0

Então, em \bar{x} a soma dos quadrados assume seu valor mínimo.

Não vale a recíproca, pois se para um certo \tilde{x} a soma dos quadrados assumir um valor positivo \epsilon, então

\sum_{i=1}^{p}\left[h_i(\bar{x})\right]^2 = \epsilon > 0 implica que existe algum i tal que h_i(\tilde{x}) > 0

Isto significa que \tilde{x} não satisfaz a equação de índice i, e portanto não pode ser uma solução exata.

Analisando a função objetivo do problema de minimização anterior, tem-se:

f(x) = \frac{1}{2}\|Ax + b\|^2
= \frac{1}{2} \langle Ax + b, Ax + b\rangle
= \frac{1}{2}x^\top A^\top A x + b^\top A x + \frac{1}{2}b^\top b

Logo, como B = A^\top A é simétrica e semi-definida positiva, tem-se f convexa. Isso implica que a condição necessária de primeira ordem é também suficiente. Assim, qualquer ponto x \in \mathbb{R}^n tal que \nabla f(x) = 0 é solução do problema aproximado.

Calculando o gradiente da função objetivo tem-se:

0 = \nabla f(x) = A^\top A x + A^\top b

Deste modo, a solução do problema é obtida resolvendo o sistema

A^\top A x = -A^\top b

Observe também que isso implica A^\top (A x + b) = 0, ou seja, A x + b \in \ker(A^\top).


[editar] Exemplo de aplicação

Considere dado um conjunto de pontos do plano \mathbb{R}^2, por exemplo \{(x_i,y_i)\}_{i=1}^{p}, representando dados obtidos experimentalmente.

Perguntas:

1. Qual é a função afim linear que melhor se aproxima dos dados experimentais?
Resolução
As funções afim lineares no plano são aquelas que possuem a forma mx + b. Então os parâmetros a identificar são m e b.

Como são as funções h_i? Lembrando que tais funções teriam imagem igual a zero em uma solução exata, pode-se tomar h_i(m,b) = mx_i + b-y_i

Assim,

H(x)
= \begin{bmatrix} h_1(x)\\ \vdots \\ h_p(x)\end{bmatrix}
= \begin{bmatrix} mx_1 + b-y_1\\ \vdots \\mx_p + b-y_p\end{bmatrix}
= \begin{bmatrix} mx_1 + b\\ \vdots \\ mx_p + b\end{bmatrix} + \begin{bmatrix} -y_1\\ \vdots \\ -y_p\end{bmatrix}
= \begin{bmatrix} x_1 & 1 \\ \vdots & \vdots \\ x_p & 1 \end{bmatrix} \begin{bmatrix} m \\ b\end{bmatrix} + \begin{bmatrix} -y_1\\ \vdots \\ -y_p\end{bmatrix}

A expressão final pode ser simplificada adotando a seguinte notação


A = \begin{bmatrix} x_1 & 1 \\ \vdots & \vdots\\ x_p & 1 \end{bmatrix},\quad
x = \begin{bmatrix} m\\ c\end{bmatrix},\quad
c = \begin{bmatrix} -y_1\\ \vdots \\ -y_p\end{bmatrix}

De fato, nesses termos tem-se

H(x) = Ax + c
2. Qual é a função quadrática que melhor se aproxima dos dados experimentais?
Resolução
As funções quadráticas possuem a forma ax^2 + bx + c. Então os parâmetros a identificar são a, b e c.

Pode-se tomar h_i(a,b,c) = ax_i^2 + bx_i + c - y_i.

Procedendo como antes, tem-se

H(x) = Ax + c

onde


A = \begin{bmatrix} x_1^2 & x_1 & 1 \\ \vdots & \vdots & \vdots\\ x_p^2 & x_p & 1 \end{bmatrix},\quad
x = \begin{bmatrix} a\\ b \\ c\end{bmatrix},\quad
c = \begin{bmatrix} -y_1\\ \vdots \\ -y_p\end{bmatrix}
Exercício

Explorar o exemplo anterior com dados concretos, implementando um dos três métodos de gradientes conjugados.

Resolução
A resolução deste exercício é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a na versão online deste material.
Observações

Conforme se aumenta o grau do polinômio que faz a aproximação dos dados, as colunas de A têm elementos elevados a potências cada vez maiores, fazendo com que os autovalores de A^\top A sejam cada vez mais dispersos. Com isso, A^\top A torna-se mal condicionada.

Exercício

É possível, usando o problema de mínimos quadrados, verificar o postulado de Euclides que diz que "por dois pontos distintos passa uma única reta"?

Resolução
A resolução deste exercício é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a na versão online deste material.
Exercício

Usando o problema de mínimos quadrados, é possível inferir que "por três pontos distintos passa uma única curva quadrática"?

Resolução
A resolução deste exercício é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a na versão online deste material.

[editar] O caso não linear

Para esse tipo de problemas, há dois métodos:

  1. Gauss-Newton
  2. Levemberg-Marquardt

Ambos são métodos do tipo Newton. Então, para entender cada um deles é preciso entender o Método de Newton.

[editar] O método de Newton

Para entender a essência do método de Newton, primeiro considere que o problema a ser resolvido é

(P)\left\{\begin{matrix}
\min f(x)\\
x \in \mathbb{R}^n
\end{matrix}\right.

sendo h_j:\mathbb{R}^n\rightarrow \mathbb{R}, e portanto f=\frac{1}{2}\sum_{j=1}^p[h_j]^2 de classe \mathcal{C}^2 \left(\mathbb{R}^n\right). A idéia de Newton é usar o desenvolvimento até segunda ordem da série de Taylor da função f em cada ponto iterado. Isto é, se o iterado é \bar{x}, então:

Q(x) = f(\bar{x}) + \nabla f(\bar{x})^\top (x-\bar{x}) + \frac{1}{2}(x-\bar{x})\nabla^2 f(\bar{x}) (x-\bar{x})


Então a condição de Newton é que em cada iteração a Hessiana \nabla^2 f deve ser definida positiva.

Calculando o gradiente da função Q, segue:

\nabla Q(x) = \nabla f(\bar{x}) + \nabla^2 f(\bar{x}) (x-\bar{x})

Se x^* é o (único) minimizador de Q, então

0 = \nabla Q(x^*) = \nabla f(\bar{x}) + \nabla^2 f(\bar{x}) (x^*-\bar{x})

donde

\nabla^2 f(\bar{x}) (x^*-\bar{x}) = - \nabla f(\bar{x})

Sendo \nabla^2 f definida positiva, tal matriz é também inversível. Portanto:

x^* = \bar{x} - \left[\nabla^2 f(\bar{x})\right]^{-1} \nabla f(\bar{x})

Assim, pode-se usar a seguinte iteração:

x_{k+1} = x_k - \left[\nabla^2 f(x_k)\right]^{-1} \nabla f(x_k)

[editar] Algoritmo de Newton (puro)

Início: Tome x_0 \in \mathbb{R}^n
  Se \nabla f(x_0) = 0, pare: x_0 é ponto crítico.
  Senão, Calcule d_0, solução de
    \nabla^2 f(x_0) d_0 = -\nabla f(x_0)
    Faça x_1 = x_0 + d_0 e k=1
Iteração: Se \nabla f(x_k) = 0, pare: x_k é ponto crítico.
  Senão, calcular d_k, solução de
    \nabla^2 f(x_k) d_k = -\nabla f(x_k)
    Faça x_{k+1} = x_k + d_k e k=k+1

Voltando ao problema original, de mínimos quadrados, se tinha:

f(x) = \frac{1}{2}\sum_{i=1}^{p}\left[h_i(x)\right]^2

Calculando o gradiente desta função, resulta:

\nabla f(x) = \sum_{i=1}^{p}h_i(x)\nabla h_i(x)

Considera-se H:\mathbb{R}^n \mapsto \mathbb{R}^p definida por

H(x) = \begin{bmatrix} h_1(x)\\ \vdots \\ h_p(x)\end{bmatrix}

Deste modo, a Jacobiana de H verifica:

\nabla f(x) = \left[J_H(x)\right]H(x)

pois o produto de uma matriz por um vetor tem como resultado um vetor que é a combinação linear das colunas da matriz, com coeficientes que são as coordenadas do vetor.

Além disso, tem-se

\nabla^2 f(x) = \sum_{i=1}^{p} h_i(x)\nabla^2 h_i(x) + \sum_{i=1}^{p} \nabla h_i(x) \nabla h_i(x)^\top

Seja

B(x) = \sum_{i=1}^{p} \nabla h_i(x) \nabla h_i(x)^{\top} = J_H(x) J_H(x)^{\top}

Sabe-se que uma matriz D é definida positiva se x^tDx > 0 para qualquer x \not = 0. Fazendo D = \nabla h_i(x) \nabla h_i(x)^\top, tem-se:

x^\top \nabla h_i(x) \nabla h_i(x)^\top x  = \left[\nabla h_i(x)^\top x\right]^2 \ge 0

Para que D seja definida positiva, é necessário que S_P \left(\nabla h_i(x)\right) = \mathbb{R}^n ( deve gerar todo o espaço), neste caso, se diz que J_H(x) é de posto máximo.

[editar] Algoritmo de Gauss-Newton

Início: Tome x_0 \in \mathbb{R}^n
  Se \nabla f(x_0) = 0, pare: x_0 é ponto crítico.
  Senão, Calcule d_0, solução de
    B(x_0) d_0 = -\nabla f(x_0), onde B(x) = \sum_{i=1}^{p} \nabla h_i(x) \nabla h_i(x)^{\top}
    Faça x_1 = x_0 + d_0 e k=1
Iteração: Se \nabla f(x_k) = 0, pare: x_k é ponto crítico.
  Senão, calcular d_k, solução de
    B(x_k) d_k = -\nabla f(x_k)
    Faça x_{k+1} = x_k + d_k e k=k+1

[editar] Algoritmo de Levemberg-Marquardt

A idéia de Levemberg-Marquardt foi perturbar a matriz B(x), considerando B(x) + \rho I, para algum \rho>0 pequeno.

Início: Tome x_0 \in \mathbb{R}^n
  Se \nabla f(x_0) = 0, pare: x_0 é ponto crítico.
  Senão, Calcule d_0, solução de
    \left(B(x_0) + \rho I\right) d_0 = -\nabla f(x_0), onde B(x) = \sum_{i=1}^{p} \nabla h_i(x) \nabla h_i(x)^{\top}
    Faça x_1 = x_0 + d_0 e k=1
Iteração: Se \nabla f(x_k) = 0, pare: x_k é ponto crítico.
  Senão, calcular d_k, solução de
    \left(B(x_0) + \rho I\right) d_k = -\nabla f(x_k)
    Faça x_{k+1} = x_k + d_k e k=k+1
Exercício

Enumere as diferenças que existem nos seguintes algoritmos:

  • Newton
  • Gauss-Newton
  • Levenberg-Marquardt
Resolução
Newton

Desenvolvendo f em série de Taylor e exigindo que a matriz hessiana de f seja definida positiva o método de Newton usa d_k=-[\nabla^2f(x_k)]^{-1}\nabla f(x_k).

Gauss-Newton

Em alguns casos a matriz hessiana de f pode ser bem complicada. Nesses casos então a idéia do método de Gauss-Newton é usar uma boa aproximação para essa matriz. Podemos escrever

f(x)=\frac{1}{2}\sum_{i=1}^p [h_i(x)]^2=\frac{1}{2}\|h(x)\|^2=\frac{1}{2}h(x)^{\top}h(x).

Usando Taylor em h temos:

h(\bar{x})=h(x)+J(x)(\bar{x}-x) onde J(x) é a jacobiana de h.

Substituindo h(\bar{x}) na f temos


\begin{align}
f(\bar{x})  & =\frac{1}{2}h(\bar{x})^{\top}h(\bar{x})\\

            & =\frac{1}{2}(h(x)+J(x)(\bar{x}-x))^{\top}(h(x)+J(x)(\bar{x}-x))\\

            & =\frac{1}{2}[(h(x)^{\top}h(x)+2h(x)^{\top}(J(x)(\bar{x}-x))+(J(x)(\bar{x}-x))^{\top}(J(x)(\bar{x}-x))]\\

            & =f(x)+h(x)^{\top}(J(x)(\bar{x}-x))+ (J(x)(\bar{x}-x))^{\top}(J(x)(\bar{x}-x))
\end{align}

Calculando o gradiente e depois a segunda derivada temos:

\nabla f(\bar{x})=J(x)^{\top}h(x)+J(x)^{\top}(J(x)(\bar{x}-x)) e

\nabla^2f(\bar{x})=J(x)^{\top}J(x)

Denotando B(x)=J(x)^{\top}J(x), então o método de Gauss-Newton escolhe d_k tal que

B(x_k)d_k=-\nabla f(x_k).


Levenberg-Marquardt

A idéia do método de Levenberg-Marquardt é pertubar a matriz B(x) um pouco, considerando B(x)+\rho I com \rho>0 pequeno.

Portanto o método de Levenberg-Marquardt escolhe d_k tal que

(B(x_k)+\rho I)d_k=-\nabla f(x_k).

KKT


Neste capítulo o objetivo é desenvolver algumas ideias e provar o teorema de Karush–Kuhn–Tucker (também chamado simplesmente de teorema KKT) que será utilizado no capítulo seguinte para explorar os métodos duais. O teorema KKT é bem útil para resolver problemas do tipo

(P)\left\{\begin{matrix}
\min f(x)\\
g_i(x) \le 0; i = 1, \ldots, p\\
h_j(x) = 0; j = 1, \ldots, q
\end{matrix}\right.

[editar] Cones

Definição

Um conjunto C \in \mathbb{R}^n é um cone quando

d \in C \Rightarrow td \in C, \forall t \in \mathbb R_{+}
Exemplo de um cone no \mathbb{R}^2.

Em outras palavras, a propriedade que caracteriza um cone é que este tipo de conjunto contém todos os múltiplos não nulos de qualquer de seus elementos.

Definição

Dado um subconjunto C \subset \mathbb{R}^n, o cone polar de C é o conjunto definido por

C^{*}=\{ p\in \mathbb{R}^n : p^{\top}x\leq 0,\ \forall x \in C \}.

Observações:

  • C^{*} é um cone: Se d \in C^{*} tem-se que d^{\top}x\leq 0,\ \forall x\in C. Logo, para qualquer t\in \mathbb R_+, vale (td)^{\top}x=t d^{\top}x\leq 0,\ \forall x\in C . Disto segue que td\in C^{*}, mostrando que C^{*} é um cone.
  • Sempre se tem que C\subseteq (C^{*})^{*} (Verifique).

Na segunda propriedade a igualdade pode não ocorrer (exemplo?). Para o objetivo deste texto, o ideal seria que a igualdade valesse. Mas será que isso ocorre para algum conjunto? A resposta é sim e, conforme o próximo lema, basta que C seja um cone convexo fechado.


Crystal Clear app kaddressbook.png Este módulo tem a seguinte tarefa pendente: Incluir a definição de projeção antes deste ponto, pois ela será usada durante a demonstração
Lema  (Farkas)

Se C\subset \mathbb{R}^n um cone convexo fechado, então C=(C^{*})^{*}.

Demonstração
Seja y\in (C^{*})^{*} e w=\text{proj}_{C}(y). Sabendo que a projeção de um ponto sobre um conjunto convexo é única, será mostrado que w=y e então ficará provada a inclusão C \supset (C^{*})^{*}. Disto seguirá a igualdade entre os dois conjuntos, já que C é sempre um subconjunto de (C^{*})^{*}.

Pelo teorema da projeção (ver Izmailov & Solodov (2007)), tem-se que (y-w)^{\top}(x-w)\leq 0,\ \forall x\in C. Usando o fato de que C é cone, segue que 0\in C e 2w\in C e substituindo isto na equação acima obtem-se

(y-w)^{\top}(-w)\leq 0 e (y-w)^{\top}w\leq 0.

Dessas desigualdades, conclui-se que (y-w)^{\top}w=0.

De (y-w)^{\top}(x-w)=(y-w)^{\top}x -(y-w)^{\top}w\leq 0, tem-se que (y-w)^{\top}x\leq 0,\ \forall x\in C.

Usando a definição de cone polar, isso implica que y-w\in C^{*}.

Uma vez que y\in (C^{*})^{*}, significa que (y-w)^{\top}y\leq 0.

Desses fatos acima se conclui que

\|y-w\|^2=(y-w)^{\top}(y-w)=(y-w)^{\top}y-(y-w)^{\top}w\leq 0

Isso mostra que y=w.


Definição

Dado x\in C, se diz que d\in \mathbb{R}^n é uma direção viável em x, com respeito a C, quando existe \epsilon > 0 tal que

x+td \in C,\ \forall t\in [0, \epsilon].
O conjunto de todas as direções viáveis em x, com respeito ao conjunto C, será denotado por V_C(x).

Esse conjunto V_C(x) é o cone das direções viáveis em x, com respeito a C.

Definição

Uma direção d\in \mathbb{R}^n é uma direção de descida da função f em x, se existe \epsilon >0 tal que

f(x+td) < f(x),\ \forall t\in (0, \epsilon].
O conjunto das direções de descida será denotado por D(x).

[editar] Caracterização das direções de descida

Lema

Seja f: \mathbb{R}^n \rightarrow \mathbb{R} uma função diferenciável em x\in \mathbb{R}^n. Então

  1. \nabla f(x)^{\top}d \leq 0,\ \forall d\in D(x).
  2. d\in \mathbb{R}^n satisfaz \nabla f(x)^{\top}d < 0 \ \Rightarrow d\in D(x).


Demonstração
1) Seja d\in D(x). Então, \forall t \in (0, \epsilon] tem-se f(x+td)<f(x).

Usando a série de Taylor, tem-se

f(x)+t\nabla f(x)^{\top}d + o(t)< f(x)

Sendo t\neq 0,

\nabla f(x)^{\top}d + \frac{o(t)}{t}< 0.

Passando o limite com t\rightarrow 0 tem-se que \nabla f(x)^{\top}d\leq 0.

2) Novamente aplicando Taylor em f(x+td)-f(x) tem-se

f(x+td)-f(x)=t\nabla f(x)^{\top}d + o(t).

Como t\neq 0, tem-se

f(x+td)-f(x)=t(\nabla f(x)^{\top}d +\frac{o(t)}{t} ).

Pela hipótese \nabla f(x)^{\top}d<0 , com isto

\lim_{t\rightarrow 0}{(\nabla f(x)^{\top}d +\frac{o(t)}{t})}=\nabla f(x)^{\top}d<0.

Pelo teorema da conservação do sinal, existe \epsilon >0 tal que \nabla f(x)^{\top}d +\frac{o(t)}{t}<0,\ \forall t\in (0, \epsilon].

Portanto,

t(\nabla f(x)^{\top}d +\frac{o(t)}{t})<0
f(x+td)<f(x)\ \forall t\in (0, \epsilon].

Conclui-se então que d\in D(x).

[editar] O cone viável linearizado

Definição

Dado x \in C = \{ x \in \mathbb{R}^n; g_i(x) \le 0 \text{ e } h_j(x) = 0\}, a desigualdade g_i(x) \le 0 é uma restrição ativa em x se g_i(x)=0.

Observações
  • O conjunto formado pelos índices das restrições de desigualdade ativas é denotado por I(x). Assim,
I(x)=\{i:g_{i}(x)=0\}
Definição

Dado um ponto x\in C e o conjunto I(x), se define o cone viável linearizado de C a partir de x como

L(x, C)=\{d\in \mathbb{R}^n: \nabla g_{j}(x)^{\top}d\leq 0,\ \forall j\in I(x) \ \text{e} \ \nabla h_{i}(x)^{\top}d=0,\ \forall i=1, \dots, q\}.

L(x, C) é um cone não-vazio convexo e fechado pois, 0\in L(x, C). E se y, w \in L(x, C), tem-se

\nabla h_{i}(x)^{\top}(\alpha y+(1-\alpha)w)= \alpha \nabla h_{i}(x)^{\top}y +(1-\alpha)\nabla h_{i}(x)^{\top}w=\alpha 0 +(1-\alpha)0=0 e
\nabla g_{j}(x)^{\top}(\alpha y+(1-\alpha)w)= \alpha \nabla g_{j}(x)^{\top}y +(1-\alpha)\nabla g_{j}(x)^{\top}w\leq \alpha 0 +(1-\alpha)0\leq 0.

Portanto \alpha y+(1-\alpha)w \in L(x, C) mostrando que L(x, C) é convexo.

Para mostrar que L(x, C) é fechado, pode-se pegar uma sequência convergente (d^{k})\in L(x, C) e mostrar que o ponto de acumulação dela esta em L(x, C).

Tem-se que \nabla h_{i}(x)^{\top}d^{k}=0 \ \text{e}\ \nabla g_{j}(x)^{\top}d^{k}\leq 0,\ \forall k\in \mathbb{N}.

Passando o limite com k \rightarrow \infty, obtem-se

0=\lim_{k \rightarrow \infty}{\nabla h_{i}(x)^{\top}d^{k}}=\nabla h_{i}(x)^{\top}\lim_{k \rightarrow \infty}{d^{k}}=\nabla h_{i}(x)^{\top}d e
0\ge \lim_{k \rightarrow \infty}{\nabla g_{j}(x)^{\top}d^{k}}= \nabla g_{j}(x)^{\top}\lim_{k \rightarrow \infty}{d^{k}}=\nabla g_{j}(x)^{\top}d.

Isso mostra que L(x, C) é fechado.


Lema  (Caratheodory)

Sejam y_1, \dots, y_m, w_1, \dots, w_p \in \mathbb{R}^n. Seja x\in \mathbb{R}^n com x\neq 0 e \alpha_1, \dots, \alpha_m, \beta_1, \dots, \beta_p escalares tais que\beta_j\ge 0\ \forall j=1, \dots, p e

x=\sum_{i=1}^m \alpha_{i}y_{i}+\sum_{j=1}^p \beta_{j}w_j.
Então existem subconjuntos I\subset \{1, \dots, m\}\text{, } J\subset\{1, \dots, p\} e escalares \alpha_{i}^{*} com i\in I e \beta_{j}^{*}\ \forall j\in J tais que
x=\sum_{i\in I} \alpha_{i}^{*}y_{i}+\sum_{j\in J} \beta_{j}^{*}w_j e os vetores \{y_i\}_{i\in I}\cup \{w_j\}_{j\in J} são linearmente independentes.
Demonstração
Sem perda de generalidade, suponha que \alpha_{i}\neq 0\ \forall i=1, \dots, m e \beta_j>0,\ \forall j=1, \dots, p. Considere que \{y_1, \dots, y_m, w_1, \dots, w_p\} sejam linearmente dependentes.

Portanto existem escalares \lambda_i com i=1, \dots, m e \delta_j com j=1, \dots, p não todos nulos tais que

0=\sum_{i=1}^m \lambda_{i}y_{i}+\sum_{j=1}^p \delta_{j}w_j

Multiplicando a igualdade acima por t e subtraindo de

x=\sum_{i=1}^m \alpha_{i}y_{i}+\sum_{j=1}^p \beta_{j}w_j tem-se
x=\sum_{i=1}^m (\alpha_{i}-t\lambda_i)y_{i}+\sum_{j=1}^p (\beta_{j}-t\delta_j)w_j

Para t=0 certamente nenhum dos coeficientes acima se anula.

Seja \bar{t} o t de menor módulo que anula pelo menos um dos coeficientes \alpha_{i}-t\lambda_i ou \beta_{j}-t\delta_j. Então

x=\sum_{i=1}^m (\alpha_{i}-\bar{t}\lambda_i)y_{i}+\sum_{j=1}^p (\beta_{j}-\bar{t}\delta_j)w_j

Assim, se escreve x como combinação linear de no máximo m+p-1 vetores já que \beta_{j}-\bar{t}\delta_j\ge 0.

Repetindo esse processo obtem-se uma combinação linearmente independente.


Definição

Dado um ponto x\in C, se define o cone G(x) por

G(x)=\{\sum_{i=1}^q \alpha_{i}\nabla h_{i}(x) +\sum_{j\in I(x)} \beta_{j}\nabla g_{j}(x):\beta_{j}\ge 0,\ \forall j\in I(x)\}.

A seguir, serão mostradas algumas propriedades deste cone.

Lema

Para qualquer x\in C, G(x) é um cone convexo e fechado.

Demonstração
Primeiro será mostrado que G(x) é de fato um cone. Seja d\in G(x) e t\ge 0. Então tem-se
td=\sum_{i=1}^q t\alpha_{i}\nabla h_{i}(x) +\sum_{j\in I(x)} t\beta_{j}\nabla g_{j}(x).

Como t\beta_{j}\ge 0 tem-se que td\in G(x).

Agora, será provado que G(x) é convexo. Para isso seja y, w\in G(x), isto é,

y=\sum_{i=1}^q \alpha_{i}\nabla h_{i}(x) +\sum_{j\in I(x)} \beta_{j}\nabla g_{j}(x) e
w=\sum_{i=1}^q \lambda_{i}\nabla h_{i}(x) +\sum_{j\in I(x)} \delta_{j}\nabla g_{j}(x) e t\in [0, 1].

Logo tem-se,

ty+(1-t)w=\sum_{i=1}^q (t\alpha_{i}+(1-t)\lambda_{i})\nabla h_{i}(x) +\sum_{j\in I(x)} (t\beta_{j}+(1-t)\delta_{j})\nabla g_{j}(x).

Como t\beta_{j}+(1-t)\delta_{j}\ge 0 visto que \beta_{j}\ge 0 e \delta_{j}\ge 0. Com isso concluímos que ty+(1-t)w\in G(x) mostrando que G(x) é convexo.

Para mostrar que G(x) é fechado, toma-se uma sequência convergente em G(x) e se mostra que o ponto de acumulação dela pertence a G(x).

Para isso seja (z^k)\subset G(x) com z^k\rightarrow z \in \mathbb{R}^n. Será mostrado que z\in G(x).

Escrevendo G(x) em forma matricial tem-se G(x)=\{A\Delta+B\Omega:\Omega\ge 0 \}.

Pelo Lema de Caratheodory podemos assumir que C=(A\ B) tem colunas linearmente independentes, e portanto C^{\top}C é não singular.

Uma vez que (z^k)\subset G(x), existem \Gamma^{k}=(\Delta^{k}\ \Omega^{k})^t com \Omega^{k}\ge 0 tais que z^{k}=C\Gamma^{k}.

Uma vez que C^{\top}C é não singular, \Gamma^{k}=(C^{\top}C)^{-1}C^{\top}z^k.

Passando o limite obtem-se,

(\Delta \ \Omega)^t=\Gamma=\lim_{k\rightarrow \infty}{\Gamma^{k}}=(C^{\top}C)^{-1}C^{\top}z com \Omega\ge 0.

Isso mostra que C\Omega \in G(x).

Agora passando o limite em z^k=C\Omega^k obtém-se z=C\Omega, mostrando que z\in G(x).


Lema

Para qualquer x\in C, G(x)=L(x, C)^{*}.

Demonstração
Como L(x, C) e G(x) são convexos e fechados, tem-se que L(x, C)=(L(x, C)^{*})^{*} e G(x)=(G(x)^{*})^{*}. Será mostrado então que L(x, C)=G(x)^{*}.

Seja d\in L(x, C). Assim, dado y\in G(x) tem-se

d^{\top}y=d^{\top}(\sum_{i=1}^q \alpha_{i}\nabla h_{i}(x) +\sum_{j\in I(x)} \beta_{j}\nabla g_{j}(x))
d^{\top}y=\sum_{i=1}^q \alpha_{i}d^{\top}\nabla h_{i}(x) +\sum_{j\in I(x)} \beta_{j}d^{\top}\nabla g_{j}(x)

Mas \beta\ge 0 e d^{\top}\nabla h_i(x)=0 e d^{\top}\nabla g_j(x)\leq 0.

Conclui-se então que d^{\top}y\leq 0. Como y é arbitrário, d\in G(x)^{*}.

Agora a volta, seja d\in G(x)^{*}, isto é, d^{\top}y\leq 0\ \forall y\in G(x).

Em particular, uma vez que \nabla h_i(x) e -\nabla h_i(x)\in G(x)\ \forall i=1, \dots, q, tem-se que d^{\top}\nabla h_i(x)=0.

Além disso, uma vez que \nabla g_j(x)\in G(x)\ \forall j\in I(x), tem-se que d^{\top}\nabla g_j(x)\leq 0.

Logo d\in L(x, C).

[editar] O cone tangente

Definição

Um vetor d\in \mathbb{R}^n é chamado direção tangente em C a partir de x\in C quando ou d=0 ou \exist (x^k)\subset C tal que

x^k\rightarrow x e \frac{x^k-x}{\|x^k-x\|}\rightarrow \frac{d}{\|d\|}.
Observações
  • O conjunto de todas as direções tangentes no ponto x \in C, é denominado cone tangente, e denotado por T(x, C).
  • Se a \in C, então T(a, C) também pode ser descrito como
T(a, C) = \left\{d \in; \exists \{d_k\} \text{ com } d_k \to d;\quad \exists \{t_k\} \text{ com } t_k \to 0 \text{ tais que } x = a + t_k d_k \in C, \, \forall k \right\}
Exercício

Verifique que T(a, C) é de fato um cone (e portanto merece ser chamado de "cone tangente").

Resolução
A resolução deste exercício é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a na versão online deste material.

[editar] Exemplo de cone tangente

Determinar o cone tangente ao ponto a = (0, 0) do quadrado unitário com vértices (0, 0), (0, 1), (-1, 1) e (-1, 0).

Resolução
Cone tangente a um quadrado unitário com vértice na origem.

Dado qualquer ponto d = (d_0, d_1) do 2º quadrante (formado pelos pontos (x, y) tais que x<0 e y>0), pode-se definir:

t_k = \left(\frac{1}{2}\right)^k
d_k = d

Com essas escolhas, tem-se:

t_k \to 0
d_k \to d

Logo, a + t_k d_k = a + \left(\frac{1}{2}\right)^k d = (0, 0) + \left(\frac{d_0}{2^k}, \frac{d_1}{2^k}\right) = \left(\frac{d_0}{2^k}, \frac{d_1}{2^k}\right) \in C.


[editar] Propriedades do cone tangente

Wikipedia
A Wikipédia tem mais sobre este assunto:
Cone tangente

O cone tangente definido anteriormente tem as seguintes propriedades:

  1. T(a, C) é fechado e 0 \in T(a, C)
  2. Se C \subset D então T(a, C) \subset T(a, D)
  3. Se V é uma vizinhança de a, então T(a, C) = T(a, V \cap C)
Observação

A terceira propriedade indica que o cone tangente só depende do que ocorre bem perto de a, no conjunto C.


Lema

Para qualquer x\in C, T(x, C) é fechado.

Demonstração
Seja (d^k)\subset T(x, C) com d^k\rightarrow d\in \mathbb{R}^n. Será mostrado que d\in T(x, C).

Caso d=0, d\in T(x, C). Então, suponha-se que d\neq 0.

Neste caso, sem perda de generalidade pode-se considerar que d^k\neq 0,\ \forall k\in \mathbb{N}, pois d^k\rightarrow d.

Fixando k\in \mathbb{N} tem-se que d^k\in T(x, C). Portanto, existe (x^{k, j})_{j\in \mathbb{N}}\subset C tal que x^{k, j}\rightarrow x e \frac{x^{k, j}-x}{\|x^{k, j}-x\|}\rightarrow \frac{d^k}{\|d^k\|} quando j\rightarrow \infty.

Assim para \epsilon=\frac{1}{k} existe j_k\in \mathbb{N} tal que para j\ge j_k, tal que \|x^{k, j}-x\|<\frac{1}{k} e \bigg|\frac{x^{k, j}-x}{\|x^{k, j}-x\|}- \frac{d^k}{\|d^k\|}\bigg|<\frac{1}{k}.

Em particular, tomando j=j_k tem-se

\|x^{k, j_k}-x\|<\frac{1}{k} e \bigg|\frac{x^{k, j_k}-x}{\|x^{k, j_k}-x\|}- \frac{d^k}{\|d^k\|}\bigg|<\frac{1}{k}.

Tomando o limite quando k\rightarrow \infty, obtem-se que x^k\rightarrow x e

\bigg|\frac{x^{k, j_k}-x}{\|x^{k, j_k}-x\|}- \frac{d}{\|d\|}\bigg|\leq \bigg|\frac{x^{k, j_k}-x}{\|x^{k, j_k}-x\|}- \frac{d^k}{\|d^k\|}\bigg|+\bigg|\frac{d^k}{\|d^k\|}- \frac{d}{\|d\|}\bigg|\rightarrow 0.

Logo \frac{x^k-x}{\|x^k-x\|}\rightarrow \frac{d}{\|d\|}.

Isso mostra que d\in T(x, C).



Exercício

Verificar que:

  1. T(a, C) \subset L(a, C).
  2. Se C = \{(x, y) \in \mathbb{R}^2; x^2 + y \le 0;\quad x^2 - y \le 0\} e a = (0, 0), então T(a, C) \not = L(a, C).
Demonstração
1) Seja d\in T(a, C), d\neq 0. Logo \exist (x^k)\subset C tal que x^k\neq a, x^k\rightarrow a e \frac{x^k-a}{\|x^k-a\|}\rightarrow \frac{d}{\|d\|}.

Usando Taylor em torno de a tem-se

0=h_j(x^k)=h_j(a)+\nabla h_j(x)^{\top}(x^k-a)+o(\|x^k-a\|).

Já que x^k\neq a, então \|x^k-a\|\neq 0 logo pode-se dividir e obtem-se

\nabla h_j(a)^{\top}\frac{(x^k-a)}{\|x^k-a\|}+\frac{o(\|x^k-a\|)}{\|x^k-a\|}=0.

Passando o limite quando k\rightarrow \infty, tem-se \nabla h_j(a)^{\top}\frac{d}{\|d\|}=0.

Novamente usando Taylor em torno de a para i\in I(x) tem-se

g_i(a)+\nabla g_i(a)^{\top}(x^k-a)+o(\|x^k-a\|)\leq 0
\nabla g_i(a)^{\top}\frac{(x^k-a)}{\|x^k-a\|}+\frac{o(\|x^k-a\|)}{\|x^k-a\|}\leq 0

Passando o limite quando k\rightarrow \infty tem-se \nabla g_i(a)^{\top}\frac{d}{\|d\|}=\leq 0.

Donde se conclui que d\in L(a, C).

2)

Crystal Clear app kaddressbook.png Este módulo tem a seguinte tarefa pendente: Colocar figura
Lema

Se a\in C é um mínimo local do problema (P), então \nabla f(a)^{\top}d\ge 0,\ \forall d\in T(a, C).

Demonstração
Por Taylor tem-se
0\ge f(x^k)-f(a)=\nabla f(a)^{\top}(x^k-a)+o(\|x^k-a\|)
0\ge \nabla f(a)^{\top}\frac{(x^k-a)}{\|x^k-a\|}+\frac{o(\|x^k-a\|)}{\|x^k-a\|}

Passando o limite quando k\rightarrow \infty obtem-se

0\ge \nabla f(a)^{\top}\frac{d}{\|d\|}

Donde \nabla f(a)^{\top}d\leq 0\ \forall d\in T(a, C).


[editar] Teorema KKT

Teorema (Condições de KKT)

Seja C = \{ x \in \mathbb{R}^n; g_i(x) \le 0 \text{ e } h_j(x) = 0\} e considere a\in C um minimizador local do problema

(P)\left\{\begin{matrix}
\min f(x)\\
x \in C
\end{matrix}\right.
Se T(a, C)^{*}=L(a, C)^{*}, então existem u \in \mathbb{R}^p e v \in \mathbb{R}^q tais que:
  1. -\nabla f(a) = \sum_{i=1}^p u_i\nabla g_i(a) + \sum_{j=1}^q v_j\nabla h_j(a)
  2. u_i \ge 0,\ \forall i=1, \dots, p
  3. u_ig_i(a)=0,\ \forall i=1, \dots, p.
Demonstração
Considere a um minimizador local do problema (P). Então (-\nabla f(a))^{\top}d\leq 0,\ \forall d\in T(a, C). Pela definição de cone polar isso significa que -\nabla f(a)\in T(a, C)^{*}.

Pela hipotése tem-se -\nabla f(a)\in L(a, C)^{*}. Como L(a, C)=G(a)^{*} obtem-se que -\nabla f(a)\in (G(a)^{*})^{*}.

Como foi visto acima G(a) é um cone convexo e fechado. Portanto usando o Lema de Farkas obtem-se que -\nabla f(a)\in G(a).

Pela definição de G(a), existem escalares \delta_i com i\in I(a) e \lambda_j com j=1, \dots, q tais que

-\nabla f(a) = \sum_{i\in I(a)} \delta_i\nabla g_i(a) + \sum_{j=1}^q \lambda_j\nabla h_j(a) com \delta_i\ge 0\ \forall i\in I(a).

Como \text{card}I(a)\leq p, define-se v_j = \lambda_j,\ \forall j=1, \dots, q e u_i=\begin{cases} 
 \delta_i & \forall i\in I(a) \\
 0 & \forall i\notin I(a) 
\end{cases}

Como g_i(a)=0,\ \forall i\in I(a) obtem-se u_ig_i(a)=0\ \forall i=1, \dots, p.

Com isso fica provado o Teorema de KKT.


Métodos duais


[editar] Lagrangiana

O conceito de lagrangiana está sempre relacionado ao seguinte problema:

(P)\left\{\begin{matrix}
\min f(x)\\
g_i(x) \le 0; i = 1, \ldots, p\\
h_j(x) = 0; j = 1, \ldots, q
\end{matrix}\right.
Definição

A função lagrangiana associada ao problema (P) é

l:\mathbb{R}^n \times \mathbb{R}^p \times \mathbb{R}^q \mapsto \mathbb{R} definida por
l(x, u, v) = f(x) + \sum_{i=1}^{p} u_i g_i (x) + \sum_{j=1}^{q} v_j h_j (x).
Wikipedia
A Wikipédia tem mais sobre este assunto:
Multiplicadores de Lagrange

Em alguns livros é usada a seguinte notação:

G:\mathbb{R}^n \mapsto \mathbb{R}^p definida por
G(x) = \begin{bmatrix} g_1(x)\\ \vdots \\ g_p(x)\end{bmatrix}

e

H:\mathbb{R}^n \mapsto \mathbb{R}^q definida por
H(x) = \begin{bmatrix} h_1(x)\\ \vdots \\ h_q(x)\end{bmatrix}

Deste modo, a lagrangiana fica expressa por

l(x, u, v) = f(x) + u^\top G(x) + v^\top H(x)

que é uma representação bem mais compacta.

[editar] Condições de otimalidade

Para permitir a compreensão da importância da função lagrangiana em otimização, é preciso ter em mente os principais resultados que garantem a otimalidade de uma solução para o problema (P). Nas próximas seções será apresentado um breve resumo das condições de otimalidade, dividindo-as em dois casos:

  • Caso particular: quando C = \{ x \in \mathbb{R}^n; g_i(x) \le 0 \text{ e } h_j(x) = 0\} é convexo e fechado.
  • Caso geral: quando C é arbitrário.

[editar] Caso particular

Proposição

Seja f função de classe \mathcal{C}^1 no conjunto C. Se \bar{x} é um minimizador local de f no conjunto convexo e fechado C, então:

\langle \nabla f(x), x - \bar{x} \rangle \ge 0, \, \forall x \in C
Demonstração
Seja \bar{x} \in C uma solução de (P), e fixe um ponto arbitrário x \in C. Denote por \theta a restrição da função f sobre o segmento de reta que passa por x e por \bar{x}, ou seja,
\theta (t) = f(\bar{x} + t (x - \bar{x})), com t \in (0, 1).

Note que \theta(0) = f(\bar{x}).

Como \bar{x} é um minimizador de f em C, tem-se:

\theta(t) = f(\bar{x} + t (x - \bar{x})) \ge f(\bar{x}) = \theta(0), para todo t \in (0, 1)

Logo,

\frac{\theta(t) - \theta(0)}{t - 0} \ge 0, ou seja, \theta'(0) = \lim_{t \to 0} \frac{\theta(t) - \theta(0)}{t - 0} \ge 0

Mas pela regra da cadeia,

\theta'(0) = \langle \nabla f(\bar{x}), x - \bar{x}\rangle, então \langle \nabla f(\bar{x}), x - \bar{x}\rangle \ge 0.

Como x \in C era arbitrário, o resultado está demonstrado.

Observação

A condição \langle u, v \rangle \ge 0 significa que os vetores u e v formam um ângulo reto ou agudo (menor ou igual a 90 graus), conforme indicado na figura a seguir:

Em um ponto de mínimo, \nabla f sempre forma ângulo menor ou igual a 90 graus com os vetores do tipo x-\bar{x}, onde \bar{x} é um ponto de mínimo da função no conjunto viável C e x é um ponto qualquer deste conjunto.

No caso específico, com u=\nabla f e v=x-\bar{x}, \langle u, v \rangle é a derivada direcional de f na direção x-\bar{x}. Quando tal número é não negativo, intuitivamente a função cresce naquela direção.


Quando C é um conjunto convexo e fechado, a existência de uma solução para o problema (P) é garantida pela seguinte proposição:

Proposição

Se \bar{x} \in C, f é convexa e

\langle \nabla f(x), x - \bar{x} \rangle \ge 0, \, \forall x \in C
então \bar{x} é um minimizador global de f no conjunto C (isto é, \bar{x} é solução de (P)).
Demonstração
Pela convexidade de f, tem-se que:
f(x) \ge f(\bar{x}) + \langle \nabla f(x), x - \bar{x} \rangle, \, \forall x \in \mathbb{R}^n

Logo,

f(x) - f(\bar{x}) \ge \langle \nabla f(x), x - \bar{x} \rangle \ge 0, \, \forall x \in C

Portanto, f(x) \ge f(\bar{x}), \, \forall x \in C, ou seja, \bar{x} é um minimizador de f em C.

Até aqui, não é exigido que qualquer das funções g_i ou h_j sejam diferenciáveis. Isto será utilizado mais adiante, nos algoritmos.

A próxima proposição fornece uma caracterização dos minimizadores, em termos do conceito de projeção.

Lembre-se que a projeção de um ponto \bar{y} sobre o conjunto C, denotada por \bar{P} = P_C (\bar{y}), satisfaz \langle \bar{P} - \bar{y}, x - \bar{P} \rangle \ge 0, \, \forall x \in C. Na verdade, vale:

\bar{P} = P_C (\bar{y}) \Leftrightarrow \langle \bar{P} - \bar{y}, x - \bar{P} \rangle \ge 0, \, \forall x \in C
O vetor \bar{P}-\bar{y} forma ângulo menor que 90 graus com o vetor x - \bar{P}, pois \bar{P} é a projeção de \bar{y} sobre o conjunto C.
Proposição

Seja \alpha > 0 um número real fixado.

  1. Se \bar{x} é um minimizador local de f em C, então \bar{x} = P_C(\bar{x} - \alpha \nabla f (\bar{x}))
  2. Se f é convexa e \bar{x} = P_C(\bar{x} - \alpha \nabla f (\bar{x})), então \bar{x} é um minimizador global de f em C.
Demonstração
(1)

Como \bar{x} é um minimizador local de f em C, então

\langle \nabla f(x), x - \bar{x} \rangle \ge 0, \, \forall x \in C

Observe que essa desigualdade equivale à

\langle \alpha \nabla f(x), x - \bar{x} \rangle \ge 0, \, \forall x \in C

ou ainda

\langle \bar{x} - (\bar{x} - \alpha \nabla f(x)), x - \bar{x} \rangle \ge 0, \, \forall x \in C

Tomando \bar{y} = \bar{x} - \alpha \nabla f(x) e P = \bar{x}, tem-se a equivalência com:

\langle P - \bar{y}, x - P \rangle \ge 0, \, \forall x \in C

que equivale a dizer que P = P_C(\bar{y}), ou seja:

\bar{x} = P = P_C(\bar{y}) = P_C(\bar{x} - \alpha \nabla f(x))

(2)

Reciprocamente, supor que \bar{x} = P = P_C(\bar{y}) = P_C(\bar{x} - \alpha \nabla f(x)), conforme a argumentação anterior, equivale a dizer que

\langle \nabla f(x), x - \bar{x} \rangle \ge 0, \, \forall x \in C

Usando a hipótese de que f é convexa, segue da proposição anterior que \bar{x} é um minimizador global de f em C.

[editar] Caso geral

Para tratar este caso, é preciso utilizar o conceito de cone tangente. O conjunto contínua o mesmo, ou seja, C = \{ x \in \mathbb{R}^n; g_i(x) \le 0 \text{ e } h_j(x) = 0\}, embora não seja mais suposto que ele é convexo. Mesmo assim, C ainda será fechado, pois as funções g_i e h_j que o definem são contínuas.

Teorema

  1. Se \bar{x} é um minimizador local de f em C, então \langle \nabla f(\bar{x}) , d \rangle \ge 0, \, \forall d \in T(\bar{x}, C).
  2. Se C é convexo, f é convexa e \langle \nabla f(\bar{x}) , d \rangle \ge 0, \, \forall d \in T(\bar{x}, C), então \bar{x} é minimizador global de f em C.

Este teorema pode ser demonstrado de forma análoga a que foi feita anteriormente.

Demonstração
Esta demonstração é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a na versão online deste material.

Seja \tilde{C} definido como:

\tilde{C} = \tilde{C}_a = \left\{x \in \mathbb{R}^n; g_i(x) \le 0, \, \forall i \in I(a);\quad h_j(x) = 0 \right\}.

Pela segunda propriedade do cone tangente, tem-se:

T(a, V \cap \tilde{C}) = T(a, C) = T(a, V \cap \tilde{C}). Logo,
T(a, C) = T(a, \tilde{C})

Em outras palavras, se é dado um conjunto \tilde{C} e depois se restringe tal conjunto para C, através do acréscimo de restrições inativas em um ponto a, os cones tangentes aos dois conjuntos (no ponto a) coincidem.

Definição

Toda condição que implica que T(a, C) = L(a, C) é chamada de condição de qualificação das restrições.

Observação: Também se pode dizer condições de regularidade das restrições (do inglês, Regularity conditions).

[editar] Exemplos de condições de qualificação das restrições

(1) Se g_i e h_i são funções afim-lineares, então para qualquer x \in C, tem-se T(x, C) = L(x, C). Prova-se isso trivialmente

Prova
Esta prova é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a na versão online deste material.

(2) Condições de Slater: Se as funções g_i são convexas e as h_i são afim lineares e, além disso, existe \hat{x} \in C tal que g_i(\hat{x}) < 0 para todo i, h_j (\hat{x}) = 0, então para qualquer x \in C, tem-se T(x, C) = L(x, C).

Prova
Esta prova é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a na versão online deste material.

(3) Condições de Mangasarian-Fromowitz: Se \{\nabla h_j(x)\} é linearmente independente e existe \tilde{d} tal que \langle \nabla h_j(x), \tilde{d}\rangle = 0 para tpdp j e \langle \nabla g_i(x), \tilde{d}\rangle < 0.

Prova
Esta prova é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a na versão online deste material.
Teorema  (Karush–Kuhn–Tucker)

Suponha que f, g_i e h_j são funções de classe \mathcal{C}^1 em uma vizinhança do ponto \bar{x} e que T(\bar{x}, C) = L(\bar{x}, C). Se \bar{x} é um minimizador local de f em C, então existem u \in \mathbb{R}^p e v \in \mathbb{R}^q tais que:

  1. \nabla f(\bar{x}) + \sum_{i=1}^p u_i \nabla g_i(\bar{x}) + \sum_{j=1}^q v_j \nabla h_j(\bar{x}) = 0
  2. u_i \ge 0, \, \forall i = 1, \ldots, p
  3. u_i g_i(\bar{x}) = 0, \, \forall i = 1, \ldots, q
Wikipedia
A Wikipédia tem mais sobre este assunto:
Condições de Karush-Kuhn-Tucker

Note que aqui aparece a lagrangiana, pois a primeira condição é equivalente a:

\nabla_x l(x, u, v) = 0

O essencial para a existência de l(u, v) é que T(\bar{x}, C) = L(\bar{x}, C).

Teorema

Suponha-se que f, g_i e h_j são funções de classe \mathcal{C}^1 em uma vizinhança de \bar{x}, que f e g_i são convexas e que h_j são afim-lineares. são funções de classe \mathcal{C}^1. Se existem u \in \mathbb{R}^p e v \in \mathbb{R}^q tais que:

  1. \nabla f(\bar{x}) + \sum_{i=1}^p u_i \nabla g_i(\bar{x}) + \sum_{j=1}^q v_i \nabla h_j(\bar{x}) = 0
  2. u_i \ge 0, \, \forall i = 1, \ldots, p
  3. u_i g_i(\bar{x}) = 0, \, \forall i = 1, \ldots, p
  4. g_i(\bar{x}) \le 0, \, \forall i = 1, \ldots, p
  5. h_j(\bar{x}) = 0, \, \forall j = 1, \ldots, q

A partir deste teorema são construídos alguns algoritmos.


Crystal Clear app kaddressbook.png Este módulo tem a seguinte tarefa pendente: Confrontar as informações presentes nestes últimos teoremas com algum livro. Parece estar precisando de pequenas correções.

Considere o seguinte problema:

\left\{\begin{matrix}
\min f(x)\\
g_i(x) \le 0; i = 1, \ldots, p\\
h_j(x) = 0; j = 1, \ldots, q
\end{matrix}\right.

Sabe-se que l:\mathbb{R}^n \times \mathbb{R}^p \times \mathbb{R}^q \mapsto \mathbb{R} dada por

l(x, u, v) = f(x) + \sum_{i=1}^{p} u_i g_i (x) + \sum_{j=1}^{q} v_j h_j (x).

Agora, tem-se as KKT:

Teorema

Supondo que f, g_i e h_j são de classe \mathcal{C}^1 em uma vizinhança de \bar{x} \in C e que T(\bar{x}, C) = L(\bar{x}, C), se \bar{x} é um minimizador local de f em C, então \exists \bar{u} \in \mathbb{R}^p tal que \exists \bar{v} \in \mathbb{R}^q de modo que:

  1. \nabla_x l (\bar{x}, \bar{u}, \bar{v}) = 0
  2. \bar{u} \ge 0
  3. \nabla_u l (\bar{x}, \bar{u}, \bar{v}) \le 0
  4. \bar{u}^\top \nabla_u l (\bar{x}, \bar{u}, \bar{v}) = 0
  5. \nabla_v l (\bar{x}, \bar{u}, \bar{v}) = 0

[editar] Uma inequação sobre ínfimos e supremos

Teorema

Sejam X \subset \mathbb{R}^n e Y \subset \mathbb{R}^m dois subconjuntos arbitrários e considere uma aplicação K: X \times Y \mapsto \mathbb{R}. Então

- \infty \le \sup_{y \in Y} \inf_{x \in X} K(x, y) \le \inf_{x \in X} \sup_{y \in Y} K(x, y) \le + \infty


Demonstração
Sabe-se que
K(x, y) \le \sup_{y \in Y} K(x, y) \le + \infty,

quaisquer que sejam x \in X, y \in Y. Então, calculando o ínfimo em ambos os membros (com relação aos valores de x), e considerando que K(x, y) pode ser ilimitado inferiormente, tem-se para qualquer y \in Y :

- \infty \le \inf_{x \in X} K(x, y) \le \inf_{x \in X} \sup_{y \in Y} K(x, y) \le + \infty

Assim, calculando o supremo (com relação aos valores de y), segue das desigualdades anteriores:

- \infty \le \sup_{y \in Y} \inf_{x \in X} K(x, y) \le \inf_{x \in X} \sup_{y \in Y} K(x, y) \le + \infty
Exercício

Verifique que:

- \infty \le \sup_{\begin{smallmatrix}u \ge 0\\v \in \mathbb{R}^q \end{smallmatrix}} \inf_{x \in \mathbb{R}^n} l(x, u, v) \le \inf_{x 

\in \mathbb{R}^n} \sup_{\begin{smallmatrix}u \ge 0\\v \in \mathbb{R}^q \end{smallmatrix}} l(x, u, v) \le + \infty
Demonstração
Esta demonstração é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a na versão online deste material.

Defina-se a função:

\alpha: \mathbb{R}^n \mapsto \mathbb{R} \cup \{+\infty\} , dada por
\alpha(x) = \sup_{\begin{smallmatrix}u \ge 0 \\v \in \mathbb{R}^q \end{smallmatrix}} l(x, u, v)

e também

\beta: [ 0, +\infty )^p \times \mathbb{R}^q \mapsto \mathbb{R} \cup \{-\infty\} , dada por
\beta(u, v) = \inf_{x \in \mathbb{R}^n} l(x, u, v)

Conforme o exercício anterior, tem-se então

- \infty \le \sup_{\begin{smallmatrix}u \ge 0 \\v \in \mathbb{R}^q\end{smallmatrix}} \beta(u, v) \le \inf_{x \in \mathbb{R}^n} 

\alpha(x) \le + \infty

Observando a relação entre essas funções, é natural considerar dois problemas de otimização:

(Pr)\left\{\begin{matrix}
\min \alpha(x)\\
x \in \mathbb{R}^n
\end{matrix}\right.

e

(D)\left\{\begin{matrix}
\max \beta(u, v)\\
u \ge 0\\
v \in \mathbb{R}^q
\end{matrix}\right.
Comentários
  1. As funções \alpha e \beta são conhecidas na literatura como funções em dualidade (ou mais frequentemente, funções duais);
  2. O problema (Pr) é conhecido como problema primal, enquanto (D) é chamado de problema dual;
  3. Fazendo \bar{\alpha} = \inf_{x \in \mathbb{R}^n} \alpha (x) e \bar{\beta} = \sup_{\begin{smallmatrix}u \ge 0 \\v \in \mathbb{R}^q\end{smallmatrix}} \beta(u, v), segue-se do exercício anterior que \bar{\beta} \le \bar{\alpha};
  4. A diferença \bar{\alpha} - \bar{\beta} é chamada de brecha de dualidade, ou salto de dualidade (do inglês, skip duality);


Exercício

Verifique que: \alpha (x) = \left\{\begin{matrix}
f(x) & , \text{ se } x \in C\\
+\infty & , \text{ se } x \not\in C
\end{matrix}\right.

Resolução
Se x\in C tem-se que g_i(x) \le 0; i = 1, \ldots, p e h_j(x) = 0; j = 1, \ldots, q .

Substituindo na lagrangeana tem-se que l(x, u, v)\leq f(x) e tomando u=0, se vê que o supremo dos valores l(x, u, v) é atingido, tendo f(x) como valor.

Por outro lado, se x\notin C existe um índice j tal que h_j(x)\neq 0 e/ou g_i(x)> 0. No primeiro caso, seja u=0 e

v_j=\begin{cases} 
 th_i(x),  & i=j \\
 0, & i\ne j 
\end{cases}
, onde t > 0 .

Com esta escolha, a lagrangiana será

l(x, u, v)= f(x)+t(h_j(x))^2.

Tomando o limite quando t \rightarrow \infty tem-se que l(x, u, v)\rightarrow \infty.

Para o outro caso a prova é análoga.

Exercício

Verifique que (D) consiste de maximizar uma função côncava em um poliedro. Lembre-se:

  1. Uma função f é côncava quando -f for convexa.
  2. Um poliedro é qualquer intersecção finita de semiespaços fechados.
Resolução
Primeiro será mostrado que :\beta: [ 0, +\infty )^p \times \mathbb{R}^q \mapsto \mathbb{R} \cup \{-\infty\} , dada por \beta(u, v) = \inf_{x \in \mathbb{R}^n} l(x, u, v) é uma função côncava. Para isso, basta mostrar que para cada t\in [0, 1] vale
\beta(ta+(1-t)b, tc+(1-t)d)\ge t\beta(a, c)+(1-t)\beta(b, d).

Chamando \bar{u}=ta+(1-t)b e \bar{v}=tc+(1-t)d tem-se:

\begin{align}
\beta(\bar{u}, \bar{v})& = \inf_{x \in \mathbb{R}^n} l(x, ta+(1-t)b, tc+(1-t)d)\\
             & = \inf_{x \in \mathbb{R}^n} f(x)+(ta+(1-t)b)^{\top}g(x)+(tc+(1-t)d)^{\top}h(x)\\
             & = \inf_{x \in \mathbb{R}^n} f(x)+ta^{\top}g(x)+(1-t)b^{\top}g(x)+tc^{\top}h(x)+(1-t)d^{\top}h(x)\\
             & = \inf_{x \in \mathbb{R}^n} tf(x)+(1-t)f(x) +ta^{\top}g(x)+(1-t)b^{\top}g(x)+tc^{\top}h(x)+(1-t)d^{\top}h(x)\\
             & = \inf_{x \in \mathbb{R}^n} t[f(x)+a^{\top}g(x)+c^{\top}h(x)]+(1-t)[f(x)+b^{\top}g(x)+d^{\top}h(x)]\\
             & = \inf_{x \in \mathbb{R}^n} tl(x, a, c)+(1-t)l(x, b, d)\\
             & \ge \inf_{x \in \mathbb{R}^n} tl(x, a, c)+\inf_{x \in \mathbb{R}^n} (1-t)l(x, b, d)\\
             & = t\inf_{x \in \mathbb{R}^n} l(x, a, c)+(1-t)\inf_{x \in \mathbb{R}^n} l(x, b, d)\\
             & =t\beta(a, c)+(1-t)\beta(b, d).


\end{align}

Quanto ao conjunto ser um poliedro, isso segue do fato de \mathbb{R}^q e u\ge 0 serem interseções finitas de semiespaços fechados.


[editar] Exemplo numérico: problema primal e seu dual

Considere o problema

(P) \left\{\begin{matrix}
\min x^2\\
1 \le x \le 2
\end{matrix}\right.
Resolução
Função quadrática restrita a um intervalo.svg

O problema é ilustrado na imagem a direita.

Primeiramente, é preciso identificar quais são as funções f, g_i e h_j envolvidas:

  • A função objetivo é aquela que se pretende minimizar, ou seja, f: \mathbb{R} \mapsto \mathbb{R}; f(x) = x^2;
  • Como aparecem apenas restrições de desigualdade, não há qualquer função h_j;
  • Reescrevendo as desigualdades como 1-x \le 0 e x-2 \le 0, conclúi-se que as funções g_i são:
g_1: \mathbb{R} \mapsto \mathbb{R}; g_1(x) = 1-x e g_2: \mathbb{R} \mapsto \mathbb{R}; g_2(x) = x-2.

Neste caso, a lagrangiana é:

l: \mathbb{R} \times \mathbb{R}^2 \mapsto \mathbb{R}, dado por
l(x, u_1, u_2) = x^2 + u_1 (1-x) + u_2 (x-2)

Em seguida, é preciso calcular \alpha e \beta:

\alpha: \mathbb{R}^n \mapsto \mathbb{R} \cup \{+\infty\} , é dada por
\alpha(x) = \sup_{\begin{smallmatrix}u_1 \ge 0\\u_2 \ge 0\end{smallmatrix}} [x^2 + u_1 (1-x) + u_2 (x-2) ] = x^2 + \sup_{u_1 \ge 0} [ u_1 (1-x) ] + \sup_{u_2 \ge 0} [ u_2 (x-2) ]

E, conforme um dos exercícios anteriores (ou através de uma breve análise dos supremos envolvidos na expressão anterior), tem-se:

\alpha(x) = \left\{\begin{matrix}
x^2 & , \text{ se } x \in [1, 2]\\
+\infty & , \text{ se } x \not\in [1, 2]
\end{matrix}\right.

Analogamente, tem-se:

\beta: \mathbb{R}^2 \mapsto \mathbb{R} \cup \{-\infty\} , dada por
\beta(u_1, u_2) = \inf_{x \in \mathbb{R}^n} [x^2 + u_1 (1-x) + u_2 (x-2) ]

Observe que em relação a x, tem-se uma função fortemente convexa, que portanto possui um único minimizador. Além disso, l é diferenciável, donde o seu único minimizador x_0 é dado pela equação:

\nabla_x l (x_0, u_1, u_2) = 0

ou seja,

2 x_0 + u_1 - u_2 = 0

donde

x_0 = \frac{u_2 - u_1}{2}

Portanto,

\beta(u_1, u_2) =\left( \frac{u_1 - u_2}{2}\right)^2 + u_1 \left(1- \frac{u_1 - u_2}{2}\right) + u_2 \left( \frac{u_1 - u_2}{2}-2\right)

Assim, os problemas em dualidade são:

(Pr)\left\{\begin{matrix}
\min \alpha(x)\\
x \in \mathbb{R}^n
\end{matrix}\right.

e

(D)\left\{\begin{matrix}
\max \beta(u_1, u_2)\\
u_1, u_2 \ge 0
\end{matrix}\right.

Pelo primeiro item do exercício, tem-se que a minimização do problema (Pr) é equivalente à minimização do problema original.

Por outro lado, através da lagrangiana, obtem-se além do problema primal (Pr), um problema dual (D), cuja estrutura é muito melhor que a do original (pois pelo outro exercício, consiste da maximização de uma função côncava \beta em um poliedro), e que servirá para encontrar o mínimo do problema primal (e consequentemente, do original).

[editar] Ponto de sela

Este é um conceito muito importante relacionado à lagrangiana.

Definição

Dado o problema (P) e a lagrangiana l associada à esse problema, se diz que (\bar{x}, \bar{u}, \bar{v}) é um ponto de sela de l se \bar{u} \ge 0 e:

l(\bar{x}, u, v) \le l(\bar{x}, \bar{u}, \bar{v}) \le l(x, \bar{u}, \bar{v}), \forall x \in \mathbb{R}^n, \forall u \ge 0, \forall v \in \mathbb{R}^q .
Teorema

O ponto (\bar{x}, \bar{u}, \bar{v}) é um ponto de sela de l se, e somente se:

  1. \bar{x} é solução de (Pr);
  2. (\bar{u}, \bar{v}) é solução de (D);
  3. \bar{\alpha} = \bar{\beta}.
Demonstração
Se (\bar{x}, \bar{u}, \bar{v}) é um ponto de sela de l, então \bar{u} \ge 0 e se tem:
l(\bar{x}, u, v) \le l(\bar{x}, \bar{u}, \bar{v}) \le l(x, \bar{u}, \bar{v}), \forall x \in \mathbb{R}^n, \forall u \ge 0, \forall v \in \mathbb{R}^q .

Em particular, como

\alpha (\bar{x}) = \sup_{\begin{smallmatrix}u \ge 0\\v \in \mathbb{R}^q \end{smallmatrix}} l(\bar{x}, u, v) \le l(\bar{x}, \bar{u}, \bar{v}).

segue da segunda desigualdade que

\bar{\alpha} \le \alpha(\bar{x}) \le l(\bar{x}, \bar{u}, \bar{v}) \le l(x, \bar{u}, \bar{v}), \forall x \in \mathbb{R}^n

e calculando o ínfimo sobre todos os x \in \mathbb{R} tem-se:

\bar{\alpha} \le \alpha(\bar{x}) \le l(\bar{x}, \bar{u}, \bar{v}) \le \inf_{x \in \mathbb{R}^n} l(x, \bar{u}, \bar{v}) = \beta(\bar{u}, \bar{v}) \le \bar{\beta}

Por outro lado, da inequação sobre ínfimos e supremos, tem-se \bar{\alpha} \ge \bar{\beta}, então  \bar{\alpha} = \bar{\beta}.

Logo, todos os membros das desigualdades acima coincidem, ou seja,

\bar{\alpha}(\bar{x}) =l(\bar{x}, \bar{u}, \bar{v}) = \beta(\bar{u}, \bar{v}) = \bar{\beta}

Mas da definição de \bar{\beta} tem-se:

\bar{\beta} = \sup_{\begin{smallmatrix}u \ge 0\\v \in \mathbb{R}^q\end{smallmatrix}} \beta(u, v).

Portanto, (\bar{u}, \bar{v}) é solução do problema (D) enquanto \bar{\alpha} = \inf_{x \in \mathbb{R}^n} \alpha (x) significa que \bar{x} é solução de (Pr)


Reciprocamente, supondo que \bar{x} é solução do problema (Pr), tem-se:

\alpha(\bar{x}) = \inf_{x \in \mathbb{R}^n} \alpha(x)

e admitindo que (\bar{u}, \bar{v}) é solução do problema (D), segue:

\beta(\bar{u}, \bar{v}) = \sup_{\begin{smallmatrix}u \ge 0 \\v \in \mathbb{R}^q \end{smallmatrix}} \beta(u, v)

Além disso, usando como hipótese que \bar{\alpha} = \bar{\beta}, segue da definição que:

l(\bar{x}, u, v) \le \alpha(\bar{x}) = \bar{\alpha} = \bar{\beta} = \beta(\bar{u}, \bar{v}) \le l(x, \bar{u}, \bar{v}), quaisquer que sejam x \in \mathbb{R}^n, u \ge 0, v \in \mathbb{R}^q.

Em particular, tomando u = \bar{u}, v = \bar{v} e x = \bar{x}, resulta:

l(\bar{x}, \bar{u}, \bar{v}) \le \bar{\alpha} = \bar{\beta} \le l(\bar{x}, \bar{u}, \bar{v}), ou seja,

l(\bar{x}, \bar{u}, \bar{v}) = \bar{\alpha} = \bar{\beta}

Logo, pela definição de \alpha, \beta tem-se:

l(\bar{x}, u, v) \le \alpha(\bar{x})= l(\bar{x}, \bar{u}, \bar{v}) = \beta(\bar{x}) = l(x, \bar{u}, \bar{v}), quaisquer que sejam x \in \mathbb{R}^n, u \ge 0, v \in \mathbb{R}^q.

Portanto, (\bar{x}, \bar{u}, \bar{v}) é um ponto de sela da lagrangiana l.


Wikipedia
A Wikipédia tem mais sobre este assunto:
Ponto de sela

[editar] Exemplos

Considere novamente o problema de otimização dado por:

\left\{\begin{matrix}
\min x^2\\
1 \le x \le 2
\end{matrix}\right.

Verificar se a lagrangiana associada à (P) tem ponto de sela.


Resolução
Função quadrática restrita a um intervalo.svg

O problema é ilustrado na imagem a direita, e conforme foi deduzido anteriormente, tem-se:

  • f: \mathbb{R} \mapsto \mathbb{R}; f(x) = x^2
  • g_1: \mathbb{R} \mapsto \mathbb{R}; g_1(x) = 1-x
  • g_2: \mathbb{R} \mapsto \mathbb{R}; g_2(x) = x-2
  • l: \mathbb{R} \times \mathbb{R}^2 \mapsto \mathbb{R}; l(x, u_1, u_2) = x^2 + u_1 (1-x) + u_2 (x-2)

Além disso, como foi mostrado anteriormente, vale

\alpha (x) = \left\{\begin{matrix}
x^2 & , \text{ se } x \in [1, 2]\\
+\infty & , \text{ se } x \not\in [1, 2]
\end{matrix}\right.

e também

\beta(u, v) = \left( \frac{u_1 - u_2}{2}\right)^2 + u_1 \left(1- \frac{u_1 - u_2}{2}\right) + u_2 \left( \frac{u_1 - u_2}{2}-2\right)

ou seja,

\begin{align}
\beta(u, v) & = \frac{1}{4} ( [u_1^2 + u_2^2 - 2u_1u_2] + [4u_1 - 2u_1^2 + 2 u_1 u_2] + [2u_1u_2 - 2u_2^2 - 8u_2])\\
      & = \frac{-1}{4} ( u_1^2 + u_2^2 - 2u_1u_2 - 4u_1 + 8u_2)\\
      & = \frac{-1}{4} ( (u_1 - u_2)^2 - 4(u_1 - u_2)) - u_2\\
      & = \frac{-1}{4} ( (u_1 - u_2)^2 - 4(u_1 - u_2) + 4) + 1 - u_2\\
      & = \frac{-1}{4} (u_1 - u_2 - 2)^2 + 1 - u_2
\end{align}

Agora, consideram-se os seguintes problemas em dualidade:

(Pr)\left\{\begin{matrix}
\min \alpha(x)\\
x \in \mathbb{R}^n
\end{matrix}\right.

e

(D)\left\{\begin{matrix}
\max \beta(u_1, u_2)\\
u_1, u_2 \ge 0
\end{matrix}\right.

Donde \bar{x} = 1 é solução do problema primal, e \alpha(\bar{x}) = 1.

Intuitivamente, nos pontos em que \beta(u_1, u_2) = \frac{-1}{4} (u_1 - u_2 - 2)^2 + 1 - u_2 assumir seu valor máximo, deve-se ter u_2=0, pois conforme u_2 aumenta, \beta(u_1, u_2) decresce.

Mas pela desigualdade a respeito de supremos e ínfimos, tem-se \beta(u_1, u_2) \le \alpha (\bar{x}) = 1, quaisquer que sejam u_1 \ge 0, u_2 \ge 0.

Como \beta(u_1, u_2) \le 1, e ao tomar u_1 = 2 e u_2 = 0 tem-se \beta(2, 0) = 1, segue-se que (\bar{u_1}, \bar{u_2}) = (2, 0) é uma solução do problema dual (D).

Finalmente, como \bar{\alpha}=\alpha(\bar{x}) = 1, e \bar{\beta} = \beta(\bar{u_1}, \bar{u_2}) = 1, segue que \bar{\alpha} = \bar{\beta}. Portanto, pelo teorema anterior, o ponto (\bar{x}, \bar{u_1}, \bar{u_2}) é ponto de sela da lagrangiana l.


[editar] Análise do problema

Considerando (P), (Pr), (D) e l, se (\bar{u}, \bar{v}) é uma solução do problema dual, considere o seguinte problema:

(P_{ \bar{u}, \bar{v} })\left\{\begin{matrix}
\min l(x, \bar{u}, \bar{v})\\
x \in \mathbb{R}
\end{matrix}\right.

que no caso do exemplo reduz-se a

(P_{ \bar{u}, \bar{v} })\left\{\begin{matrix}
\min x^2 + 2(1-x)\\
x \in \mathbb{R}
\end{matrix}\right.

A solução deste problema é obtida resolvendo 2x-2=0, sendo portanto x=1. Note que esta é a solução do problema primal (Pr) (e consequentemente do problema original (P)).

Será apenas uma coincidência? Ou ainda, em que situações a solução deste último problema será também solução do problema original?

A resposta será dada pelo próximo teorema. Acompanhe.

Teorema  (dualidade lagrangiana convexa)

Suponha-se que f, g_i: \mathbb{R}^n \mapsto \mathbb{R} são de classe \mathcal{C}^1 e convexas, e que as funções h_j são afim lineares e T(x, C) = L(x, C), \forall x \in C. Nestas condições:

  1. Se o problema (P) tem solução, então \bar{\alpha} = \bar{\beta} e os problemas (Pr) e (D) têm solução.
  2. Se (P) tem solução, e (\bar{u}, \bar{v}) é solução de (D), então as soluções de (Pr) são as soluções de (P_{\bar{u}, \bar{v}}), onde
(P_{ \bar{u}, \bar{v} })\left\{\begin{matrix}
\min l(x, \bar{u}, \bar{v})\\
x \in \mathbb{R}
\end{matrix}\right.
que também são as soluções de (P).

Antes da demonstração, vale a pena imaginar a seguinte aplicação da segunda parte do teorema: Se por alguma razão não se sabe resolver o problema original (P), pode-se optar por resolver (D) (um problema concavo em um poliedro), e depois resolver (P_{ \bar{u}, \bar{v} }) (que é um problema convexo sem restrições). Em outras palavras, é possível trocar um problema difícil por dois problemas mais fáceis, (D) e (P_{ \bar{u}, \bar{v} }).

Agora a demonstração do teorema:

Demonstração
(1) Como (P) tem solução, seja \bar{x} uma solução de (P). Pelo teorema de KKT (que se aplica pois as funções são de classe \mathcal{C}^1 e se tem a qualificação das restrições), segue a existência de \bar{u} \in \mathbb{R}^p e \bar{v} \in \mathbb{R}^q, tais que:
  1. \nabla_x l (\bar{x}, \bar{u}, \bar{v}) = 0
  2. \bar{u} \ge 0
  3. \nabla_u l (\bar{x}, \bar{u}, \bar{v}) \le 0
  4. \bar{u}^\top \nabla_u l (\bar{x}, \bar{u}, \bar{v}) = 0
  5. \nabla_v l (\bar{x}, \bar{u}, \bar{v}) = 0

Como f, g_i são funções convexas e h_j afim lineares, então a função lagrangiana l(x, \bar{u}, \bar{v}) é convexa em x (pela forma de l). Isto significa que:

l(x, \bar{u}, \bar{v}) \ge l(\bar{x}, \bar{u}, \bar{v}) + < \nabla_x l(\bar{x}, \bar{u}, \bar{v}), x-\bar{x}>, \forall x \in \mathbb{R}^n.

Usando (1), conclúi-se uma das duas desigualdades que definem um ponto de sela em (\bar{x}, \bar{u}, \bar{v}):

l(x, \bar{u}, \bar{v}) \ge l(\bar{x}, \bar{u}, \bar{v}).

Considerando que \bar{x} é solução de (P), tal ponto satisfaz as restrições g_i(x) \le 0; i = 1, \ldots, p e h_j(x) = 0; j = 1, \ldots, q. Então, usando l(\bar{x}, u, v) = f(\bar{x}) + \sum_{i=1}^{p} u_i g_i (\bar{x}) + \sum_{j=1}^{q} v_j h_j (\bar{x}), segue que:

l(\bar{x}, u, v) \le f(\bar{x})

pois h_j (\bar{x}) = 0 e u_i \ge 0.

Mas, por KKT, 0 = \bar{u}^\top \nabla_u l (\bar{x}, \bar{u}, \bar{v}) = \sum_{i=1}^{p} \bar{u}_i g_i (\bar{x}). Além disso, \sum_{j=1}^{q} \bar{v}_j h_j (\bar{x}) = 0.

Logo:

l(\bar{x}, \bar{u}, \bar{v}) = f(\bar{x}) + \sum_{i=1}^{p} \bar{u}_i g_i (\bar{x}) + \sum_{j=1}^{q} \bar{v}_j h_j (\bar{x}) = f(\bar{x})

Assim,

l(x, \bar{u}, \bar{v}) \ge l(\bar{x}, \bar{u}, \bar{v}) \ge l(\bar{x}, u, v)

quaisquer que sejam x \in \mathbb{R}^n, u \ge 0 e v \in \mathbb{R}^q. Mas isso implica, por definição, que (\bar{x}, \bar{u}, \bar{v}) é um ponto de sela de l.

Logo, \bar{\alpha} = \bar{\beta}, \bar{x} é solução do problema primal (Pr) e (\bar{u}, \bar{v}) é solução do problema dual (D)

(2)Seja (\bar{u}, \bar{v}) solução de (D) (que existe, pelo item 1). Sabe-se pelo item anterior que o problema primal (Pr) tem solução e que \bar{\alpha} = \bar{\beta}.

Seja \bar{x} solução de (Pr). Então, (\bar{x}, \bar{u}, \bar{v}) é ponto de sela de l, logo valem as desigualdades:

l(\bar{x}, u, v) \le l(\bar{x}, \bar{u}, \bar{v}) \le l(x, \bar{u}, \bar{v}), \forall x \in \mathbb{R}^n, \forall u \ge 0, \forall v \in \mathbb{R}^q .

Uma vez que ( \bar{u}, \bar{v} ) está fixo e a segunda desigualdade vale para qualquer x \in \mathbb{R}^n, conclui-se que \bar{x} é solução do problema

(P_{ \bar{u}, \bar{v} })\left\{\begin{matrix}
\min l(x, \bar{u}, \bar{v})\\
x \in \mathbb{R}
\end{matrix}\right.

As soluções deste problema são soluções de (Pr) e, consequentemente, do problema original.

Exercício

Verificar se existe salto de dualidade nos problemas em dualidade para o seguinte problema de minimização:

(P_1)\left\{\begin{matrix}
\min x-y^2\\
x^2 + y^2 \le 1\\
x + y \le 1
\end{matrix}\right.
Exercício

Juntamente com o problema (P_1) do exercício anterior, considere o seguinte problema:

(P_2)\left\{\begin{matrix}
\min x-y^2\\
x \ge 0\\
y \ge 0\\
x + y \le 1
\end{matrix}\right.
Os problemas são equivalentes? (no sentido de que têm as mesmas funções objetivo e o mesmo conjunto de restrições) O que acontece com relação às condições de KKT? Apesar de f(x, y) = x-y^2 não ser convexa, é válido o resultado do teorema? Para que serve KKT? É possível resolver o problema dual e usar a resposta para resolver o primal?
Exercício

Experimente escolher x \in \mathbb{R}^5 e uma transformação linear, e um poliedro (intersecção finita de semi-espaços). É difícil de resolver o problema primal. Tente resolver o dual, usando os métodos conhecidos.

A seguir será apresentada uma proposição que responde a uma pergunta deixada anteriormente:

Proposição

Considere o problema (P) a lagrangiana l e os problemas em dualidade (Pr) e (D) e (\bar{u}, \bar{v}) soluções de (D). Se \bar{x} é solução do problema

(P_{ \bar{u}, \bar{v} })\left\{\begin{matrix}
\min l(x, \bar{u}, \bar{v})\\
x \in \mathbb{R}
\end{matrix}\right. ,
\bar{x} \in C (o conjunto dos pontos que satisfazem as restrições) e g_i(\bar{x})\bar{u}_i = 0, \forall i, então \bar{x} é também uma solução do problema (P).

Tal proposição fornece um roteiro para quem precisa resolver o problema (P) relativamente difícil:

  1. Primeiramente, resolve-se (D);
  2. Depois, constrói-se o problema (P_{ \bar{u}, \bar{v} }) e encontra-se uma solução \bar{x} para o mesmo;
  3. Finalmente, se \bar{x} satisfaz as últimas condições da proposição, ele é também uma solução de (P).
Demonstração
Se (\bar{u}, \bar{v}) é uma solução de (D), \bar{\beta} = \beta (\bar{u}, \bar{v}). Então, pela definição de \beta (\bar{u}, \bar{v}), tem-se:
\bar{\beta} = \sup_{\begin{smallmatrix}u \ge 0 \\v \in \mathbb{R}^q\end{smallmatrix}} \beta(u, v)

Como \bar{x} é solução de (P_{ \bar{u}, \bar{v} }), então:

l(\bar{x}, \bar{u}, \bar{v}) \le l(x, \bar{u}, \bar{v}), \forall x \in \mathbb{R}^n

ou seja,

l(\bar{x}, \bar{u}, \bar{v}) = \min_{x \in \mathbb{R}^n} l(x, \bar{u}, \bar{v}) = \inf_{x \in \mathbb{R}^n} l(x, \bar{u}, \bar{v}) = \bar{\beta}

Portanto,

l(\bar{x}, \bar{u}, \bar{v}) = \bar{\beta} = \beta (\bar{u}, \bar{v})

Por outro lado, como x \in C, tem-se g_i(x) \le 0; i = 1, \ldots, p e h_j(x) = 0; j = 1, \ldots, q.

Por hipótese, g_i(\bar{x})\bar{u}_i = 0, então

\sum_{i=1}^{p} \bar{u}_i g_i (\bar{x}) = 0 = \sum_{j=1}^{q} \bar{v}_j h_j (\bar{x})

Logo,

l(\bar{x}, \bar{u}, \bar{v}) = f(\bar{x}) + \sum_{i=1}^{p} \bar{u}_i g_i (\bar{x}) + \sum_{j=1}^{q} \bar{v}_j h_j (\bar{x}) = f(\bar{x}).

Então para u \ge 0 tem-se que u_i g_i(\bar{x}) \le 0, \forall i e v_j h_j(\bar{x}) = 0, \forall j. Consequentemente,

l(\bar{x}, u, v) = f(\bar{x}) + \sum_{i=1}^{p} u_i g_i (\bar{x}) + \sum_{j=1}^{q} v_j h_j (\bar{x}) \le f(\bar{x})

quaisquer que sejam u \ge 0, v \in \mathbb{R}^q, ou seja,

l(\bar{x}, u, v) = f(\bar{x}) = l(\bar{x}, \bar{u}, \bar{v}) \le l(x, \bar{u}, \bar{v}) ,

quaisquer que sejam u \ge 0, v \in \mathbb{R}^q.

Portanto, (\bar{x}, \bar{u}, \bar{v}) é ponto de sela de l, e assim, \bar{x} é solução primal (solução de (Pr)), e em consequência, solução de (P).


[editar] Resumo do esquema de dualidade

Neste ponto, pode-se sintetizar a estratégia geral para a resolução de problemas de otimização utilizando esquemas de dualidade.

Considere o problema

\left\{\begin{matrix}
\min f(x)\\
g_i(x) \le 0; i = 1, \ldots, p\\
h_j(x) = 0; j = 1, \ldots, q
\end{matrix}\right.

onde f, g_i, h_j:\mathbb{R}^n \mapsto \mathbb{R} são funções de classe \mathcal{C}^1 (\mathbb{R}^n), e seja a lagrangiana l:\mathbb{R}^n \times \mathbb{R}^p \times \mathbb{R}^q \mapsto \mathbb{R} definido por

l(x, u, v) = f(x) + \sum_{i=1}^{p} u_i g_i (x) + \sum_{j=1}^{q} v_j h_j (x).

Convenciona-se que:

  • x é a variável primal;
  • (u, v) é a variável dual;

Nesse sentido, "o \mathbb{R}^p \times \mathbb{R}^q é o dual de \mathbb{R}^n" (não confundir com o significado que essa expressão teria na análise funcional. Ver nota sobre a terminologia), ou seja:

  • \mathbb{R}^n é onde "mora" a variável primal x;
  • \mathbb{R}^p \times \mathbb{R}^q é onde "mora" a variável dual (u, v);

Definem-se então as funções \alpha e \beta da seguinte maneira:

\alpha: \mathbb{R}^n \mapsto \mathbb{R} \cup \{+\infty\} , dada por
\alpha(x) = \sup_{\begin{smallmatrix}u \ge 0 \\v \in \mathbb{R}^q \end{smallmatrix}} l(x, u, v)

e

\beta: \mathbb{R}^p \times \mathbb{R}^q \mapsto \mathbb{R} \cup \{-\infty\} , dada por
\beta(u, v) = \inf_{x \in \mathbb{R}^n} l(x, u, v)

A partir destas duas funções, formulam-se os seguintes problemas duais:

(Pr)\left\{\begin{matrix}
\min \alpha(x)\\
x \in \mathbb{R}^n
\end{matrix}\right.

e

(D)\left\{\begin{matrix}
\max \beta(u, v)\\
u \ge 0\\
v \in \mathbb{R}^q
\end{matrix}\right.

É possível verificar que (Pr) equivale ao problema original (P) e que (D) consiste da maximização de uma função côncava em um poliedro (convexo).

Logo, o dual de (P) é (D).

[editar] Conclusões

Dado qualquer problema (P), seu dual (D) é um problema côncavo (isto é, a função objetivo é côncava), tal que os pontos satisfazendo o conjunto de restrições formam um poliedro convexo.

Apesar da controvérsia filosófica existente acerca do nome destes conceitos (coisa que poderia muito bem vir a ser alterada no futuro), a moral da história é que "transforma-se um problema geralmente difícil (sem estrutura) em um problema mais fácil (cheio de estrutura)".

[editar] Uma nota sobre a terminologia

Na subárea da matemática denominada Análise Funcional, quando se tem um espaço topológico X, costuma-se chamar de dual topológico ao conjunto X^* = \{ t : X \mapsto \mathbb{R}; t \text{ é linear e continua} \}.

Aparentemente, os conceitos de dual da otimização e da análise funcional não tem relação um com o outro.

Um dos primeiros a falar de dualidade (espaços duais) foi o francês Fenchel, mas foi fortemente criticado por Urruty e Lemaréchal, pois os dois conceitos de dualidade não estão relacionados. Também o francês Brezis concordou que há um problema a ser resolvido com a nomenclatura, e um dos conceitos deveria deixar de ser chamado assim.

Wikipedia
A Wikipédia tem mais sobre este assunto:
Análise funcional

Aplicações dos métodos duais

[editar] Aplicação à programação linear

Considere um problema típico da programação linear como:

\left\{\begin{matrix}
\min a^\top x + b^\top y\\
Mx + Ny \ge c\\
Px + Qy  =  d\\
x \ge 0; \, y \in \mathbb{R}^n
\end{matrix}\right.

onde são dados a \in \mathbb{R}^n, b \in \mathbb{R}^m, M_{p \times n}, N_{p \times m}, P_{q \times n}, Q_{q \times m}, d \in \mathbb{R}^p e c \in \mathbb{R}^q. Por simplicidade, pode-se ainda adotar a seguinte notação:

C = \left\{\begin{bmatrix}x\\y\end{bmatrix}; Mx + Ny \ge c,\,Px + Qy = d,\,x \ge 0 \text{ e } y \in \mathbb{R}^n\right\}

Nesta seção será mostrado como a "bonita teoria dos métodos duais" se aplica a esse tipo de problema.

Primeiramente, calcula-se a lagrangiana:

\begin{align}
l(x,y,u,v,w) & = [a^\top x + b^\top y] + u^\top[c - Mx - Ny] + v^\top[d - Px - Qy] + w^\top [-x]\\
             & = [a - M^\top u - P^\top v - w]^\top x + [b - N^\top u - Q^\top v]^\top y + [c^\top u + d^\top v]\\
\end{align}

Note que:

  • As variáveis primais são x e y;
  • As variáveis duais são u, v e w;

Agora é preciso identificar as funções \alpha e \beta correspondentes a este problema. Conforme anteriormente, tem-se:

\alpha: \mathbb{R}^n \mapsto \mathbb{R} \cup \{+\infty\} , dada por
\alpha(x)
= \sup_{\begin{smallmatrix}u \ge 0 \\v \in \mathbb{R}^q\\w \ge 0 \end{smallmatrix}} l(x,y,u,v,w)
= \left\{\begin{array}{rcl}
a^\top x + b^\top y & , & \text{ se } \begin{bmatrix}x\\y\end{bmatrix} \in C\\
            +\infty & , & \text{ em outros casos}
\end{array}\right.

e

\beta: \mathbb{R}^p \times \mathbb{R}^q \mapsto \mathbb{R} \cup \{-\infty\} , dada por
\begin{align}
\beta(u,v) & = \inf_{\begin{smallmatrix}x \in \mathbb{R}^n \\y \in \mathbb{R}^m\end{smallmatrix}} l(x,y,u,v,w)\\
           & = \inf_{\begin{smallmatrix}x \in \mathbb{R}^n \\y \in \mathbb{R}^m\end{smallmatrix}} [a - M^\top u - P^\top v - w]^\top x + [b - N^\top u - Q^\top v]^\top y + [c^\top u + d^\top v]\\
           & = [c^\top u + d^\top v] + \inf_{x \in \mathbb{R}^n} [a - M^\top u - P^\top v - w]^\top x + \inf_{y \in \mathbb{R}^m}[b - N^\top u - Q^\top v]^\top y\\
           & = \left\{\begin{array}{rcl}
               c^\top u + d^\top v & , & \text{ se } \left\{\begin{array}{rcl}
                                                   a - M^\top u - P^\top v - w & = & 0\\
                                                   b - N^\top u - Q^\top v     & = & 0\\
                                                   u\ge 0;\, w \ge 0           &   &
                                                   \end{array}\right.\\
                           -\infty & , & \text{ em outros casos}
               \end{array}\right.
\end{align}

Logo, considerando que a - M^\top u - P^\top v = w \ge 0, o problema dual consiste no seguinte:

\left\{\begin{matrix}
\max c^\top u + d^\top v\\
M^\top u + P^\top v \le a\\
N^\top u + Q^\top v  =  b\\
u \ge 0
\end{matrix}\right.
Exercício

Verificar que \bar{x} é uma solução de

(P)\left\{\begin{matrix}
\min f(x)\\
x \in C
\end{matrix}\right.
se, e somente se, \bar{x} é uma solução de
(\bar{P})\left\{\begin{matrix}
\max (-f(x))\\
x \in C
\end{matrix}\right.
Resolução
Seja \bar{x} uma solução de (P). Então, por definição, tem-se para todo x \in C que
f(\bar{x}) \le f(x)

mas isto equivale a

-f(x) \le -f(\bar{x}),

ou seja, \bar{x} é uma solução de (\bar{P}).

Exercício

Verificar que:

\inf_{x \in C} f(x) = - \sup_{x \in C} [-f(x)]
Resolução
A resolução deste exercício é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a na versão online deste material.

[editar] Exemplificando com um problema de programação linear

O seguinte problema é chamado de problema standard (padrão) de programação linear:

(PL)\left\{\begin{matrix}
\min c^\top x\\
     Ax  =  b\\
      x \ge 0\\
\end{matrix}\right.

onde são dados A_{m \times n}, b \in \mathbb{R}^m e c \in \mathbb{R}^n.

[editar] Calculando o dual de (PL)

Primeiramente,

l(x,u,v) = c^\top x + u^\top(-x) + v^\top(b-Ax)

A função \alpha não precisa ser calculada, pois já se mostrou que

\alpha (x) = \left\{\begin{matrix}
f(x) & , \text{ se } x \in C\\
+\infty & , \text{ se } x \not\in C
\end{matrix}\right.

Por outro lado, quanto à função \beta tem-se:

\beta(u,v) = \inf_{x \in \mathbb{R}^n} l(x,u,v)
\begin{align}
\beta(u,v) & = \inf_{x \in \mathbb{R}^n} l(x,u,v)\\
           & = \inf_{x \in \mathbb{R}^n} [c^\top x - u^\top x + v^\top(b-Ax)]\\
           & = b^\top v + \inf_{x \in \mathbb{R}^n} [c - u - A^\top v]^\top x\\
           & = \left\{\begin{array}{rcl}
               b^\top v & , & \text{ se } c - u - A^\top v = 0\\
                -\infty & , & \text{ em outros casos}
               \end{array}\right.
\end{align}

Logo, o problema dual é:

(D)\left\{\begin{matrix}
         \max b^\top v\\
c - u - A^\top v  =  0\\
u \ge 0;\, v \in \mathbb{R}^m\\
\end{matrix}\right.

ou ainda

(D)\left\{\begin{matrix}
     \max b^\top v\\
    A^\top v \le c\\
v \in \mathbb{R}^m\\
\end{matrix}\right.

[editar] Calculando o dual do dual de (PL)

Considere o seguinte problema:

(D)\left\{\begin{matrix}
     \max b^\top v\\
    A^\top v \le c\\
v \in \mathbb{R}^m\\
\end{matrix}\right.

que, conforme já foi mostrado em um exercício anteriormente, equivale a

(D)\left\{\begin{matrix}
     \min -b^\top v\\
    A^\top v \le c\\
v \in \mathbb{R}^m\\
\end{matrix}\right.

A lagrangiana é dado por:

l(v,y) = -b^\top v + y^\top (A^\top v - c)

Logo,

\begin{align}
\beta(y) & = \inf_{v \in \mathbb{R}^m} l(v,y)\\
         & = \inf_{v \in \mathbb{R}^m} [-b^\top v + y^\top (A^\top v - c)]\\
         & = -c^\top y + \inf_{v \in \mathbb{R}^m} [(Ay-b)^\top v]\\
         & = \left\{\begin{array}{rcl}
              -c^\top y & , & \text{ se } Ay-b = 0\\
                -\infty & , & \text{ em outros casos}
               \end{array}\right.
\end{align}

Logo, o dual de (D) é:

(DD)\left\{\begin{matrix}
\max \beta(y)\\
      y \ge 0\\
\end{matrix}\right. ,

ou seja,

(DD)\left\{\begin{matrix}
\max -c^\top y\\
      Ay  =  b\\
       y \ge 0\\
\end{matrix}\right.

que equivale a

(P)\left\{\begin{matrix}
\min c^\top x\\
     Ax  =  b\\
      x \ge 0\\
\end{matrix}\right.

[editar] Um exemplo numérico contextualizado

Considere a seguinte situação:

Um empresário que produz cerveja dispões de 240 kg de milho, 5 kg de lúpulo e 596 kg de Malta. Para produzir um barril de cerveja preta requer 2,5 kg de milho, 0,125 kg de lúpulo e 17,5 kg de malta. Enquanto que para produzir um barril de cerveja branca, se precisa de 7,5 kg de milho, 0,125 kg de lúpulo e 10 kg de malta. Por barril de cerveja branca vendido, o empresário recebe 130 reais, enquanto por um barril de cerveja preta, recebe 230 reais. Achar o modelo matemático para otimizar o ganho do empresário.
Resolução
Primeiramente, é preciso identificar quais são as variáveis do problema, e quais são os dados. Pode-se adotar a seguinte notação para as quantidades dos dois tipos de cerveja:
  • p: Indica a quantidade (em litros) de cerveja preta;
  • b: Indica a quantidade (em litros) de cerveja branca;

O ganho do empresário, que é o que se pretende maximizar, pode ser obtido pela fórmula:

g(p,b) = 230 p + 130 b

Por hora, considere que são aceitáveis os valores de p e b serem reais (não apenas inteiros positivos, mas também "números quebrados" como 3,5, 7,11, etc).

Como o estoque de cada ingrediente é limitado, tem-se restrições que devem ser consideradas. Matematicamente tais restrições podem ser expressas assim:

\left\{\begin{matrix}
  2,5 p & + &   7,5 b & \le & 240\\
0,125 p & + & 0,125 b & \le &   5\\
 17,5 p & + &    10 b & \le & 595\\
p \ge 0 &   &         &     &    \\
b \ge 0 &   &         &     &    \\
\end{matrix}\right.

Pode-se simplificar a escrita das restrições e também da função objetivo utilizando a seguinte notação:

A = \begin{bmatrix}2,5 & 7,5 \\0,125 & 0,125\\17,5 & 10\end{bmatrix},\quad x = \begin{bmatrix}p\\b\end{bmatrix},\quad b = \begin{bmatrix}240\\5\\595\end{bmatrix} \quad \text{e}\quad c = \begin{bmatrix}230\\130\end{bmatrix}

Nesses termos, o problema de otimização a ser resolvido é:

\left\{\begin{matrix}
\max c^\top x\\
     Ax \le b\\
      x \ge 0
\end{matrix}\right.

ou apenas,

\left\{\begin{matrix}
\max c^\top x\\
x \in C
\end{matrix}\right.

onde C é o conjunto definido pelo conjunto de restrições:

C = \left\{(p,b);Ax \le b \text{  e  } x \ge 0\right\}

Tal problema tem solução pois a função objetivo g(p,b) é linear (contínua) e o conjunto de restrições forma um conjunto compacto.

Para este problema, a lagrangiana é dado por

l(x,u_1,u_2) = -c^\top x + u_1^\top (Ax-b) + u_2^\top (-x)

Donde

\begin{align}
\beta(u_1,u_2) & = \inf_{x \in \mathbb{R}^2} l(x,u_1,u_2)\\
               & = \inf_{x \in \mathbb{R}^2} [-c^\top x + u_1^\top (Ax-b) + u_2^\top (-x)]\\
               & = -b^\top u_1 + \inf_{x \in \mathbb{R}^2} [(A^\top u_1 - c - u_2)^\top x]\\
               & = \left\{\begin{array}{rcl}
                   -b^\top u_1 & , & \text{ se } A^\top u_1 - c - u_2 = 0\\
                       -\infty & , & \text{ em outros casos}
                   \end{array}\right.
\end{align}

Então, o problema dual é dado por

(D)\left\{\begin{matrix}
\max -b^\top u_1\\
A^\top u_1 \ge 0\\
u_1 \ge 0
\end{matrix}\right.

que é equivalente a

(D)\left\{\begin{matrix}
\min b^\top u_1\\
A^\top u_1 \ge 0\\
u_1 \ge 0
\end{matrix}\right.

ou, escrevendo novamente em termos dos valores numéricos,

(D)\left\{\begin{matrix}
\min 240u_1 + 5u_2 + 595u_3\\
  7,5u_1 + 0,125u_2 + 10u_3 \le 130\\
2,5u_1 + 0,125u_2 + 17,5u_3 \le 230\\
u_1 \ge 0; u_2 \ge 0\\
\end{matrix}\right.


Na década de 30, 40 e 50 havia diversos livros que tratavam cada problema de programação linear individualmente, deduzindo vez após vez os seus duais, e disso extraindo certas "regras" que eram então sugeridas ao leitor na forma "se o problema for desse tipo, use tal regra, se for daquele tipo, use esta outra, e se for deste outro tipo, use esta regra". Um dos primeiros autores que começou a trabalhar os problemas sob um novo ponto de vista, mais generalizado, foi Werner Oettio (grafia?) . Seguindo-se por George Dantzig (conhecido como inventor do método simplex), Eugen Blumb (grafia?) e Jean-Pierre Crouzeix.

Crystal Clear app kaddressbook.png Este módulo tem a seguinte tarefa pendente: Identificar e corrigir a grafia correta dos nomes dos pesquisadores mostrados acima; Encontrar fontes para comprovar a informação deste parágrafo.

[editar] Aplicação à programação quadrática

Agora, o problema a considerar passa a ser

\left\{\begin{matrix}
\min \frac{1}{2}x^\top Q x + q^\top x + \alpha\\
x \in C
\end{matrix}\right.

onde C é um poliedro (interseção finita de semi-espaços), q \in \mathbb{R}^n, \alpha \in \mathbb{R} e Q_{n \times n} é uma matriz simétrica positiva definida.

Note que este problema tem solução, uma vez que o problema irrestrito correspondente tem solução (já que Q é uma matriz simétrica positiva definida, a função é limitada inferiormente, e como C é fechado, a função objetivo assume seu valor mínimo em C, por Wolfe).

Mesmo para n = 5, os problemas de programação linear já são difíceis de resolver "à mão". É preciso utilizar alguma técnica mais sofisticada.

Para dar continuidade ao exemplo, considere que o poliedro C é dado por

C = \left\{x \in \mathbb{R}^n; x \ge 0 \text{  e  } Ax = b\right\}

com A_{m\times n} e b \in \mathbb{R}^m.

Agora será aplicado o esquema de dualidade. A lagrangiana é

l: \mathbb{R}^n \times \mathbb{R}^n \times \mathbb{R}^m \mapsto \mathbb{R}
l(x,u,v) = \frac{1}{2}x^\top Q x + q^\top x + \alpha \quad + \quad u^\top (b - Ax) \quad + \quad v^\top (-x)

além disso,

\beta: \mathbb{R}^n \times \mathbb{R}^m \mapsto \mathbb{R} \cup \{-\infty\}
\beta(u,v) = \inf_{x \in \mathbb{R}^n} l(x,u,v) = \min_{x \in \mathbb{R}^n} l(x,u,v)

e a última igualdade vale pois a função é fortemente convexa.

\begin{align}
\beta(u,v) & = \inf_{x \in \mathbb{R}^n} l(x,u,v)\\
           & = \min_{x \in \mathbb{R}^n} l(x,u,v)\\
           & = \min_{x \in \mathbb{R}^n} [\frac{1}{2}x^\top Q x + q^\top x + \alpha + u^\top (b - Ax) + v^\top (-x)]\\
           & = \min_{x \in \mathbb{R}^n} [\frac{1}{2}x^\top Q x + (q - v - A^\top u)^\top x + u^\top b + \alpha]\\
\end{align}

Considerando \nabla_x l(\bar{x}, u, v) = Q \bar{x} + q - v - A^\top u = 0, se deduz que

\bar{x} = Q^{-1} (A^\top u + v - q).

Logo,

\begin{align}
\beta(u,v) & = \frac{1}{2}\left[Q^{-1} (A^\top u + v - q)\right]^\top Q \left[Q^{-1} (A^\top u + v - q)\right] + (q - v - A^\top u)^\top \left[Q^{-1} (A^\top u + v - q)\right] + u^\top b + \alpha\\
           & = \frac{1}{2}(A^\top u + v - q)^\top Q^{-1} (A^\top u + v - q) - (A^\top u + v - q)^\top Q^{-1} (A^\top u + v - q) + u^\top b + \alpha\\
           & = \frac{-1}{2}(A^\top u + v - q)^\top Q^{-1} (A^\top u + v - q) + u^\top b + \alpha
\end{align}

Observe que, sendo os autovalores de Q positivos, o mesmo vale obrigatoriamente para Q^{-1}. Assim, como a expressão de \beta envolve (-Q^{-1}), tal função é fortemente côncava (conforme já era esperado para tal função).

Baseado nestas deduções, o problema dual é

(D)\left\{\begin{matrix}
         \max \frac{-1}{2}(A^\top u + v - q)^\top Q^{-1} (A^\top u + v - q) + u^\top b + \alpha\\
u \ge 0;\, v \in \mathbb{R}^m\\
\end{matrix}\right.

ou seja,

(D)\left\{\begin{matrix}
         \min \frac{1}{2}(A^\top u + v - q)^\top Q^{-1} (A^\top u + v - q) - (u^\top b + \alpha)\\
u \ge 0;\, v \in \mathbb{R}^m\\
\end{matrix}\right.

Usualmente este tipo de problema (D) é resolvido por meio do método do gradiente projetado.

[editar] Revisão do método do gradiente projetado

O método baseia-se na seguinte proposição:

Proposição

Seja f uma função convexa em C, um conjunto convexo e fechado. Se o ponto \bar{x} \in C é tal que \bar{x} = P_C (\bar{x} - \alpha \nabla f(\bar{x})), então f(\bar{x}) \le f(x), \forall x \in C.

Prova
Esta prova é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a na versão online deste material.

[editar] Um algoritmo para o método do gradiente projetado

Este algoritmo é bastante simples.

Primeiro passo:
  Escolha x_0 \in \mathbb{R}^n e fixe \alpha > 0.  
Passo iterativo k: Enquanto x_k \not = P_C(x_k - \alpha \nabla f(x_k))
  x_{k+1} = P_C(x_k - \alpha \nabla f(x_k))

Agora, é interessante observar como se faz para projetar um ponto em C = [0, +\infty )^n \times \mathbb{R}^m.

Exercício

Dado C = [0, +\infty )^n \times \mathbb{R}^m, mostre que

P_C\left(\begin{bmatrix}x\\y\end{bmatrix}\right) = \begin{bmatrix}\max \{x,0\}\\0\end{bmatrix}
Resolução
A resolução deste exercício é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a na versão online deste material.

[editar] Exemplificando a projeção

Seja u = \begin{bmatrix}1 & -1 & 2 & -2 & 3 & -3 & 4 & -4 & 5 & -5\end{bmatrix}^\top. Então a projeção de u sobre C = [0, +\infty )^6 \times \mathbb{R}^4 é :

P_C\left(\begin{bmatrix}1 & -1 & 2 & -2 & 3 & -3 & 4 & -4 & 5 & -5\end{bmatrix}^\top\right) = \begin{bmatrix}1 & 0 & 2 & 0 & 3 & 0 & 4 & -4 & 5 & -5\end{bmatrix}^\top

Devido a essa simplicidade ao se fazer a projeção de um ponto, o método do gradiente projetado é muito eficiente para resolver o problema (D).

[editar] Exemplo concreto

Seja C = \{ (x,y); \quad x + y = 1,\quad x \ge 0, \quad y \ge 0\}.

Crystal Clear app kaddressbook.png Este módulo tem a seguinte tarefa pendente: Incluir uma ilustração deste conjunto.

Como calcular a projeção do ponto (1,2) sobre o conjunto C, P_C (1,2)?

Resolução
Seja \bar{P} a projeção de (1,2) sobre C. Então:
\|\bar{P} - (1,2)\|^2 \le \|(x,y) - (1,2)\|^2, \quad \forall (x,y) \in C

Então:

\frac{1}{2}\|\bar{P} - (1,2)\|^2 = \min_{(x,y) \in C} \frac{1}{2} \|(x,y) - (1,2)\|^2, \quad \forall (x,y) \in C

Nota: O leitor deve observar que a inclusão de \frac{1}{2} é válida, e foi feita apenas para simplificar as contas.

Logo o modelo matemático para resolver este problema é

(P)\left\{\begin{matrix}
\min \frac{1}{2}\left[(x-1)^2 + (y-2)^2\right]\\
x \ge 0; y \ge 0\\
x + y = 1
\end{matrix}\right.

A lagrangiana é dada por

l(x,y,u,v,w) = \frac{1}{2}\left[(x-1)^2 + (y-2)^2\right] - ux - vy + w(1-x-y)

Logo,

\begin{align}
\beta(u,v,w) & = \inf_{\begin{smallmatrix}x \in \mathbb{R}\\y \in \mathbb{R} \end{smallmatrix}} l(x,y,u,v,w)\\
             & = \inf_{\begin{smallmatrix}x \in \mathbb{R}\\y \in \mathbb{R} \end{smallmatrix}} \left[\frac{1}{2}\left[(x-1)^2 + (y-2)^2\right] - ux - vy + w(1-x-y)\right]\\
             & = \inf_{\begin{smallmatrix}x \in \mathbb{R}\\y \in \mathbb{R} \end{smallmatrix}} \left[\frac{1}{2}\left[(x-1)^2 + (y-2)^2\right] - (u+w)x - (v+w)y + w\right]
\end{align}

Fazendo \nabla_{(x,y)} l(\bar{x},\bar{y},u,v,w) = 0 tem-se:

\begin{bmatrix}\bar{x} - 1 - (u + w)\\\bar{y} - 2 - (v + w)\end{bmatrix}
= \begin{bmatrix}0\\0\end{bmatrix}

Donde

\begin{bmatrix}\bar{x} \\ \bar{y} \end{bmatrix} = \begin{bmatrix}u+w+1\\v+w+2\end{bmatrix}

Logo,

\begin{align}
\beta(u,v,w) & = \frac{1}{2}\left[\left((u+w+1)-1\right)^2 + \left((v+w+2)-2\right)^2\right] - (u+w)(u+w+1) - (v+w)(v+w+2) + w\\
             & = \frac{1}{2}\left[(u+w)^2 + (v+w)^2\right] - (u+w)(u+w+1) - (v+w)(v+w+2) + w\\
             & = \frac{-1}{2}\left[(u+w)^2 + (v+w)^2\right] - (u+w) - 2(v+w) + w\\
             & = \frac{-1}{2}\left[(u+w)^2 + (v+w)^2\right] - u - 2v - 2w
\end{align}

Logo, o problema dual é

(D)\left\{\begin{matrix}
         \max \beta(u,v,w)\\
u \ge 0;\, v \ge 0;\, w \in \mathbb{R}^m\\
\end{matrix}\right.

ou seja,

(D)\left\{\begin{matrix}
         \min \frac{1}{2}\left[(u+w)^2 + (v+w)^2\right] + u + 2(v + w)\\
u \ge 0;\, v \ge 0;\, w \in \mathbb{R}^m\\
\end{matrix}\right.

Agora, para a resolução deste problema dual, pode-se usar o método do gradiente projetado.

Para isso, note que o gradiente de \beta é:

\nabla \beta (u,v,w) = \begin{bmatrix}u+w+1\\v+w+2\\u+v+2w+2\end{bmatrix}
Passo 0

Seja \alpha = \frac{1}{2} e escolha x_0 = \begin{bmatrix}u_0\\v_0\\w_0\end{bmatrix}= \begin{bmatrix}0\\0\\0\end{bmatrix}.

Passo 1

x_1
= \begin{bmatrix}u_1\\v_1\\w_1\end{bmatrix}
= P_{[0, +\infty )^2 \times \mathbb{R} }\left(\begin{bmatrix}0\\0\\0\end{bmatrix} - \frac{1}{2}\begin{bmatrix}1\\2\\2\end{bmatrix}\right)
= P_{[0, +\infty )^2 \times \mathbb{R} }\left(\begin{bmatrix}-1/2\\-1\\-1\end{bmatrix}\right)
= \begin{bmatrix}0\\0\\-1\end{bmatrix}

Passo 2

x_2
= \begin{bmatrix}u_2\\v_2\\w_2\end{bmatrix}
= P_{[0, +\infty )^2 \times \mathbb{R} }\left(\begin{bmatrix}0\\0\\-1\end{bmatrix} - \frac{1}{2}\begin{bmatrix}0\\1\\0\end{bmatrix}\right)
= P_{[0, +\infty )^2 \times \mathbb{R} }\left(\begin{bmatrix}0\\-1/2\\-1\end{bmatrix}\right)
= \begin{bmatrix}0\\0\\-1\end{bmatrix}
= \begin{bmatrix}u_1\\v_1\\w_1\end{bmatrix}

Logo, a solução dual é

\begin{bmatrix}\bar{u}\\\bar{v}\\\bar{w}\end{bmatrix} = \begin{bmatrix}0\\0\\-1\end{bmatrix}

Agora, substituindo tal solução na lagrangiana, obtem-se o problema:

(P_{ \bar{u},\bar{v},\bar{w} })\left\{\begin{matrix}
\min \frac{1}{2}\left[(x-1)^2 + (y-2)^2 + x + y - 1\right]\\
x \in \mathbb{R};\,y \in \mathbb{R}
\end{matrix}\right.

que é um problema quadrático sem restrições. Neste caso, basta igualar o gradiente a zero para determinar uma solução: \begin{bmatrix}0\\0\end{bmatrix}
= \nabla_{(x,y)} l(\bar{x},\bar{y},\bar{u},\bar{v},\bar{w})
= \begin{bmatrix}(\bar{x} - 1) + 1\\(\bar{y} - 2) + 1\end{bmatrix}
= \begin{bmatrix}\bar{x}\\\bar{y} - 1\end{bmatrix}

Donde \bar{x} = 0 e \bar{y} = 1. De acordo com a teoria desenvolvida, a solução \begin{bmatrix}0\\1\end{bmatrix} do problema é também solução do problema original, pois o produto é igual a zero (ver condições da proposição).

[editar] Exercícios resolvidos

Exercício

Encontrar a solução do seguinte problema de programação linear:

\left\{\begin{matrix}
\min c^\top x\\
Ax = b\\
x \ge 0
\end{matrix}\right.
sendo
c^\top = \begin{bmatrix}2 & 5 & 2 & 1 & 1\end{bmatrix}
A = \begin{bmatrix}-1 & 2 & 1 & -1 &  0\\
                           1 & 1 & 0 &  0 & -1\end{bmatrix}
b^\top = \begin{bmatrix}1 & 1\end{bmatrix}
Resolução
Primeiramente observe que x \in \mathbb{R}^5, então não é possível obter uma interpretação geométrica do problema. Será usado o esquema de dualidade lagrangiana:
Identificar a lagrangiana l
l: \mathbb{R}^5 \times \mathbb{R}^5 \times \mathbb{R}^2 \mapsto \mathbb{R}
l(x,u,v) = c^\top x - u^\top x + v^\top (b-Ax)
Identificar a função \beta
\beta: \mathbb{R}^5 \times \mathbb{R}^5 \mapsto \mathbb{R} \cup \{-\infty\}
\begin{align}
\beta(u,v) & = \inf_{x \in \mathbb{R}^5} l(x,u,v)\\
           & = \inf_{x \in \mathbb{R}^5} [c^\top x - u^\top x + v^\top (b-Ax)]\\
           & = \inf_{x \in \mathbb{R}^5} b^\top v + (c - u - A^\top v)^\top x\\
           & = \left\{\begin{array}{rcl}
                       b^\top v & , & \text{ se } c - u - A^\top v = 0\\
                        -\infty & , & \text{ em outros casos}
                      \end{array}\right.
\end{align}
Identificar o problema dual (D)
(D)\left\{\begin{matrix}
              \max \beta(u,v)\\
u \ge 0;\, v \in \mathbb{R}^2
\end{matrix}\right.

que equivale a

(D)\left\{\begin{matrix}
              \max b^\top v\\
     c - u - A^\top v  =  0\\
                    u \ge 0
\end{matrix}\right.

ou simplesmente

(D)\left\{\begin{matrix}
 \max b^\top v\\
A^\top v \le c\\
\end{matrix}\right.

Agora, se tem um problema em \mathbb{R}^2, e portanto, pode-se ter uma interpretação geométrica para o mesmo:

Crystal Clear app kaddressbook.png Este módulo tem a seguinte tarefa pendente: Incluir ilustração do conjunto dos pontos que verificam as restrições do problema dual, bem como algumas curvas de nível da função objetivo (as mais relevantes).

As inequações que definem o conjunto viável são:

\left\{\begin{matrix}
-v_1 & + & v_2 & \le 2\\
2v_1 & + & v_2 & \le 5\\
 v_1 &   &     & \le 2\\
-v_1 &   &     & \le 1\\
-v_2 &   &     & \le 1
\end{matrix}\right.

Tal conjunto é um pentágono (irregular, mas convexo), logo o valor máximo de b^\top v é assumido quando v = \begin{bmatrix}v_1\\v_2\end{bmatrix} for um dos cinco vértices. Através de simples verificação, comprova-se que o ponto de máximo é v = \begin{bmatrix}1\\3\end{bmatrix}. Conforme a teoria, este ponto é uma solução do problema dual (D).

Como (D) é um problema diferenciável e convexo (as condições de KKT são necessárias e suficientes), o dual terá solução e as duas juntas compõe um ponto de sela de l.

Considere as condições de KKT:

  1. \nabla_x l (\bar{x}, \bar{u}, \bar{v}) = c - I\bar{u} - A^\top \bar{v} = 0
  2. \bar{u} \ge 0
  3. \nabla_u l (\bar{x}, \bar{u}, \bar{v}) = - \bar{x} \le 0, ou seja, \bar{x} \ge 0
  4. \bar{u}^\top \nabla_u l (\bar{x}, \bar{u}, \bar{v}) = \bar{u}^\top  (- \bar{x}) = 0, ou seja, \bar{x}^\top \bar{u} = 0
  5. \nabla_v l (\bar{x}, \bar{u}, \bar{v}) = b - A\bar{x} = 0, ou seja, A\bar{x} = b

Note que, a partir de 2 e 3, conclui-se que 4 só é possível quando \bar{x}_i^\top \bar{u}_i = 0,\, \forall i.

Para se obter \bar{u}, basta lembrar que em (D) se tem:

c - \bar{u} - A^\top \bar{v} = 0 (comprovando 1)

Logo, \bar{u}
= c - A^\top \bar{v}
= \begin{bmatrix}2 \\ 5 \\ 2 \\ 1 \\ 1\end{bmatrix} - \begin{bmatrix}-1 & 1\\2 & 1\\ 1 & 0\\ -1 & 0\\ 0 & -1\end{bmatrix}\begin{bmatrix}1\\3\end{bmatrix} = \begin{bmatrix}2 \\ 5 \\ 2 \\ 1 \\ 1\end{bmatrix} - \begin{bmatrix}2 \\ 5 \\ 1 \\ -1 \\ -3\end{bmatrix} = \begin{bmatrix}0 \\ 0 \\ 1 \\ 2 \\ 4\end{bmatrix} \ge 0

comprovando a propriedade (2).

Como valem as propriedades (3) e (5), tem-se:

0 = x_1u_1 + x_2u_2 + x_3u_3 + x_4u_4 + x_5u_5 = x_3 + 2x_4 + 4x_5

Donde x_3 = x_4 = x_5 = 0. Resta usar a informação da propriedade (4) para obter x_1 e x_2. O sistema A\bar{x} = b pode ser escrito como:

\left\{\begin{matrix}
-x_1 & + & 2x_2 & = 1\\
 x_1 & + &  x_2 & = 1\\
     & + & 3x_2 & = 2
\end{matrix}\right.

Logo, x_2 = \frac{2}{3} e x_1 = \frac{1}{3}.

Se \bar{x}, \bar{u} e \bar{v} satisfazem todas as propriedades, então \bar{x} = \begin{bmatrix} \frac{1}{3} & \frac{2}{3} & 0 & 0 & 0\end{bmatrix}^\top é solução do problema (P), pois todas as condições de KKT são necessárias e suficientes.

Ao resolver o problema (P), poderia ter sido escolhido Ax-b em vez de b-Ax. Será que isso influenciaria o resultado final?

Acompanhe como ficaria a resolução desta maneira:

Resolução
;Identificar a lagrangiana l
l: \mathbb{R}^5 \times \mathbb{R}^5 \times \mathbb{R}^2 \mapsto \mathbb{R}
l(x,u,v) = c^\top x - u^\top x + v^\top (Ax-b)
Identificar a função \beta
\beta: \mathbb{R}^5 \times \mathbb{R}^5 \mapsto \mathbb{R} \cup \{-\infty\}
\begin{align}
\beta(u,v) & = \inf_{x \in \mathbb{R}^5} l(x,u,v)\\
           & = \inf_{x \in \mathbb{R}^5} [c^\top x - u^\top x + v^\top (Ax-b)]\\
           & = \inf_{x \in \mathbb{R}^5} -b^\top v + (c - u + A^\top v)^\top x\\
           & = \left\{\begin{array}{rcl}
                       -b^\top v & , & \text{ se } c - u + A^\top v = 0\\
                         -\infty & , & \text{ em outros casos}
                      \end{array}\right.
\end{align}
Identificar o problema dual (D)
(D)\left\{\begin{matrix}
              \max \beta(u,v)\\
u \ge 0;\, v \in \mathbb{R}^2
\end{matrix}\right.

que equivale a

(D)\left\{\begin{matrix}
              \max -b^\top v\\
      c - u + A^\top v  =  0\\
                     u \ge 0
\end{matrix}\right.

ou simplesmente

(D)\left\{\begin{matrix}
 \min b^\top v\\
-A^\top v \le c\\
\end{matrix}\right.

Novamente, chega-se até um problema em \mathbb{R}^2, que pode ser interpretado de forma geométrica.

Crystal Clear app kaddressbook.png Este módulo tem a seguinte tarefa pendente: Incluir ilustração do conjunto dos pontos que verificam as restrições deste problema dual, bem como algumas curvas de nível da função objetivo (as mais relevantes).

As inequações que definem o conjunto viável são:

\left\{\begin{matrix}
  v_1 & - & v_2 & \le 2\\
-2v_1 & - & v_2 & \le 5\\
 -v_1 &   &     & \le 2\\
  v_1 &   &     & \le 1\\
  v_2 &   &     & \le 1
\end{matrix}\right.

Este conjunto também é um pentágono, de modo que o valor mínimo de b^\top v é assumido quando v = \begin{bmatrix}v_1\\v_2\end{bmatrix} for um dos cinco vértices. Através de simples verificação, comprova-se que o ponto de mínimo, ou seja uma solução do problema dual (D), é v = \begin{bmatrix}-1\\-3\end{bmatrix}.

Nas condições de KKT, a única mudança é na propriedade (1), que se torna:

c - \bar{u} + A^\top \bar{v} = 0

Resta ainda saber quem é \bar{u} e quem é \bar{x}:

\bar{u}
= c + A^\top \bar{v}
= \begin{bmatrix}2 \\ 5 \\ 2 \\ 1 \\ 1\end{bmatrix} + \begin{bmatrix}-1 & 1\\2 & 1\\ 1 & 0\\ -1 & 0\\ 0 & -1\end{bmatrix}\begin{bmatrix}-1\\-3\end{bmatrix}
= \begin{bmatrix}2 \\ 5 \\ 2 \\ 1 \\ 1\end{bmatrix} + \begin{bmatrix}-2 \\ -5 \\ -1 \\ 1 \\ 3\end{bmatrix} = \begin{bmatrix}0 \\ 0 \\ 1 \\ 2 \\ 4\end{bmatrix} \ge 0

como antes.

Como ao obter \bar{u} e \bar{v} chegou-se ao mesmo resultado de antes, o ponto \bar{x} continuará sendo o mesmo.

Exercício

Formule como um problema de minimização com restrições o problema de projetar ortogonalmente o ponto (-5,2) sobre o conjunto C = \{(x,y): x \ge 0, y \ge 0\}. Depois, calcule explicitamente a função lagrangiana e o problema dual.

Resolução
Matematicamente resolver esse problema é resolver


\begin{cases} 
  \text{min} \frac{1}{2}[(x+5)^2+(y-2)^2]   \\
   x\ge 0 \ y\ge 0 
\end{cases}

Definimos a lagrangeana como l(x,y,u,v)=\frac{1}{2}[(x+5)^2+(y-2)^2] +u(-x)+v(-y)

\beta(u,v)= \inf_{\begin{smallmatrix}x \in \mathbb{R} \\y \in \mathbb{R}\end{smallmatrix}} 
l(x,y,u,v)

\beta(u,v)= \inf_{\begin{smallmatrix}x \in \mathbb{R} \\y \in \mathbb{R}\end{smallmatrix}} 
[\frac{1}{2}[(x+5)^2+(y-2)^2]-ux-vy]

Como a função é quadrática , podemos calcular o gradiente e igualar a zero:

\nabla_{x,y}l(x,y,u,v)= \begin{bmatrix}
\bar{x}+5-u \\ \bar{y}-2 -v
\end{bmatrix} =\begin{bmatrix}
0 \\ 0
\end{bmatrix}

Donde, \bar{x}=u-5 e \bar{y}=v+2

Substituindo na função \beta temos:

\beta(u,v)= \frac{1}{2}[(u-5+5)^2+(v+2-2)^2]-u(u-5)-v(v+2)

\beta(u,v)= \frac{1}{2}[u^2+v^2]-u^2+5u-v^2-2v)

\beta(u,v)= -[\frac{1}{2}(u^2+v^2)-5u+2v)]

O problema dual fica então:


\begin{cases} 
  \text{max} -[\frac{1}{2}(u^2+v^2)-5u+2v]   \\
   u\ge 0 \ v\ge 0 
\end{cases}

Que é equivalente a:


\begin{cases} 
  \text{min}\ \frac{1}{2}(u^2+v^2)-5u+2v   \\
   u\ge 0 \ v\ge 0 
\end{cases}


Exercício

No exercício anterior, se u é a variável dual relacionada a variável primal x e v é a variável dual relacionada a variável primal y, então verifique se (u,v) = (10,0) é solução dual e se a lagrangiana tem pontos de sela. Em caso afirmativo, calcule um ponto de sela, caso contrário, argumente porque a lagrangiana não tem pontos de sela.

Resolução
Olhando o problema


\begin{cases} 
  \text{min}\ \frac{1}{2}(u^2+v^2)-5u+2v   \\
   u\ge 0 \ v\ge 0 
\end{cases}

observamos que é um problema com restrições.

Podemos usar as equações de KKT para resolve-lo.

Definimos a lagrangeana como

L(u,v,a,b)=\frac{1}{2}(u^2+v^2)-5u+2v-au-bv

Agora usando o teorema de KKT devemos ter:

  1. \nabla_{u,v}L(u,v,a,b)=0
  2. u\ge 0
  3. v\ge 0
  4. a\ge 0
  5. b\ge 0
  6. -au=0
  7. -bv=0

Calculando o gradiente temos

\nabla_{u,v}L(u,v,a,b)= \begin{bmatrix}
\bar{u}-5-a \\ \bar{v}+2 -b
\end{bmatrix} =\begin{bmatrix}
0 \\ 0
\end{bmatrix}

E portanto \bar{u}=a+5 e \bar{v}=b-2.

Olhando nas equações acima obtemos \bar{u}=5 , \bar{v}=0, b=2 e a=0.

Portanto (5,0) é a solução dual.

Voltando na lagrangiana do problema original e substituindo os multiplicadores de Lagrange obtemos um problema sem restrições:


\begin{cases} 
  \text{min}\ \frac{1}{2}[(x+5)^2+(y-2)^2] +5(-x)   \\
   x,y\in \mathbb{R} 
\end{cases}

Calculando o gradiente temos:

\nabla l(x,y,5,0)= \begin{bmatrix}
\bar{x}+5-5 \\ \bar{y}-2 

\end{bmatrix} =\begin{bmatrix}
0 \\ 0
\end{bmatrix}

Portanto (0,2) é a solução primal.

Logo (0,2,5,0) é um candidato a ponto de sela.

Para verificar que ele é ponto de sela, basta ver se não há salto de dualidade. Mas isso não ocorre já que o valor ótimo do primal é \frac{25}{2} e o valor ótimo do dual também é \frac{25}{2}.

Método da lagrangiana aumentada


O problema a ser resolvido é:

\left\{\begin{matrix}
\min f(x)\\
h_j(x)  =  0; j = 1, \ldots, q
\end{matrix}\right.

Sabe-se que se for aplicado o método da lagrangiana, será considerada a função:

l:\mathbb{R}^n \times \mathbb{R}^q \mapsto \mathbb{R} dada por
l(x,v) = f(x) + \sum_{j=1}^{q} v_j h_j (x).

e também

\beta:\mathbb{R}^q \mapsto \mathbb{R} \cup \{-\infty\} , dada por
\beta(v) = \inf_{x \in \mathbb{R}^n} l(x,v)

A grande dificuldade seria saber quando o valor de \beta é finito. Uma idéia seria modificar um pouco a lagrangiana (aumentando-a, com um termo extra), da seguinte maneira:

l(x,v) = f(x) + \sum_{j=1}^{q} v_j h_j (x) + r \sum_{j=1}^{q} \left(h_j (x)\right)^2

Com isso, seria necessário garantir que a idéia de fato resolve o problema. Por este motivo, é preciso desenvolver alguns resultados teóricos. Para fazer a análise deste método, um primeiro resultado importante é o seguinte:

Lema  (Finsler-Debreu)

Seja A_{n\times n} uma matriz simétrica e B_{p\times n}. As seguintes afirmações são equivalentes:

  1. Se Bx = 0, com x \not = 0, então \langle Ax, x \rangle > 0;
  2. Existe r>0 tal que a matriz A + rB^\top B é definida positiva;
  3. Existe s>0 tal que a matriz A + rB^\top B é definida positiva para todo r>s.
Demonstração
Primeiramente, será mostrado que a afirmação (2) é equivalente à afirmação (3).

Obviamente, (3) implica em (2).

Por outro lado, supondo que se tem (2), existe algum s > 0, tal que a matriz A + sB^\top B é definida positiva. Logo, para concluir a outra implicação, basta notar que para qualquer x \not = 0 vale:

\begin{align} \langle (A + rB^\top B)x, x \rangle
                 & = \langle A x, x \rangle + r\langle B^\top B x, x \rangle\\
                 & = \langle A x, x \rangle + r\langle B x, Bx \rangle\\
                 & = \langle A x, x \rangle + r\|Bx\|^2\\
                 & > \langle A x, x \rangle + s\|Bx\|^2,\, \forall r > s\\
                 & = \langle (A + sB^\top B)x, x \rangle\\
                 & > 0
\end{align}

Assim, (2) também implica (3).

Agora, será mostrado que de (2) se conclui (1). De fato, se Bx = 0, para algum x \not = 0, então:

\langle A x, x \rangle = \langle A x, x \rangle + r\|Bx\|^2 = \langle (A + rB^\top B)x, x \rangle > 0

Finalmente, será garantido que se vale (1) então vale (3) (ver também página 25, de Izmailov & Solodov (2007)). Isso será mostrado por contradição, ou seja, supondo que vale (1) e que mesmo assim, não seja válido (3). Se disso seguir uma contradição, a implicação desejada é verdadeira.

Supondo que para todo x \not = 0 tal que Bx = 0, tem-se \langle Ax, x \rangle > 0, e que fosse falsa a afirmação (3), existiria para cada s>0 algum r>s e algum x \not =0 de modo que

\langle A x, x \rangle + r\|Bx\|^2 \le 0

Neste caso, dividindo por \|x\|^2, e denotando u = \frac{x}{\|x\|}, segue que para cada k \in \mathbb{N}, existe algum r_k > k e algum u_k \in \mathbb{R}^n, com \|u_k\| = 1, tais que

\langle A u_k, u_k \rangle + r_k\|Bu_k\|^2 \le 0

Como todos os u_k estão na esfera unitária, que é compacta, a sequência \{u_k\} tem algum ponto de acumulação, por exemplo \bar{u}. Passando o limite em k, e considerando apenas a subsequência de \{u_k\} que converge para \bar{u}, tem-se

  • r_k \to +\infty (pois r_k > k)
  • \langle A u_k, u_k \rangle \to \langle A \bar{u}, \bar{u} \rangle e
  • \|Bu_k\|^2 \to \|B\bar{u}\|^2 \not = 0

Neste caso,

\lim_{k \to +\infty} \langle A u_k, u_k \rangle + r_k\|Bu_k\|^2 \le 0

só é possível se \|B\bar{u}\|^2 = 0, ou seja, se B\bar{u}=0. Logo,

\langle A \bar{u}, \bar{u} \rangle \le 0

contradizendo a hipótese (1).


Exercício

Prove a seguinte variante do lema anterior:

Seja A_{n\times n} uma matriz simétrica e semi-definida positiva, e B_{p\times n}. As seguintes afirmações são equivalentes:
  1. Se Bx = 0, com x \not = 0, então \langle Ax, x \rangle > 0;
  2. Existe r>0 tal que a matriz A + rB^\top B é definida positiva;
  3. Para todo r>0, a matriz A + rB^\top B é definida positiva.

A partir de agora, o problema será:

\left\{\begin{matrix}
\min f(x)\\
h_j(x)  =  0; j = 1, \ldots, q
\end{matrix}\right.

onde se supõe que f,\, h_i : \mathbb{R}^n \mapsto \mathbb{R} são funções de classe \mathcal{C}^2\left( \mathbb{R}^n\right) e que para todo \lambda \in \mathbb{R}, o conjunto \{x \in \mathbb{R}^n; f(x) \le \lambda\} é compacto (em inglês costuma-se usar a expressão inf-compact para descrever tais funções).

Sabe-se que a lagrangiana associada ao problema (P) é:

l:\mathbb{R}^n \times \mathbb{R}^q \mapsto \mathbb{R} dada por
l(x,v) = f(x) + \sum_{j=1}^{q} v_j h_j (x).

e ainda, em uma notação mais sintética, considerando a função H dada por:

H:\mathbb{R}^n \mapsto \mathbb{R}^q definida por
H(x) = \begin{bmatrix} h_1(x)\\ \vdots \\ h_q(x)\end{bmatrix}

tem-se a lagrangiana expressa da seguinte maneira:

l(x,v) = f(x) + H(x)^\top v

Para o método da lagrangiana aumentada serão assumidas as seguintes hipóteses:

  1. Se \bar{x} é solução, então existe \bar{u} \in \mathbb{R}^q tal que \nabla_x l(\bar{x}, \bar{u}) = 0;
  2. Para todo d \not = 0, o jacobiano de H satizfaz:
J_H(x) d = 0 \Rightarrow d^\top H_x^2 l (\bar{x}, \bar{u}) d > 0

Note que a segunda hipótese tem exatamente a mesma forma de uma das condições que aparece no lema de Finsler-Debreu.

Definição

Dado \rho > 0, se define a lagrangiana aumentada para o problema (P) como:

l_{\rho}:\mathbb{R}^n \times \mathbb{R}^q \mapsto \mathbb{R} dada por
l_{\rho}(x,v) = f(x) + H(x)^\top v + \frac{\rho}{2}H(x)^\top H(x) = l(x,v) + \frac{\rho}{2}H(x)^\top H(x)

Observe que é justamente a aparição do termo \frac{\rho}{2}H(x)^\top H(x) sendo somado à lagrangiana que justifica o nome lagrangiana aumentada.

Esse conceito possui algumas interpretações:

Exercício

Verifique que a lagrangiana aumentada l_\rho é justamente a lagrangiana do problema (P) penalizado, ou seja, de

\left\{\begin{matrix}
\min f(x) + \frac{\rho}{2}H(x)^\top H(x)\\
H(x)  =  0
\end{matrix}\right.
Demonstração
Basta notar que dentro do conjunto viável as funções objetivos são iguais.
Exercício

Verifique que a função dual \beta_\rho (o ínfimo da lagrangiana aumentada em relação à primeira componente), dada por

\beta_\rho:\mathbb{R}^q \mapsto \mathbb{R} \cup \{-\infty\} , com
\beta_\rho(v) = \inf_{x \in \mathbb{R}^n} l_\rho(x,v)
para 0 < \rho < \delta satisfaz
\beta_\rho (v) \le \beta_\delta (v) \le f(\bar{x})
onde \bar{x} é solução de (P)
Demonstração
\begin{align}
\beta_\rho(v) & = \inf_{x \in \mathbb{R}^n} l_\rho(x,v)\\

              & = \inf_{x \in \mathbb{R}^n} [f(x)+H(x)^{\top}v+\frac{\rho}{2}H(x)^{\top}H(x)]\\

              & \leq \inf_{x \in \mathbb{R}^n} [f(x)+H(x)^{\top}v+\frac{\delta}{2}H(x)^{\top}H(x)]\\

              & = \inf_{x \in \mathbb{R}^n} l_\delta (x,v)\\

              & = \beta_\delta(v) 
\end{align}

E a segunda desigualdade

  \beta_\delta(v)= \inf_{x \in \mathbb{R}^n} l_\delta (x,v)\leq l_\delta (\bar{x},v)=f(\bar{x}) .


Exercício

Verifique que \{x \in \mathbb{R}^n; l_\rho(x,v) \le \lambda \} é compacto, qualquer que seja \lambda \in \mathbb{R} e v \in \mathbb{R}^q.

Resolução
A resolução deste exercício é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a na versão online deste material.
Proposição

Se \bar{x} é solução de (P), então existe algum \rho_0 > 0, algum \delta_0 > 0 e alguma vizinhança X_0 de \bar{x} tais que l_{\rho_0} ( \cdot, \bar{v}) é fortemente convexa com parâmetro \delta_0.

Demonstração
Conforme o lema de Finsler-Debreu, a hipótese (2), segundo a qual para todo d \not = 0, o jacobiano de H satizfaz:
J_H(x) d = 0 \Rightarrow d^\top H_x^2 l (\bar{x}, \bar{u}) d > 0

equivale a dizer que existe algum \rho_0 > 0 tal que \nabla_x l(\bar{x},\bar{v}) + \rho_0 \left(J_H(x)\right)^\top \left(J_H(x)\right) é definida positiva.

Mas a Hessiana de l_{\rho_0} é dada por \nabla_x^2 l_{\rho_0}(\bar{x},\bar{v}) = \nabla_x^2 l(\bar{x},\bar{v}) + \left(J_H(x)\right)^\top \left(J_H(x)\right). Logo, a hipótese equivale a dizer que \nabla_x^2 l_{\rho_0} (\bar{x}, \bar{v}) é definida positiva.

Considere agora a função \delta: \mathbb{R}^n \times \mathbb{R}^n \mapsto \mathbb{R} definida por

\delta(d,x) = d^\top \nabla_x^2 l_\rho (x,\bar{u}) d

Logo, para qualquer d \not = 0, tomando x = \bar{x}, tem-se

\delta(d,x) = \delta(d,\bar{x}) = d^\top \nabla_x^2 l_\rho (\bar{x},\bar{u}) d > 0

Pela continuidade da função \delta, existe uma vizinhança X_0 de \bar{x} tal que \delta(d,x) > 0 para todo ponto x \in X_0, e qualquer que seja d\not = 0. Além disso, tal vizinhança pode ser tomada aberta, convexa e suficientemente pequena para que \delta(d,x) \ge \epsilon > 0.

Assim, tomando o ínfimo, pode-se definir \delta_0 como

\delta_0 = \inf_{\begin{smallmatrix}x \in X_0\\\|d\| = 1\end{smallmatrix}} \delta(d,x)

Então

\delta\left(\frac{d}{\|d\|}, x\right) = \left(\frac{d}{\|d\|}\right)^\top \nabla_x^2 l_{\rho_0} (x,\bar{u}) \left(\frac{d}{\|d\|}\right) \ge \delta_0

quaisquer que sejam d\not=0 e x \in X_0.

Portanto,

d^\top \left[\nabla_x^2 l_{\rho_0} (x,\bar{u})\right] d \ge \delta_0 \|d\|^2

ou seja, os auto-valores de \nabla_x^2 l_{\rho_0} (x,\bar{u}) são todos positivos.

Para concluir que l_{\rho_0}(\cdot, \bar{u}) é fortemente convexa, basta recordar-se de dois fatos:

  1. Uma função f de classe \mathcal{C}^2 é convexa se, e somente se, \nabla^2 f é semidefinida positiva;
  2. Uma função f é fortemente convexa se, e somente se, existe uma constante \alpha > 0 tal que f(x) - \frac{\alpha}{2}\|x\|^2 é convexa, ou seja, d^\top \nabla^2 f(x) d - \alpha \|x\|^2 \ge 0.

Com isso, l_{\rho_0}(\cdot, \bar{u}) é fortemente convexa pois

d^\top \left[\nabla_x^2 l_{\rho_0} (x,\bar{u})\right] d - \delta_0 \|d\|^2 \ge 0

Isso significa que há um único mínimo local para tal função, e que consequentemente ele é um mínimo global. Das hipóteses 1 e 2 colocadas no início da discussão sobre a lagrangiana aumentada, segue que l_{\rho_0} é fortemente convexa em X_0.

Com essas condições, mostrou-se que em um ponto que seja solução, a lagrangiana aumentada é fortemente convexa.

Antes de apresentar o algoritmo, será fixada mais uma notação:

(P_{\rho, u})\left\{\begin{matrix}
\min l_{\rho} (x,u)\\
x \in \mathbb{R}^n
\end{matrix}\right.

[editar] Algoritmo da lagrangiana aumentada

Dados \rho > 0 e \epsilon \in (0,\rho).

Início: Tome u_0 \in \mathbb{R}^q e k=0.

Iteração: Calcule x_k, solução de (P_{\rho, u_k}).
     Se H(x_k) = 0, pare: \bar{x} = x_k.
     Senão, faça
           u_{k+1} = u_k + \epsilon H(x_k)
           k = k+1


Este é um dos algoritmos mais usados e mais eficientes para problemas de programação não linear. A garantia de convergência segue dos próximos teoremas.

Teorema

Sejam \rho_0, \delta_0 e X_0 como na proposição anterior. Se U é uma vizinhança de \bar{u}, existe algum \rho_U > \rho_0 tal que:

  1. X(\rho, u) \subset X_0, \forall u \in U
  2. X(\rho, \bar{u}) = \{ \bar{x} \} e \beta_\rho (\bar{u}) = f(\bar{x})

Observações:

  • X(\rho, u) denota as soluções do problema;
  • A igualdade \beta_\rho (\bar{u}) = f(\bar{x}) significa que não há salto de dualidade.
  • Já foi mostrado que a função é fortemente convexa em uma vizinhança. Logo, os minimizadores devem estar em tal vizinhança.
  • A prova é um pouco técnica, e usa as condições de KKT, mostrando que o cone linearizado é igual ao cone tangente.
Prova
Esta prova é deixada a cargo do leitor. Sinta-se livre para melhorar a qualidade deste texto, incluindo-a na versão online deste material.

O segundo teorema é:

Teorema

Se \rho é suficientemente grande e \delta suficientemente pequeno, então:

  1. x_k \to \bar{x}
  2. l_\rho (x_k,u_k) \to f(\bar{x})
  3. \{u_k\} é limitada

Observações:

  • A propriedade 2 praticamente segue do fato de não haver salto de dualidade.
Justificativa
Fica a cargo do leitor justificar este fato. Sinta-se livre para melhorar a qualidade deste texto, incluindo a justificativa na versão online deste material.

Com esses resultados, tem-se a garantia de que o algoritmo realmente converge para uma solução, desde que os parâmetros sejam tomados adequadamente. A questão que ainda permanece é como identificar os valores adequados de \rho e de \delta para que tal convergência ocorra.

Exercício

Argumente porque as hipóteses 1, 2 e 3 garantem que a iteração do algoritmo da Lagrangiana aumentada para problemas de minimização com restrições de igualdade tem uma única solução, sendo:

  1. Todas as funções são de classe \mathcal{C}^2 e a função objetivo tem todos os seus subníveis compactos.
  2. Se \bar{x} é solução do problema, então existe \bar{u} tal que o gradiente da Lagrangiana não aumentada com respeito a variável primal se anula em (\bar{x},\bar{u}).
  3. A hessiana da Lagrangiana não aumentada com respeito a variável primal é definida positiva sob a variedade ortogonal de todos os gradientes no ponto \bar{x} das restrições.

Estilo


Nesta página estarão indicadas as convenções adotadas neste wikilivro sobre otimização, no que diz respeito a sua formatação. Recomenda-se a leitura do mesmo, por todos que pretendem contribuir com a melhoria desde texto.

[editar] Definições

Definição

As definições aparecem em caixas com essa aparência.

Observações
  • Geralmente um termo importante aparece pela primeira vez em uma definição.
  • Devido à sua importância, é bom destacar a definição do restante do texto.
  • No momento, a forma de destacar uma definição neste wikilivro é a inclusão de bordas duplas em torno do texto da definição, como no exemplo acima. Para facilitar essa tarefa, utiliza-se a predefinição {{Definição}}.
  • O conceito que está sendo definido costuma ser colocado em negrito, sendo que o texto da explicação tem sido alinhado a esquerda.

[editar] Propriedades

Teorema

Sempre que uma propriedade importante dos objetos tratados no texto precisa ser destacada, isto deve ser feito em uma caixa como essa.

Observações
  • Os principais tipos de propriedades a ser destacados são: teoremas, proposições, corolários e pequenos lemas.
  • Para conseguir a formatação acima, utiliza-se a predefinição {{Teorema}}.

[editar] Exercícios

Exercício

Aqui vai o enunciado de um exercício.

  1. Pode-se colocar vários itens;
  2. Cada um com a numeração adequada;

[editar] Demonstrações e resolução de exercícios

Resolução
É interessante destacar a resolução de um exercício, ou a demonstração de alguma propriedade em um espaço como este.


[editar] Tarefas pendentes

Crystal Clear app kaddressbook.png Um dos autores deste material sugeriu que as tarefas pendentes, como a inclusão de figuras, sejam indicadas em caixas como esta.

Lista de símbolos

Otimização/Lista de símbolos

Índice remissivo



Nuvola apps korganizer.svg

Faltam capítulos neste índice.
Ampliando-o você ajudará a melhorar o Wikilivros.

Nesta página estão listados os conceitos abordados neste livro em ordem alfabética.

O nome de cada conceito possui um link para a página onde o mesmo é definido. Outras ocorrências importantes do conceito são indicadas pelos links numerados, logo após o link principal.


[editar] A

[editar] B

  • Brecha de dualidade
  • Busca linear

[editar] C

  • Condição
suficiente de segunda ordem
  • Condições de otimalidade
  • Condições de otimalidade
de Karush-Kuhn-Tucker
necessárias
de primeira ordem
de segunda ordem
suficientes
de segunda ordem
  • Condições de qualificação das restrições
de independência linear dos gradientes das restrições ativas
de linearidade
de Mangasarian-Fromovitz
de Slater
  • Cone
tangente
tangente linearizado
  • Conjunto
convexo
de nível
viável
  • Convexidade

[editar] D

  • Direção
de descida
viável
  • Direções
conjugadas
  • Dualidade

[editar] E

[editar] F

  • Função
coerciva
côncava
convexa
de barreira
de penalidade
exterior
interior
fortemente
côncava
convexa
objetivo
  • Funções em dualidade

[editar] G

[editar] H

  • Hessiana

[editar] I

[editar] J

[editar] K

[editar] L

  • Lagrangiana
aumentada
  • Lema
de Finsler-Debreu

[editar] M

  • Método
dos gradientes conjugados
da Lagrangiana aumentada
de barreira
de descida
de direções conjugadas
de Newton
de penalidade
interior
exterior
das regiões de confiança
do gradiente
conjugado
projetado
  • Minimizador
global
local
  • Multiplicadores de Lagrange

[editar] N

[editar] O

[editar] P

  • Penalidade
exterior
interior
  • Ponto
crítico
de mínimo
  • Problema
de minimização
de programação
linear
quadrática
dual
irrestrito
primal
  • Projeção

[editar] Q

[editar] R

  • Região de confiança
  • Regra
de Armijo
  • Restrições ativas

[editar] S

  • Sistema
de Lagrange
  • Solução

[editar] T

  • Teorema
de Wolfe

[editar] U

[editar] V

[editar] W

  • Wolfe, teorema

[editar] X

[editar] Y

[editar] Z

Bibliografia


[editar] Livros e artigos

  • Polak, E.; Ribière, G..Note sur la convergence de directions conjugées. Rev. Francaise Informat Recherche Operationelle, 3e Année 16 (1969) 35-43.

[editar] Páginas da internet

  • João Antônio de Vasconcelos. Apresentações de slides (em pdf):

[editar] Material complementar

Ferramentas pessoais
Espaços nominais

Variantes
Acções
Navegação
Projecto
Imprimir/exportar
Ferramentas