Prof. Dr. Emílio Carlos Nelli Silva 4. MÉTODOS NUMÉRICOS ...sites.poli.usp.br/d/pmr5215/a2-5215.pdf · Os métodos de programação matemática são classificados em métodos de

Otimização em Engenharia Mecânica Prof. Dr. Emílio Carlos Nelli Silva

4. MÉTODOS NUMÉRICOS APLICADOS EM PROBLEMAS DE OTIMIZAÇÃO DE ENGENHARIA

Os métodos numéricos discutidos nesse curso são baseados na chamada teoria de programação matemática como já comentados acima. 4.1. Classificação dos Algoritimos Numéricos em Otimização

Os algoritimos numéricos para solução de problemas de otimização são essencialmente classificados em métodos de programação matemática e métodos probabilísticos. Os métodos de programação matemática são classificados em métodos de programação linear, programação não-linear e métodos baseados em teoria de aproximações como programação linear seqüencial (PLS ou "SLP" em inglês) e programação quadrática seqüencial (PQS ou "SQP" em inglês). Por sua vez, os métodos de programação não-linear são classificados em métodos para solução de problemas de otimização sem restrição e com restrição. Entre os métodos probabilísticos temos os algoritimos genéticos e o "Simulated Annealing". A diferença essencial dos métodos de programação matemática para os métodos probabilísticos é que os últimos procuram encontrar o mínimo global do problema de otimização evitando os mínimos locais. Já os métodos de programação matemática fornecem um mínimo local1. Os métodos probabilísticos, como o próprio nome sugere, se utilizam de um processo de busca randômica guiados por decisões probabilísticas para obter o mínimo global. Além disso, os métodos probabilísticos são também ferramentas poderosas para problemas com variáveis discretas.

4.2. Programação Linear O método de programação linear (PL ou "LP" em inglês) se destina a solução de problemas de otimização lineares, ou seja, problemas em que a função objetivo e as restrições são funções lineares em relação às variáveis de projeto, ou seja:

Minimizar ;xaxc...xcxc)(fn

1iiinn2211 ∑

=

=+++=x

x tal que gnn2211k n1,..,j ,0xb...xbxb)(h ==+++=x

enn2211j n1,...,j ,0xd...xdxd)(g =≥+++=x

Existem uma grande gama de problemas de otimização nas diversas áreas do conhecimento (engenharia, economia, biologia, etc.) que se encaixam nessa formulação, e por essa razão existe uma vasta literatura disponível sobre programação linear (PL), bem como softwares disponíveis para a sua solução.

1 a menos que o problema somente possua um único mínimo em que nesse caso será o mínimo global.

44


Calculando as derivadas, por exemplo da função f, em relação às variáveis de projeto, temos:

∑=

==∂∂

⇒=n

1ii

iii .ctea

xfxaf

ou seja, os gradientes (ou derivadas) são constantes que não são necessariamente zeros. Isso implica que o extremo de um problema de programação linear não pode ser encontrado no interior do domínio viável (não há gradientes nulos no interior) e portanto está localizado nas fronteiras do domínio definido pelas equações de restrições. A Fig.4.2.1 ilustra esse conceito num espaço tridimensional de variáveis (x1, x2, x3), onde são representados o domínio viável limitado por um contorno denominado polítope, e as curvas de nível da função f.

x2

x3

x1

o

Ponto Ótimo

↓f

0)(gi =xpolítope

Fig.4.2.1: Solução de um problema de otimização linear no espaço tridimensional.

Já Fig.4.2.2 ilustra uma situação que pode ocorrer devido a natureza da solução do problema de otimização linear. A curva de nível de f é paralela a uma das restrições causando uma solução degenerada, ou seja, inúmeras soluções.

x3

x2

x1

o

Solução Degenerada

↓f

0)(gi =x

Fig.4.2.2: Caso em que ocorre uma solução degenerada num problema de otimização linear.

Apesar de somente resolver problemas de otimização lineares, uma aplicação importante da PL é nos algoritimos de solução de problemas de otimização não-linear baseados em métodos seqüenciais. Nesse caso, o problema de otimização não-linear é dividido numa seqüência de subproblemas de PL, sendo a solução obtida após um número suficiente de iterações. Esse método é denominado programação linear seqüencial (PLS ou "SLP") e será apresentado em detalhe na seção 4.4.

45


Para ilustrar os conceitos, vejamos a solução gráfica de um exemplo de problema de PL. Considere o problema no espaço bidimensional formulado abaixo:

Min x1, x2

tal que

2121 xx2)x,f(x +=

3x2x43x4x2

1x21xx

1x2x2 1x4

21

21

1

21

21

2

≥+≥+

≥≥+

≥+≥

onde as variáveis são x1 e x2. Na Fig.4.2.3 são plotados as curvas de nível da função objetivo e as equações de restrição do problema.

x1

x2

o

0,5

0,5

f=3,5

f=3,0f=2,5

f=2,0

Fig.4.2.3: Solução gráfica do problema de PL.

Nesse caso a solução ótima (obtida graficamente) vale x1=x2=1/2 com fmin=3/2. Note que a solução se encontra num vértice do domínio em que todas as restrições estão ativas com exceção das duas primeiras. Para verificar se um certo ponto (x1,x2) está no domínio viável, basta substituí-lo na equação e verificar se as restrições são respeitadas. Essa solução foi obtida simplesmente calculando o valor da função objetivo nos vértices do contorno do domínio viável, uma vez que já se sabe, como discutido acima, que a solução ótima será um desses pontos. Os vértices do domínio são determinados resolvendo-se a intersecção das equações de restrição. Se modificarmos a equação da função objetivo para f=c(2x1+4x2), a solução se torna degenerada (inúmeras soluções), pois essa equação é proporcional a uma das restrições do problema. A solução numérica de um problema de PL é obtida usando um método eficiente e confiável chamado Simplex. Essencialmente, esse método resolve um problema de programação linear na sua forma padrão, ou seja:

46


Min xtal que

xcTf =

0xbAx

≥=

onde c é um vetor n×1, A é uma matriz m×n e b é um vetor m×1. A princípio essa formulação parece não representar o exemplo resolvido acima, no entanto fazendo uso de variáveis auxiliares (ou "slack", em inglês), qualquer problema pode ser representado. Assim, por exemplo, considere as inegualdades abaixo:

1xx1x2x2

1x4

21

21

2

≥+≥+

≥

, elas podem ser facilmente transformadas em igualdades através da subtração de três variáveis auxiliares não negativas x3, x4, x5, assim:

1,...,5i 0,x1x-xx

1x-2x2x 1x-4x

i

521

421

32

=≥=+

=+=

, onde a condição de não negatividade está prevista na formulação padrão acima. Se as restrições fossem opostas (≤0), então as variáveis auxiliares seriam adicionadas e não subtraídas. Já restrições do tipo xi≤0 podem ser convertidas na forma padrão através de expressões do tipo:

0v,v,u,u :e v-u xe v-ux 2121222111 ≥== ou adicionando uma constante positiva de valor suficientemente alto:

xMx 11 +=

de forma que a nova variável nunca fique negativa durante a otimização. Esse tipo de artifício tem a vantagem de não aumentar a dimensão do problema, pois não cria novas variáveis, no entanto, pode ser difícil conhecer com antecedência o valor da constante M que fará com que a variável fique positiva na solução do problema. A utilização de valores muito altos de M pode resultar em mau condicionamento numérico do problema. Uma vez transformado o problema inicial de PL que se deseja resolver na forma padrão acima, o algoritimo inicia a solução do problema através do método Simplex. Note no problema da forma padrão que se m=n e todas as equações são linearmente independentes, há uma única solução para o sistema; se m>n não há solução (a menos que as equações sejam linearmente dependentes), pois o sistema de equações é inconsistente; é somente quando m<n que temos várias soluções possíveis para o problema e entre todas essas soluções o Simplex procura a solução que satisfaz a restrição de não negatividade (x≥0) e minimize a função objetivo f. Para isso, o Simplex segue o seguinte procedimento:

Inicia com uma solução viável básica; Seleciona-se m colunas linearmente independentes entre as n colunas da matriz A (lembre-se que m<n) como ilustrado abaixo.

•

47


[ ][ ]Km x mDA =

m colunas linearmenteindependentes

n

Essa matriz m×m é denominada D e é não singular. Assim podemos resolver o sistema:

1 x mm x m1 x mD

1DDD bDxbDx −=⇒=

onde xD é um vetor de variáveis independentes. Verifica-se que

⎪⎭

⎪⎬

⎫

⎪⎩

⎪⎨

⎧=

0

D

L

xx

é solução do sistema Ax=b. Tal solução é denominada solução básica, e xD é chamado de vetor de variáveis básicas. Uma solução básica, no entanto, não necessariamente satisfaz a restrição x≥0. As soluções básicas que satisfazem essa restrição são chamadas soluções básicas viáveis e demonstra-se que representam os "cantos" da polítope (ver definição acima e Fig.4.2.4). O número total de soluções básicas do sistema Ax=b é dado pelo número de possibilidades de se selecionar arbitrariamente m colunas de um número de n, ou seja, da teoria das permutações e combinações:

( )!m-nm!n!

mn

=⎟⎟⎠

⎞⎜⎜⎝

⎛

onde somente algumas soluções serão viáveis. Muda o conjunto de variáveis básicas (mas o novo conjunto deve ser viável) melhorando a função objetivo ao mesmo tempo;

•

• Repete essa operação até que não consiga mais reduzir o valor da função objetivo; Assim a idéia do algoritimo Simplex é continuamente reduzir o valor da função objetivo indo de uma solução básica viável para outra, e portanto percorrendo os extremos (ou "cantos") da polítope como mostrado na Fig.4.2.4.

x3

x1

x2o

1

0

2

3 Ponto Ótimo

polítope

Soluções básicas viáveis

Fig.4.2.4: Caminho seguido pelo método Simplex para a solução.

48


Uma descrição mais detalhada do algoritimo numérico do Simplex encontra-se nas referências da bibliografia dessa apostila. Nota-se que quanto maior o número de variáveis de projeto (n) ou de restrições (m), maior será o número de soluções básicas viáveis ( e portanto "cantos" da polítope), e portanto maior será o custo computacional (e portanto tempo) para o Simplex encontrar a solução. De fato, o número de operações aritméticas cresce exponencialmente com o número de variáveis, o que torna a utilização do Simplex inviável para problemas com grande número de restrições ou variáveis (ordem de milhares). Nesse caso, são utilizados algoritimos denominados métodos interiores, uma vez que ao invés de percorrerem os extremos da polítope, movem-se no seu interior. Essencialmente, a cada passo, esses algoritimos calculam uma direção a ser seguida no interior da polítope de forma a minimizar a função objetivo. Utilizam conceitos similares aos que serão descritos adiante nos algoritimos para a solução de problemas de otimização não-lineares com ou sem restrição. Isso torna o tempo de solução mais independente do número de variáveis e restrições. Em alguns casos a solução pode ser obtida em apenas uma iteração. Esses algoritimos chegam a ser até 50 vezes mais rápido do que o Simplex. Um método interior clássico é o algoritimo de Kamarkar (AT&T Bell Labs), descrito em detalhe na referência [1] da bibliografia. A Fig.4.2.4 ilustra o caminho seguido pelo algoritimo de Kamarkar na solução do problema anterior.

x2

x3

x1

o

Ponto Ótimo

↓f

Fig.4.2.4: Caminho seguido pelo algoritimo de Kamarkar para a solução.

No aplicativo SCILAB a rotina que resolve um problema de PL usando o método Simplex é denominada LINPRO. Essa rotina resolve um problema de PL da forma:

Min xtal que

xpTf =

( )

( )si

d

22

i

11

nm

nm

cxc

bxC

bxC

≤≤×

≤×

=

49


[ ]( )

[ ]( )di

21

di

21

mm1 ,

nmm ,

+×=

×+=

bbb

CCC

A sua sintaxe é dada por: [x,lagr,f]=linpro(p,C,b,ci,cs,mi,x0)

onde: x - vetor da solução ótima; •

•

• • • •

lagr - vetor de multiplicadores de Lagrange (dimensão n+mi+md) associado a cada uma das restrições. Se o valor de um multiplicador é zero, a restrição correspondente está inativa, caso contrário está ativa no ponto ótimo; f - valor ótimo da função; mi - número de restrições de igualdade; md - número de restrições de inegualdade; x0 - chute inicial de x; ou se “x0=v” o cálculo é iniciado num vértice do domínio; ou se “x0=g” o valor inicial é arbitrário;

A rotina que resolve um problema de PL usando o algoritimo de Kamarkar é denominada KAMARKAR. Para solução de problemas de PL com variáveis discretas é utilizado uma técnica denominada programação linear inteira que não será abordado nesse curso.

4.3. Métodos Para Solução de Problemas Não-Lineares em Otimização Nesse capítulo são descritos os métodos de solução de problemas em que a função objetivo e restrições não são lineares em relação às variáveis de projeto. Esses métodos estão divididos em métodos para a solução de problemas sem restrição e com restrição. Considerando que a grande maior parte dos problemas de otimização apresentam restrições é natural se perguntar por que estudar métodos para solução de problemas sem restrições. Basicamente existem duas razões: primeiro porque os métodos para a solução de problemas com restrições baseiam-se nos métodos para a solução de problemas sem restrições; segundo porque existem métodos que convertem o problema com restrição num problema sem restrição para depois resolvê-lo, como é o caso dos métodos de penalização descritos a seguir. 4.3.1. Problemas de Otimização sem Restrição

Ao contrário da programação linear, a cada iteração esses métodos procuram inicialmente encontrar uma direção (s0) a seguir que reduza a função objetivo. Uma vez obtida essa direção, decidem o quanto "andar" nessa direção (α). Através desse procedimento, a cada passo, um problema de encontrar n variáveis (x) é reduzido a um problema de encontrar uma variável (α), como descrito na equação a seguir:

)f()f()f( 0000 ααα =+=⇒+= sxxsxx onde x0 é o ponto inicial. O problema de encontrar α pode ser resolvido agora usando técnicas de minimização de uma função com uma variável que são de implementação mais simples. Essa etapa é denominada busca unidimensional e quando a direção s coincide com a direção de um dos eixos coordenados é denominada busca univariada.

50


Assim o procedimento seguido por esses métodos em cada passo é resumido a seguir:

Encontrar x0 e s0 que reduza a função objetivo f ; • • Encontrar α na direção de s0 que minimize f (busca unidimensional) e obter

00 sxx α+= ; Verificar a convergência, se satisfeita, pare e . Caso contrário,

e volte à etapa inicial.

*xx =1ii e 1 +==+ xx i

•

Para melhor entendermos a idéia desses métodos, podemos fazer uma analogia com a subida de uma montanha. Suponha que a forma da montanha é a função em questão e queremos atingir o seu pico máximo (pico da montanha). Num ato intuitivo, a cada passo dado, inicialmente procuramos uma direção a seguir e depois decidimos o quanto andar. O quanto andamos num passo, vai depender de quão bem comportado é o terreno da montanha no local. Por exemplo, se há um buraco adiante (ponto de mínimo) não vamos querer pisar nele ao final do passo, então devemos reduzir o passo, e escolher uma nova direção a seguir no passo seguinte, e assim por diante. Esses métodos são classificados em três tipos de acordo com a informação necessária da função objetivo: métodos de ordem zero, de primeira ordem e de segunda ordem. Nos métodos de ordem zero somente o valor da função objetivo é utilizado. São usados quando a função não é diferenciável ou quando a função é altamente não-linear, sendo difícil obter as derivadas de forma precisa. Entre um dos mais importantes métodos temos o método das direções conjugadas de "Powell", descrito abaixo. Nos métodos de primeira ordem são usados os valores da função objetivo e de suas derivadas (gradientes) em relação às variáveis de projeto. Entre os métodos mais clássicos existentes temos o método "Steepest Descent" e o método dos gradientes conjugados descritos abaixo. Já nos métodos de segunda ordem são utilizados o valores da função objetivo, de suas derivadas e também da matriz Hessiana. Entre os métodos existentes, destacam-se os métodos de Newton e Quasi-Newton, descritos a seguir. A rapidez na convergência do resultado aumenta do primeiro para o último método. A utilização de um ou outro método vai depender da informação disponível (ou a facilidade de calcular) sobre a função objetivo.

4.3.1.1.Métodos de Ordem Zero

Esses métodos são usados quando o valor da função é obtido com precisão pobre, e portanto os valores das derivadas (ou gradientes) não são confiáveis e não devem ser utilizados.

51


a) Método das Direções Conjugadas de Powell Muitos métodos para a solução de problemas sem restrição são desenvolvidos para minimizar funções quadráticas, embora a maior parte dos problemas tenham funções que não sejam quadráticas. Isso porque toda função pode ser bem aproximada por uma função quadrática próxima do mínimo. Esse é o caso do método de "Powell". Nesse método, a cada passo a função a ser minimizada é aproximada localmente por uma função quadrática, dada por:

cxbQxxx ++= TT

21)(f

Considere agora um conjunto de direções si, i=1,2… Q-conjugadas linearmente independentes, ou seja:

ji para 0jTi ≠=Qss

Pode ser mostrado que: "Se f for minimizada ao longo de cada direção s definida acima, então o mínimo de f será encontrado no (ou antes) do nésimo passo independentemente do ponto inicial, dado que erros de arredondamento não sejam acumulados", onde n é o número de variáveis. É importante que as direções sejam linearmente independentes (como definido acima) caso contrário não há convergência para o mínimo. Powell sugere um método conveniente de gerar direções Q-conjugadas linearmente independentes. No entanto, às vezes esse método pode gerar direções linearmente dependentes quando nesse caso não converge para o mínimo. A sua estratégia é baseada na seguinte propriedade (ver referências para detalhe): se x1 e x2 são dois pontos e s é uma direção especificada e, x1s corresponde ao ponto mínimo de f na linha paralela à s iniciando em x1 e x2s é o ponto mínimo de f na linha paralela à s iniciando em x2, então s e a direção (x2s-x1s) serão Q conjugadas. O algoritimo de Powell é descrito abaixo: 1) Minimize f ao longo das direções coordenadas (busca univariada), iniciando em

x0k e gerando os pontos x1

k,…, xnk onde k é o número do ciclo;

2) Após encerrar a busca univariada encontre o índice m correspondente à direção em que a função f apresenta o maior decréscimo indo de xm-1

k para xmk;

3) Calcule a direção “padrão” (soma de todos os movimentos

univariados) e determine o valor de α que minimize f tal que: ;

k0

kn

kp xxs −=

kp

k0 sxx α+=

4) Se ( ) ( )( ) ( )

21

km

k1m

1k0

k0

ffff

⎥⎦

⎤⎢⎣

⎡−

−<

−

+

xxxxα então utilize as mesmas direções para a próxima

busca univariada. Se a equação NÃO for satisfeita, então substitua a mésima direção pela direção padrão sp

k; 5) Comece a nova busca univariada com as direções obtidas no passo 4, e repita os

passos 2, 3, e 4 até que a convergência seja atingida, ou seja: ε≤−+k

1k xx . A Fig.4.3.1 ilustra o caminho seguido pelo método das direções conjugadas de "Powell" durante as iterações.

52


Fig.4.3.1: Exemplo do caminho seguido pelo algoritimo de Powell.

Como exemplo de aplicação, vamos considerar a determinação da máxima deflexão e rotação da extremidade de uma viga através da minimização da sua energia potencial total modelada usando um único elemento finito de viga cúbico descrito na Fig.4.3.2.

Fig.4.3.2: Viga engastada carregada na extremidade e seu modelo de MEF.

Considerando a formulação do elemento de viga cúbico, os deslocamentos são dados por:

( ) ( ) ( ) ( )[ ]⎪⎪⎭

⎪⎪⎬

⎫

⎪⎪⎩

⎪⎪⎨

⎧

+−−+−+−=

2

2

1

1

32323232

v

v

232231)(v

θ

θξξξξξξξξξξ ll

onde l/x=ξ . A energia potencial correspondente do modelo da viga é dado por:

20

2

2

2

3 pvdd

vd2EI

+⎟⎟⎠

⎞⎜⎜⎝

⎛=Π ∫

l

lξ

ξ

Como a viga está engastada em 0v0 11 ==⇒= θξ e substituindo v(ξ) em Π obtemos:

( ) 22222

22

23 pvv122v122EI

+−+=Π lll

θθ

Agora, definindo: ll2221

3

x;v x;EI

2f θ==Π

= a função Π acima fica:

53


1212

22

1 2xx12x4x12xf +−+= Assim, o problema inicial de elementos finitos foi reduzido a minimização da função acima.

Min 1212

22

1 2xx12x4x12xf +−+=Ou seja, o problema passa a ser um problema de otimização independente onde os conceitos de MEF foram utilizados apenas para a obtenção da função acima não sendo mais necessários na continuação da solução. No caso vamos utilizar o método de Powell na solução do problema. Solução: Todo método numérico exige um "chute inicial" para iniciar as iterações. Nesse caso será adotado: . Para efeito de comparação final a solução ótima do problema pode ser obtida por cálculo diferencial sendo igual à:

.

( ) 2)f( 2,1 10

T10 =⇒−−= xx

( ) 2/1,3/1 T* −−=x Assim sendo s1 e s2 as direções coordenadas, e seguindo o procedimento de Powell descrito acima temos para a primeira iteração:

( ) ( )

( )

21

20

20

20

10

12

1p

12

12

12

T12

11

11

2211

T11

3541,19166,131972,12

4940 e

31972,1)f( e 49/83147/157

40/49

fMin )univariada (busca 8/32

12/18/312/1

21

8/312/1

-3541,1)f( e 8/13

12/133/8

fMin )univariada (busca 2

12/1310

212/13

1,0

;9166,1)f( e 212/13

-1/12 fMin )univariada (busca )1(2)2)(1(12

)2(4112)(f2

101

21

0,1

⎥⎦

⎤⎢⎣

⎡−

−<=

=⎭⎬⎫

⎩⎨⎧

−−

=⇒=

⇒⇒⎭⎬⎫

⎩⎨⎧

+−−−

=⎭⎬⎫

⎩⎨⎧−

+⎭⎬⎫

⎩⎨⎧−−

=⇒

⇒⎭⎬⎫

⎩⎨⎧−

==⇒=⎭⎬⎫

⎩⎨⎧

−−

=⇒=

⇒⇒⎭⎬⎫

⎩⎨⎧

+−−

=⎭⎬⎫

⎩⎨⎧

+⎭⎬⎫

⎩⎨⎧

−−

=⇒=

=⎭⎬⎫

⎩⎨⎧

−−

=

⇒=⇒⇒+−+−+−−

+−++−=⇒⎭⎬⎫

⎩⎨⎧

−+−

=⎭⎬⎫

⎩⎨⎧

+⎭⎬⎫

⎩⎨⎧−−

=⇒=

α

α

αα

α

α

αα

ααα

ααα

α

xx

x

xxsxx

xs

xx

xs

A tabela abaixo descreve os resultados nos ciclos de iterações seguintes.

54


Tabela 4.3.1: Ciclos de iteração. Ciclos x1 x2 f

0 -1.0 -2.0 2.0 1 -1.083334 -2.0 1.916667 1 -1.083334 -1.625 1.354167 2 -0.895834 -1.625 0.9322967 2 -0.895834 -1.34375 0.6158854 2 -0.33334 -0.499999 -0.333333

4.3.1.2.Métodos de Primeira Ordem Os métodos de primeira ordem se utilizam da informação das derivadas da função objetivo para encontrar o ponto ótimo. Esses métodos apresentam uma taxa de convergência linear ou superlinear, ou seja, uma seqüência xk é dita superlinear convergente para x* de ordem pelo menos p se:

p*kk

*1k c xxxx −≤−+

e p=1 e ck é uma seqüência que converge para zero. Se p=1 é ck é uma constante então a convergência é linear.

a) Método "Steepest Descent" É um dos métodos mais antigos, sendo proposto em 1847. Nesse método a direção s a ser seguida é obtida resolvendo-se o problema de otimização:

Min stal que

i

n

1i i

T sxff ∑

= ∂∂

=∇ s

1T =ss

ou seja, procura-se encontrar a direção paralela ao gradiente com sinal oposto. Essa é a direção que determina a maior minimização de f. Considerando o Lagrangeano do problema, temos:

( ) ( )1f,L TT −+∇= ssss λλ e impondo as condições de estacionaridade obtemos a solução:

ff

f24f12

f02fL 22T

∇∇

−=⇒

⇒∇=⇒=∇⇒=⇒∇

−=⇒=+∇=∂∂

s

sssss

λλλ

λ

onde é a norma Euclideana. Assim, segue-se para próxima etapa que é a busca unidimensional, substituindo sxx α+=+ k1k na expressão de f. A Fig.4.3.3 ilustra o caminho seguido pelo "Steepest Descent" para uma função no espaço bidimensional. Note que a cada passo, o método segue na direção do gradiente da função.

55


Fig.4.3.3: Exemplo do caminho seguido pelo "Steepest Descent".

Consideremos a minimização de uma função quadrática para ilustrar alguns conceitos importantes do "Steepest Descent". Seja a função quadrática dada por:

c21f TT ++= xbQxx

onde Q é a matriz Hessiana (no caso simétrica). Substituindo sxx α+=+ k1k na expressão de f usando s dado acima e realizando a busca unidimensional, obtém-se o seguinte valor de α (faça como exercício):

( )Qss

sbQxT

TTk* +

−=α

O desempenho do método nesse caso vai depender do número de condicionamento da matriz Q. O número de condicionamento (CN) é dado pela razão entre o maior e o menor autovalores da matriz Q, ou seja:

min

maxCNλλ

=

Se CN é alto, o método caminha de forma lenta e seguindo o padrão "zig-zag", o que é conhecido como fenômeno de "hemstitching", como mostrado nas Fig.4.3.3. Consideremos a minimização da função para ilustrar novamente o fenômeno. A Fig.4.3.4 mostra as curvas de nível da função e o caminho "zig-zag" seguido pelo "Steepest Descent".

12

2212

1 x2x4xx12x12f ++−=

56


Fig.4.3.4: Fenômeno de "hemstitching".

Fazendo as mudanças de variáveis:

( )212

22

1212

2211 y3y61yy)y,f(y

12xy ;x

21-xy +++=⇒=⎟

⎠⎞

⎜⎝⎛=

Essa nova função apresenta curvas de nível como circunferências e agora CN=1. A Fig.4.3.5 ilustra o caminho seguido pelo método nessa nova função.

Fig.4.3.5: Convergência em uma iteração.

Ou seja, a solução é obtida em apenas uma iteração e o padrão "zig-zag" não está mais presente. Assim, demonstra-se que o "Steepest Descent" apresenta apenas uma taxa linear de convergência sem a transformação de variáveis. No entanto, como é díficil na maior parte dos casos encontrar a transformação de variáveis correta foi proposto o método de Direções Conjugadas ou de "Fletcher-Reeves" descrito a seguir.

b) Método de Direções Conjugadas ou "Fletcher-Reeves" Na primeira iteração, a direção a ser seguida nesse método (s0) segue o método "Steepest Descent". A partir da segunda iteração, a direção s1 é obtida tal que seja Q conjugada a s0, e assim por diante, ou seja:

57


1,2,...k ,0k1-kT ==Qss

onde Q é a matriz Hessiana da função quadrática f. Assim, as iterações prosseguem até que ocorra a convergência. Devido ao teorema de "Powell" das direções conjugadas (veja seção 4.3.1.1) a convergência para o mínimo de uma função quadrática é teoreticamente garantida na direção correspondente à iteração n-1 (sn-1), onde n é o número de variáveis. Para funções que não sejam quadráticas, a utilização de direções conjugadas perde o sentido pelo fato da matriz Hessiana não ser mais uma matriz de constantes. Mesmo assim, é prática comum utilizar esse algoritimo para funções que não sejam quadráticas. No entanto, para essas funções a convergência raramente ocorre depois de n ou menos iterações, o que exige que o algoritimo seja reiniciado depois de n iterações. O procedimento seguido pelo algoritimo de "Fletcher-Reeves" é brevemente descrito abaixo:

1. Calcular k1kk1k sxx ++ += α onde 1k+α é determinado tal que 0d

)(df

1k

1k =+

+

αα ;

2. ( )

( )kk1-k

T1-k

kT

kk

1-kkkkkkk

f- e

com 0k se e 0k se f

xggggg

sgsxgs

∇==

>+==−∇==

β

β

3. Se 1k+g ou ( ) ( )k1k ff xx −+ é suficiente pequeno, pare;

4. Caso contrário, se k<n vá para 1 ou recomece.

A Fig.4.3.6 ilustra o caminho seguido no método de Direções Conjugadas.

Fig.4.3.5: Exemplo do caminho seguido pelo método de direções conjugadas.

58


Considere o mesmo exemplo usado anteriormente: Min e 121

22

21 2xx12x4x12xf +−+= )2,1(T

0 −−=xSolução:

;00

)f(50,0

3334,0 0,4334

0f 0178,3

76036,18077,10961,1

e 0178,3

76036,142

4275,0

3077,16154,2

4275,0)4()2(

)3077.1()6154.2()f(

3077,16154,2

)f( e 8077,10961,1

048077,0

0fonal)unidimensi (busca )21(2)42)(21(12

)41(4)21(12)(f42

21

42

)(f128

21224)(f)2,1(

2*

22

222

122

22

1011

1111

1111

221

21111

0012

2100

0

⎭⎬⎫

⎩⎨⎧

=∇⇒⎭⎬⎫

⎩⎨⎧

−−

==⇒=⇒

⇒=⇒⎭⎬⎫

⎩⎨⎧

+⎭⎬⎫

⎩⎨⎧−−

=⎭⎬⎫

⎩⎨⎧

=⎭⎬⎫

⎩⎨⎧−

+

+⎭⎬⎫

⎩⎨⎧

−−

=⇒=+−

−+−=⇒+−∇=

=⇒⎭⎬⎫

⎩⎨⎧

−−

=∇⎭⎬⎫

⎩⎨⎧−−

=⇒=⇒

⇒=⇒−−++−−−−

++−+−−=⇒⎭⎬⎫

⎩⎨⎧−

+⎭⎬⎫

⎩⎨⎧−−

=⇒

⇒⎭⎬⎫

⎩⎨⎧−

=−∇=⇒⎭⎬⎫

⎩⎨⎧

−+−

=∇⇒−−==

xxx

x

ssx

sxx

x

xsxxxx

α

αα

ββ

α

αααα

αααα

dd

dd

xxxxT

Finalmente, como

00178,3

76036,18121224

42

≅⎭⎬⎫

⎩⎨⎧

⎥⎦

⎤⎢⎣

⎡−

−

⎭⎬⎫

⎩⎨⎧−

é verificado aque s0 e s1 são Q conjugadas.

4.3.1.3.Métodos de Segunda Ordem Esses métodos utilizam os valores da função objetivo, suas derivadas e sua matriz Hessiana. Os métodos mais antigos e clássicos são os métodos de Newton e Quasi-Newton, descritos a seguir.

a) Método de Newton Nesse método a direção s a ser seguida é obtida resolvendo-se o problema de otimização:

Min stal que

i

n

1i i

T sxff ∑

= ∂∂

=∇ s

1T =Qss onde essencialmente a norma Euclideana no método "Steepest Descent" foi substituída pela norma onde Q é a matriz Hessiana da função objetivo. 1=QssT

Considerando o Lagrangeano do problema, temos:

( ) ( )1f,L TT −+∇= Qssss λλ e impondo as condições de estacionaridade obtemos a solução:

59


⇒∇

−=⇒=+∇=∂∂ −

λλ

2f02fL 1QsQs

sff1 −∇=⇒∇−= − QsQs

e a matriz Q deve ser positiva-definida. Assim, segue-se para próxima etapa que é a busca unidimensional, substituindo na expressão de f.

)x(f k1

k1kkk1kk1k ∇−=+= −+++ Qxsxx αα

Para Q=I o método se iguala ao método "Steepest Descent". Para funções quadráticas demonstra-se que a solução é obtida em uma iteração com α=1. Demonstra-se que o método de Newton tem uma taxa de convergência quadrática, mas entre as suas sérias desvantagens está a necessidade de avaliar a matriz Hessiana e resolver o sistema de equações:

f−∇=Qs para obter a direção s. Assim, para cada iteração o método de newton exige o cálculo de n(n+1)/2 elementos da matriz simétrica Q, e n3 operações matemáticas para resolver o sistema de equações acima. Isso fez com que fossem desenvolvidos métodos denominados Quasi-Newton em que a se procura usar a informação dos gradientes da função para construir uma aproximação da matriz Hessiana como descrito a seguir.

b) Método de Quasi-Newton Considerando a expansão em série de Taylor dos gradientes de f em xk+1, temos:

{ kk1kkkkk1kk1k e )()(f)(fkkk

pyBpAyxxQxxpAy

==⇒−≅∇−∇ +++ 43421444 3444 21

onde Ak é a aproximação de Q e Bk+1 é a aproximação de Q-1. A relação Bk+1yk=pk é denominada relação da secante ou Quasi-Newton. Ak e Bk devem permanecer positivo-definidos durante as iterações. Eventualmente, . IAB =+ k1k

Assim a direção s é dada por )f(- kkk xBs ∇= e a busca unidimensional é realizada substituindo k1kk1k sxx ++ += α na expressão de f. A fórmula de aproximação e atualização mais utilizada para a matriz Bk e que fornece os melhores resultados é a BFGS dada por:

kT

k

Tkk

kT

k

kT

kk

kT

k

Tkk

1k yppp

yppyIB

ypypIB +⎥

⎦

⎤⎢⎣

⎡−⎥

⎦

⎤⎢⎣

⎡−=+

onde o algoritimo é iniciado fazendo-se B=I, quando nesse caso o método é equivalente ao "Steepest Descent". Essa fórmula que a matriz Bk permaneça simétrica e positiva-definida durante as atualizações. Para funções quadráticas, após n iterações Bn é exatamente igual à Q-1. Considere o mesmo exempo anterior, mas agora sendo resolvido pelo método de Quasi-Newton com aproximação BFGS.

Min e e B0=I 1212

22

1 2xx12x4x12xf +−+= )2,1(T0 −−=x

60


Solução: O primeiro passo é igual ao do "Steepest Descent"

{ }

12

2

221212

111

1

0T

0

T00

00000

011

5.025.025.01667.0

:que se- verificae

cia)(convergên 00

)f(

e 0.5-

0.3333-0186.37608.1

1.8077-1.0961-

0186.37608.1

3077.16154.2

10385.160225.060225.037213.0

)(f10385.160225.060225.037213.0

03698.001848.001848.000923.0

96127.01

51773.025873.088754.044354.0

96127.01

1001

1001

51773.088754.025873.044354.0

96127.01

1001

51773.088754.025873.044354.0

6923.26154.40.19230.0961-

e 96127.06923.26154.4

42

3077.16154.2

42

)f(-)f(- como e 0.19230.0961-

21

1.8077-1.0961-

3077.16154.2

)f( e 1.8077-1.0961-

−=⎥⎦

⎤⎢⎣

⎡=

⎭⎬⎫

⎩⎨⎧

≅∇

⎭⎬⎫

⎩⎨⎧

=⇒⎭⎬⎫

⎩⎨⎧

+⎭⎬⎫

⎩⎨⎧

=+=⇒⎭⎬⎫

⎩⎨⎧

=

=⎭⎬⎫

⎩⎨⎧

−−

⎥⎦

⎤⎢⎣

⎡−=∇=⇒⎥

⎦

⎤⎢⎣

⎡=

=⎥⎦

⎤⎢⎣

⎡−

−+⎟⎟

⎠

⎞⎜⎜⎝

⎛⎥⎦

⎤⎢⎣

⎡−

−−⎥

⎦

⎤⎢⎣

⎡×

×⎥⎦

⎤⎢⎣

⎡⎟⎟⎠

⎞⎜⎜⎝

⎛⎥⎦

⎤⎢⎣

⎡−

−−⎥

⎦

⎤⎢⎣

⎡=⇒

⇒⎥⎦

⎤⎢⎣

⎡−

−=−

⎭⎬⎫

⎩⎨⎧

=

=⇒⎭⎬⎫

⎩⎨⎧−

=⎭⎬⎫

⎩⎨⎧−

−⎭⎬⎫

⎩⎨⎧

−−

=

=⇒⎭⎬⎫

⎩⎨⎧−

=∇=∇=⎭⎬⎫

⎩⎨⎧

=

=⎭⎬⎫

⎩⎨⎧−−

−⎭⎬⎫

⎩⎨⎧

=⇒⎭⎬⎫

⎩⎨⎧

−−

=∇⎭⎬⎫

⎩⎨⎧

=

QB

x

xsxx

xBs

B

yp

yp

yxIxBs

pxx

αα

Portanto, como a função f é quadrática B2 (n=2) é exatamente igual à Q-1. No aplicativo SCILAB estão implementados alguns dos algoritimos acima que são chamados através da rotina OPTIM. Essa rotina resolve problemas de otimização do tipo:

Min xtal que

)(costf x

si bxb ≤≤

A sintaxe da função é: [f,[xopt,[gradopt,[work]]]]=optim(costf,[contr],x0,['algo'],[df0,[mem]],[work],[stop],['in']) onde:

costf - função a ser minimizada • • • •

xopt - vetor da solução ótima; x0 - chute inicial de x; f - valor ótimo da função;

61


contr - ‘b’, bi, bs; • •

• •

• •

algo - ‘qn’ (Quasi-Newton); ‘gc’ (gradientes conjugados); ‘nd’ (não diferenciável - método de ordem zero); mem - número de variáveis usadas na aproximação da matriz Hessiana; stop - ‘ar’ (palavra chave); ‘nap’ (máximo número de chamadas da função); ‘iter’ (máximo número de iterações permitida); ‘epsg’ (valor de corte na norma do gradiente); ‘epsf’ (valor de corte na variação de costf); ‘epsx’ (vetor com valores de corte para x); gradopt - gradiente de costf é fornecido; work - vetor para ‘restart’;

Cabem agora algumas considerações extras sobre os algoritimos estudados até aqui. A princípio esses algoritimo se aplicam somente para solução de problemas de otimização, no entanto podem também ser aplicados para a solução de sistemas de equação lineares e não-lineares, uma vez que esses problemas derivam intrinsicamente de problemas de otimização. Por exemplo, em problemas de análise estrutural não-linear, a condição é a minimização da energia potencial (Π), o que considerando a condição de estacionaridade, ocorre quando seus gradientes são nulos, resultando no sistema de equações:

0)( =Π∇ x No caso em que os problemas são posados diretamente como:

Ax=b A solução pode ser obtida através da minimização da função:

( ) ( )bAxbAxbAx −−≡= T

21Min

A vantagem de resolver esse problema de minimização ao invés de resolver Ax=b numa análise estrutural não-linear é que permite obter não somente as configurações de equilíbrio estáveis da estrutura, mas também as configurações de equilíbrio instáveis, dado que não haja convergência da minimização para um mínimo local. No caso de convergência para um mínimo local, algumas técnicas numéricas são utilizadas para forçar a convergência para o mínimo global. Além disso, o métodos descritos acima, mesmo que iterativos, se mostram mais eficientes na solução de sistemas Ax=b (formulados como minimização), com grande número de variáveis, do que os métodos diretos (como Eliminação de Gauss), sendo por isso (e pela razão acima) utilizados na maior parte dos softwares comerciais de MEF atualmente.

4.3.2. Problemas de Otimização com Restrição Esses métodos resolvem problemas de otimização com restrições. Para entendermos melhor o problema enfrentado por esses métodos, consideremos novamente a analogia com a subida ao pico da montanha usada na seção 4.3.1. Imagine que agora, durante a subida encontremos obstáculos como cercas e grandes pedras que não fazem parte da topologia do terreno em si, e não nos permitem sempre caminhar no sentido de inclinação máxima, como mostrado na figura 4.3.5a. Esses obstáculos representam as restrições e vão exigir que façamos grandes desvios pra chegar ao pico da montanha, ou eventualmente eles nem nos permitam chegar ao pico da montanha, pois esse se

62


encontra no terreno do vizinho que é limitado por uma cerca. Nesse último caso, a máxima altura que conseguiríamos subir está situada em algum ponto da cerca. As diversas situações de localização do ponto ótimo estão descritas na seção 1.4. Ou seja, o ponto de mínimo (ou de máximo) não mais está associado com a condição de que o gradiente de f deva ser nulo (estacionaridade).

Fig.4.3.5a: Ilustração sobre problemas de otimização com restrição.

Existem essencialmente, dois tipos de métodos para solução desses problemas: métodos diretos e indiretos. Os métodos diretos verificam as restrições após seguirem o procedimento do algoritimo de otimização sem restrição. Já os métodos indiretos transformam o problema com restrição em um problema sem restrição e resolvem. Para entender os algoritimos desses métodos é necessário inicialmente estudar alguns conceitos que serão descritos a seguir, como condições KKT de optimalidade, problemas convexos, ponto sela (ou "Saddle Point") e dualidade.

4.3.2.1.Condições KKT de Optimalidade

As condições KKT (Kharash Kuhn-Tucker) permitem verificar se uma solução x* obtida é um mínimo local ou não, ou seja, a condição de optimalidade. O termo mínimo local ao invés de global é usado por uma questão de generalidade. Para funções convexas o mínimo local é também o global como discutido adiante. Considerando um problema com ng restrições de inegualdade, o Lagrangeano do problema é dado por:

∑ =−−= gn

1j2

jjj )t(g)(f),(L λλ xx

onde as restrições são do tipo . As condições KKT afirmam que se x* é um mínimo local então:

0)(g j ≥x

63


ativo está não g

0 se0)(g 3. ;0

;0)(g)f( .2

n1,j 0)(g viávelé .1

j

j*

jj*

j*

n

1j

*jj

**

g*

j*

g

⇒

⇒=⇒∴=

≥

=∇−∇

=≥⇒

∑=

λλ

λ

λ

x

xx

xx

onde o superescrito * indica o ponto ótimo. A condição 1 indica que a solução seja viável, a condição 2 é a condição de estacionaridade da função Lagrangeana, com a novidade que os multiplicadores de Lagrange devem ser positivos no ponto ótimo, e a condição 3 é a condição de complementaridade já discutida na seção 3.2. A razão pela qual os multiplicadores de Lagrange devam ser positivos no ponto ótimo é entendido com ajuda da Fig.4.3.6, onde ∇g1 e ∇g2 são os gradientes das restrições ativas (perpendiculares às curvas dessas restrições) e ∇f é o gradiente da função objetivo (normal a sua curva de nível).

g1=0 g2=0S

A

11 g- ∇λ

1g-∇2g-∇

22 g- ∇λ

f-∇

Fig.4.3.6: Interpretação geométrica das condições KKT. Se A é um candidato à mínimo local então pela condição de estacionaridade:

( )2211 ggf ∇+∇−=∇− λλ

Para determinar se o ponto A é um mínimo local ou não, devemos verificar se existe uma direção s saindo de A que seja viável (não desrespeite as restrições ativas) e útil (melhore o valor da função objetivo). Se essa direção existir então o resultado pode ser "melhorado" e o ponto A não é um mínimo local, caso contrário ele será. Para uma direção ser viável deve formar um ângulo obtuso com -∇g1 e -∇g2, ou seja:

0g jT ≥∇s

onde j representa as restrições ativas em A. Para a direção ser útil devemos ter: 0fT <∇s

Pré-multiplicando a condição de estacionaridade por sT temos:

gfng

1jj

Tj

T ∑=

∇=∇ ss λ

Assim nota-se que não é possível ter uma direção útil ( ) e viável ( ) se

0fT <∇s 0g j

T ≥∇s 0j ≥λ . Ou seja, não é possível encontrar uma direção s saindo de A

64


que ao mesmo tempo não desrespeite as restrições ativas e melhore o valor de f. Existe uma situação no entanto em que é possível caminhar tangente às restrições ativas e perpendicular ao gradiente de f, nesse caso:

0gf jTT =∇=∇ ss

gradiente ⊥ tangente àrestrição

Nesse caso a verificação do efeito do movimento nessa direção somente pode ser obtido investigando-se as derivadas mais altas. Dessa forma as condições KKT são condições necessárias, mas não suficientes para a optimalidade. Assim a condição de suficiência das condições KKT fica restrita à existência de solução da equação . Note que esse sistema apresenta n variáveis e ng equações (igual ao número de restrições ativas). Ou seja, se ∇gj incluir um número de direções linearmente independentes (entre os ng ∇gj disponíveis) que é igual ao número de variáveis (n), então a única solução para a equação em questão é s=0 e nesse caso as consições KKT são também condições suficientes de optimalidade. Se ∇gj incluir um número de direções linearmente independentes (entre os ng ∇gj disponíveis) que não é igual (no caso sempre menor) ao número de variáveis (n), a equação acima apresenta solução com

0g jT =∇s

0s ≠ , e nesse caso as condições KKT não são suficientes. Note que não podemos ter num espaço n dimensional mais do que n direções linearmente independentes, embora possamos ter um número maior de restrições ativas (algumas delas seriam linearmente dependentes portanto). No último caso, a condição de suficiência exige a análise da matriz Hessiana da função Lagrangeana do problema que envolve as segundas derivadas das restrições e função objetivo. Nesse caso a condição de suficiência é que essa matriz Hessiana seja positiva-definida no subespaço tangente às restrições ativas, ou seja, 0s ≠ obtido com a solução da equação . Assim, temos: 0g j

T =∇s

( ) 0) com ativas s(restriçõe 0g que tal todopara 0L jjT2T >=∇>∇ λssss

onde representa a matriz Hessiana da função Lagrangeana. Caso a matriz Hessiana seja negativa definida o ponto será de máximo.

L2∇

Uma outra condição em que as condições KKT também são condições suficientes é quando o problema é convexo, como descrito adiante. Como exemplo, consideremos a solução do problema de otimização abaixo pela análise das condições KKT. Essa certamente não é uma forma interessante de se obter a solução do problema, no entanto é didática e nos permite fixar os conceitos acima.

Min x1,x2

tal que

321

22

31 x26x10x2xf −−+−−=

0x10g0xg

0xx10g

23

12

211

≥−=≥=

≥−=

65


Solução: Note que inicialmente é importante colocar as restrições na forma , de forma a ficar compatível com as definições acima. O Lagrangeano do problema é dado por:

0)(g j ≥x

( ) ( 2312211321

22

31 x10xxx10x26x10x2xL −−−−−−−+−−= λλλ )

e as condições de estacionaridade valem:

0xx6x4xg

xg

xf

0x10x3xg

xg

xf

311222

2

33

2

11

2

22121

1

22

1

11

1

=−+−−=∂∂

−∂∂

−∂∂

=−++−=∂∂

−∂∂

−∂∂

λλλλ

λλλλ

Além disso, a matriz Hessiana da função Lagrangiana vale:

⎥⎦

⎤⎢⎣

⎡−−−−−

=∇⇒∂∂

∂=−=

∂∂∂

−−=∂∂

−=∂∂

2 1

112

12

2

121

2

222

2

121

2

x124x6

Lxx

LxxL ;x124

xL ;x6

xL

λλ

λ

Agora, a tabela abaixo descreve todos os casos possíveis que podem ocorrer com as restrições (representadas pelos seus multiplicadores de Lagrange) na solução ótima. Um valor de λj =0 indica restrição inativa e λj ≠0 indica restrição ativa.

λ1 λ2 λ3 Caso 1 0 0 0 Caso 2 ≠0 0 ≠0 Caso 3 ≠0 0 0 Caso 4 0 ≠0 0 Caso 5 0 0 ≠0 Caso 6 0 ≠0 ≠0 Caso 7 ≠0 ≠0 0 Caso 8 ≠0 ≠0 ≠0

O sétimo e oitavo casos não são possíveis, pois não há solução para os sistemas de equação. No primeiro caso a solução das equações acima (estacionaridade mais as restrições) fornece:

17,6f0 x;826,1x 21 =⇒== e 875,5f0 x;3/2x 21 =⇒=−= em ambos os casos, a princípio as condições KKT são satisfeitas ( 0j ≥λ ), mas como o número de restrições ativas (e gradientes linearmente independentes) (0) é diferente do número de variáveis (2), devemos verificar a matriz Hessiana da função L. Nesse caso ela vale:

⎥⎦

⎤⎢⎣

⎡−−−−−

=∇21

112

x124x6

Lλ

λ

e portanto é negativa-definida no primeiro ponto indicando que é um ponto de máximo e indefinida no segundo ponto indicando que é um ponto “sela”.

66


No segundo caso, 639,3 -0,7; ;10 x;1x 3121 ==== λλ . Assim as condiçõs KKT não são satisfeitas e o ponto não é de mínimo nem de máximo. No terceiro caso, temos:

08,73f24,13

599,2x847,3x

10xx0xx6x4

0x103x-

1

2

1

21

11222

2121

−====

⇒⎪⎩

⎪⎨

⎧

==+−−

=++

λλ

λ

a princípio as condições KKT são satisfeitas, mas como o número de restrições ativas (e gradientes linearmente independentes) (1) é diferente do número de variáveis (2), devemos verificar a matriz Hessiana da função L, que vale:

⎥⎦

⎤⎢⎣

⎡−

−=∇

19,3524,1324,1308,23

L2

Novamente, é negativa-definida indicando ser um ponto de máximo. No quarto caso, temos:

99,6f103/2x0x6f100x0x

0x6x4

0103x-

221

221222

221

−==−==−====

⇒⎪⎩

⎪⎨⎧

=−−

=−+λλλ

a princípio as condições KKT são satisfeitas, mas como o número de restrições ativas (e gradientes linearmente independentes) (1) é diferente do número de variáveis (2), devemos verificar a matriz Hessiana da função L, que no caso também é negativa-definida e portanto são pontos de máximo. No quinto caso, temos:

2194f64010x826,1x0x6x4

0103x-321

3222

21 −====⇒

⎪⎩

⎪⎨⎧

=−−−

=+λ

λ

a princípio as condições KKT são satisfeitas, mas como o número de restrições ativas (e gradientes linearmente independentes) (1) é diferente do número de variáveis (2), devemos novamente verificar a matriz Hessiana da função L, que no caso também é negativa-definida e portanto é novamente um ponto de máximo. No sexto caso, temos:

-2206f ;640 ;10 ;10 x;0x 3221 ===== λλ as condições KKT são satisfeitas, e como o número de restrições ativas (e gradientes linearmente independentes) (2) é igual ao número de variáveis (2) o ponto obtido é um ponto de mínimo e portanto o ponto ótimo procurado. Na prática, os softwares de otimização resolvem o problema de otimização usando um algoritimo numérico e as condições KKT são apenas verificadas no final para solução ótima. A seguir é descrito como calcular os multiplicadores de Lagrange de forma numérica para verificar as soluções KKT.

67

Documents

Prof. Dr. Emílio Carlos Nelli Silva 4. MÉTODOS NUMÉRICOS ...sites.poli.usp.br/d/pmr5215/a2-5215.pdf · Os métodos de programação matemática são classificados em métodos de